Welch’s Algorithm of Averaging Modified Periodograms

Question

I am second year Undergrad student in computer engineering . I have a homework problem where I needed to compute Power Spectral Density . I'm using MATLAB as I noticed that it has good documentation and also a function to that. I'm going through this page to understand using dsp.SpectrumAnalyzer System object™

https://www.mathworks.com/help/dsp/ug/spectral-analysis.html

I used this object to get my results , but i failed to understand how it is doing so. I vaguely understand the concept of PSD that it is computed by Fourier transform of AutoCorr of discrete sequence as FT cannot be applied directly to the random processes.

I'm stuck and confused on this point itself

If I have "m" samples and I need to compute PSD , Is this the correct procedure as described in POINT 1

• Divide m Samples into N point segments . --> this results in total of m/N segments (true ?)

• Each segment which has N points is further divided into L overlapping segment having M points .

• And these M points might have D overlapped points

  Is my understanding correct ?

Thanks for taking time to read my question.

Answer 1

You can relate the time-averaged autocorrelation sequence $r_{xx}[k]$ to the sample sequence $x[n]$ as follows. You have that the sequence $r_{xx}[k]$ does not depend on $n$ but on the lag $k$ $$ r_{xx}[k] = \frac 1N \sum_{n = 0}^{N - k - 1}x^{*}[n]x[n + k], \qquad 0\leq k\leq N - 1\tag{1} $$ The corresponding PSD, the auto-PSD, estimate is $$ P_{xx}(f) = \sum_{k = -(N-1)}^{N - 1}r_{xx}[k]\exp\left(-j2\pi f k\right)\tag{2} $$ From $(1)$ and $(2)$ you can simply compute the peridogram $P_{xx}(f)$ as function of the sample sequence $x[n]$ instead of $r_{xx}[m]$ $$ P_{xx}(f) = \frac 1N \bigg|\sum_{n = 0}^{N - 1}x[n]\exp\left(-j2\pi f n\right)\bigg|^2 = \frac 1N \big|X(f)\big|^2\tag{3} $$ The PSD in $(3)$ describes the distribution in frequency of the power of $x[n]$. Also have a look at this question and answer.

Welch's method does two additional things as follows

• Windowing of each segment. Instead of $x[n]$, you use $x_i[n]w[n]$ where $w[n]$ is the chosen window function where $x_i[n]$ are signal samples over the $i^{\rm th}$ segment. Here $i = 0, 1, \ldots, L - 1$ and $n = 0, 1, \ldots, M - 1$

• Normalization of the periodogram (over windowed data segments) by the power of the chosen window function. You normalize everything by $P_w$ where $$P_w = \frac 1M\sum_{n = 0}^{M-1} w^2[n]$$

• from the two steps above you get what's called the modified periodogram over each segment as follows: $$P_{xx}^i(f) = \frac{1}{M\cdot P_w}\bigg\lvert \sum_{n = 0}^{M - 1}\big(x_i[n]w[n]\big)\exp\left(-j2\pi f n\right)\bigg\rvert^2\tag{4}$$ With $i = 0, 1, \ldots, L - 1$ and $n = 0, 1, \ldots, M - 1$. Meaning we have $L$ segments in total, and each one of them is of length $M$.

The resulting Welch PSD estimate is then the average of the modified periodograms in $(4)$ as follows

$$P_{xx}^{\rm Welch}(f) = \frac{1}{L}\sum_{i = 0}^{L - 1}P_{xx}^i(f)$$

Note that your second bullet point is not correct. You overlap the segments in your first bullet point, and not subdivide them further for overlap.

Answer 2

Understanding the Windowing Method in PSD Calculation

Please go through the link .. Someone has posted their code snippet which matches with MATLAB output. Hope that might give you better understanding