0
$\begingroup$

Convoluting a signal (using discrete FT) for a given interval [a, b] with a Gaussian can be done by circular wrapping as shown in Numerical Recipes. I found a shortcut for circular wrapping so that the phase remains zero for an even-N data. Consequently, the peaks in the signal will not shift their positions after convolution. The shortcut is to center the Gaussian at the mean (a,b)+dt/2, where dt is a sampling interval. After that, one can do a circular shift by N/2.

Consider a hypothetical case where we have a skewed Gaussian instead of a Gaussian. Should I have the centroid at the mean (a,b)+dt/2 before circular shift by N/2?

Continuous version

enter image description here

Circularly wrapped version

enter image description here

$\endgroup$

1 Answer 1

1
$\begingroup$

The key is being even, notice that $f_{N+k} = f_{k}$, and when you do a circular wrapping of an infinite function you should have $f_{k} = \sum_{i=-\infty}^\infty f((i \cdot N + k) T_s)$, then naturally $f_{N-k} = f_{k}$, in practice you want to choose $f$ and $N$ so that $f$ vanishes before $N/2$, then you are safe to use only one term of the summation that define $f_{k}$.

This holds true for asymmetric signals, the important point is that if the convolution result has more than $N$ non-zero points then then you can only recover the samples that are less than $N$ indices from both ends of the signal.

In other words if you have a signal that is non-zero from $L \le i \le H$ then you can only recover the samples at $i$ such that $H - N < i < L + N$. If you have a response of length not greater than $N$ then you can recover all the outputs. If your response has length $N + d$ you can only recover $N - d$ samples from the output. And by circularly shifting by the inputs by different number of samples you can control which samples of the output will be usable. you can also split your operators in multiple chunks to process an infinite duration signal as in overlap-save method.

An example in python

There is a pair of useful function fftshift and ifftshift to perform these operations without unintended code obfuscation.

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import skewnorm;
a = 3

plt.figure(figsize=(6, 6))

# circularly shifted equally spaced grid
t = np.fft.ifftshift(np.linspace(-10, 10, 2048));

# filter
h = skewnorm.pdf(t, a)
# signal without circular shift
x = np.zeros((2048))
x[500:1300] = np.arange(-400, 400);
# apply filter to signal (result is circularly shifted, you can 
# recover the signal in the natural order with ifftshift)
y = np.fft.irfft(np.fft.rfft(h) * np.fft.rfft(x)) / len(t)

Filtering with a skewed gaussian

$\endgroup$
8
  • 1
    $\begingroup$ This argument applies to even function like a Gaussian. Asymmetric (skewed) Gaussian cannot be even. If I split it at the centroid, there will be more points on the right than on the left. The gist of my question is how to circularly wrap a skewed function, say, if it's centroid were centered at zero. $\endgroup$
    – AChem
    Feb 26, 2021 at 15:08
  • $\begingroup$ Yes, absolutely $\endgroup$
    – Bob
    Feb 26, 2021 at 16:13
  • $\begingroup$ So I believe the shortcut would be generating a skewed Gaussian with its centroid at mean(a,b)+dt, and then do a circular shift of N/2, with N being an even number of data points. $\endgroup$
    – AChem
    Feb 26, 2021 at 22:05
  • $\begingroup$ Yes, I added an example to illustrate what I said. $\endgroup$
    – Bob
    Feb 28, 2021 at 10:03
  • 1
    $\begingroup$ Not necessarily the centroid. You can shift it by any amount as long as you compensate the shift at the output. $\endgroup$
    – Bob
    Feb 28, 2021 at 17:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.