0
$\begingroup$

I am looking for a python code to double integrate signals in the time domain using Fast Fourier Transformation (FFT).

The idea is to apply FFT to the signal data in the time domain. Then dividing the transformed data by FFT scaling factor (time-integration property of FFT) in the frequency domain. After that, applying Inverse-FFT to return the double integrated signal to the time domain.

$\endgroup$
3
  • 3
    $\begingroup$ Divide by $j\omega$ for integration in frequency domain $\endgroup$ Nov 10, 2022 at 2:54
  • 1
    $\begingroup$ It seems to me that time-domain integration would be more computationally efficient that frequency-domain processing. $\endgroup$ Nov 10, 2022 at 9:56
  • 1
    $\begingroup$ Uhm isn't the title self-contradictory? $\endgroup$ Nov 10, 2022 at 15:43

1 Answer 1

3
$\begingroup$

Integrators are tricky. There time continuous transfer function is

$$H(\omega) = \frac{1}{j\omega} \tag{1}$$

which unfortunately means that $H(0) \rightarrow \infty$ and that an integrator is not BIBO stable.

The impulse response of an integrator is the unit step $u(t)$.

Translating this to discrete time is difficult. It follows directly from the definition of the Discrete Fourier Transform (DFT) that both time and frequency domain signals are periodic. You need to make both frequency and time windows large enough to cover "almost all" of the impulse response and transfer function. That's easy enough in the frequency domain but impossible for a unit step, so you can't sample it without a large amount of time domain aliasing.

The most common way to implement an integrator is simply

$$y[n] = x[n]+y[n-1] \tag{2}$$

That system is marginally stable, i.e. it has a pole at $z = 1$ and infinite gain at DC, so typically it's used in conjunction with a DC blocker. That's an approximation but it's "good enough" for many applications.

Another option is to design a 1rst order lowpass filter with a cutoff frequency that's lower than your lowest frequency of interest. Above the cutoff, this behaves similar to an integrator, i.e. a slope of $-6dB$ per octave and a phase of -90 degrees.

double integrate signals in the time domain using Fast Fourier Transformation (FFT).

IMO that's a really bad idea. Frequency domain filtering is complicated and the filter can't easily be sampled without major error. You can try to sample it in the frequency domain, but it's going to look ugly in the time domain and you may also have to look into Overlap add/save methods for longer signals. These require zero padding, which is also ill defined for this type of impulse response.

$\endgroup$
7
  • $\begingroup$ Thank you so much. I have accleration data collected from accelerometer sensor and my idea is to covert the accleration to displacement data. Do you have code to implement the above that you mentioned? $\endgroup$ Nov 10, 2022 at 13:25
  • 2
    $\begingroup$ What software/language are you using? Most have functions that can create the filters and perform the filtering for you. Custom code tailored to a particular question/problem is outside the scope of SE.DSP $\endgroup$
    – Jdip
    Nov 10, 2022 at 13:56
  • 2
    $\begingroup$ @YazanAlatoom: writing the code is trivial: just implement equation $(2)$. The tricky part is tuning the integrator to your requirements and the properties of the data. There is no one-size-fits-all solution for that. Integration type algorithm amplify low frequencies A LOT, so this needs to be carefully managed depending on the signal to noise ratio of your data $\endgroup$
    – Hilmar
    Nov 10, 2022 at 14:09
  • 1
    $\begingroup$ This is a standard time domain difference equation. $x[n]$ is your input signal and $y[n]$ is your output signal. $\endgroup$
    – Hilmar
    Nov 10, 2022 at 14:48
  • 2
    $\begingroup$ @YazanAlatoom: The information at dsprelated.com/blogs-1/nf/Rick_Lyons.php?searchfor=integrators might be of some interest to you. $\endgroup$ Nov 10, 2022 at 16:26

Your Answer

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

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