I have been tasked with computing the power spectrum of a noisy signal. Specifically, I am asked to do so through first attaining the autocorrelation function. From what I have read, the PSD is simply the Fourier transform of the biased autocorrelation sequence. Given that this function is symmetric, I only need to compute the cosine transform of it. I know Mathematica has functions for the autocorrelation, but I need to do it manually.

(*---Define Parameters---*)
tau = 1; (* Duration of Sampling Period in Seconds *)
sf = 2^10; (* Sampling Frequency *)
npts = tau*sf; (* Number of Samples *) 
tau0 = (1/sf); (* Time Interval Between Samples *)
nyq = sf/2; (* Nyquist Frequency *)
mean = 0; (* Expectation Value of the data *)
sd = 10; (* Standard Deviation *)

(*---Generate Data---*)
xlst = RandomReal[NormalDistribution[mean, sd], npts];
μ = Mean[xlst];
σ = StandardDeviation[xlst];

(*---Biased AutoCorrelation---*)
bacf = 
    (1/npts)*Sum[(xlst[[n]] - μ)*(xlst[[n + lag]] - μ), {n, 1, npts - lag}], 
    {lag, 0, npts - 1}];

(*---Power Spectrum---*)
p1 = 
    PlotRange -> All, SampleRate -> sf, Frame -> True, ImageSize -> Medium, 
    PlotStyle -> Thick, PlotLabel -> Style["Periodogram", FontSize -> 18], 
    ScalingFunctions -> "dB"];

p2 = 
   ListLinePlot[FourierDCT[bacf, 2], 
     PlotRange -> All, Frame -> True, ImageSize -> Medium, FontSize -> 18, 
     PlotLabel -> Style["PSD (ACF Transform)"]];

GraphicsRow[{p1, p2}, ImageSize -> Full]

enter image description here

Not only do these two plots not resemble each other, but I also know that the power spectrum of white noise should be flat (it isn't in plot 2).

I am new to signal processing. What am I missing to make this work?

  • $\begingroup$ As you are using an uncorrelated time history I would expect the autocorrelation to be approximately zero. Try using a sine function. Also I would use Abs Fourier with appropriate parameters. $\endgroup$ – Hugh Oct 5 '15 at 5:47
  • $\begingroup$ A further thought. The autocorrelation of a random time history is a delta function at lag zero. The Fourier transform of a delta function is a flat spectrum. One approach to develop your method would be to low pass filter the time history to give some correlation between points. $\endgroup$ – Hugh Oct 5 '15 at 7:48

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.