0
$\begingroup$

I am trying to compute the autocorrelation function of a signal for which I only know the power-spectrum.

In order to test my approach I wanted to try it out on the spectrum of $1/f^2$ noise for which the autocorrelation is known (decreasing slowly). here is my test code

import scipy as sp
import colorednoise as cn

def plot_acf(length, avg):
    fig, ax1 = plt.subplots(1,1)

# get 1/f^2 noise.
    data = cn.powerlaw_psd_gaussian(2, length)

# getting the power-spectrum. scipy.welch returns the power-spectrum, with 
# positive frequencies first and then negative frequencies.
    f, pxx = sp.welch(data, nperseg = avg, window = "hanning", return_onesided = False)

# remove the spurious imaginary part of the inverse fourier transform.
    acf = np.real_if_close(np.fft.ifft(pxx))

    ax1.plot(acf[:int(avg/2)])

But it produces plots that look like that :

enter image description here

Which is obviously wrong : the strengh of the correlation is way too weak and it decreases too fast for a $1/f^2$ noise signal.

Do you know where the poblem might come from ?

$\endgroup$
  • $\begingroup$ Welch's power spectrum estimation is not the power spectrum. I think that here lays your bug... $\endgroup$ – Gideon Genadi Kogan Feb 4 at 5:41
  • $\begingroup$ Welch's method provides an consistent estimation, so given enough sample and averaging one should be able to get a pretty good approximation of the power-spectrum. $\endgroup$ – Johncowk Feb 5 at 14:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.