I'm trying to understand PSD, and, it's functions.
So far I have carried out the following:
- Segmented a signal with
NFFT = 256
andhop_size = 128
- Windowed each of these segments (Hamming)
- DFT of each Window
- Taken the first half of the resulting matrix
(n/2+1)
From what I have read so far, the algorithm can be best described as the following:
We take the squared magnitude of the LENGTH NFFT
, we calculate this by: mag = (bin.re * bin.re + bin.im * bin.im)
. This gives the magnitudes of each of the points in the DFT bin.
Q1: These are "averaged" to form Pxx
Pxx
is a bin of size NFFT/2+1
which I assume store the averages, but, do these store the averages for each of the bins? Because I don't quite understand:
If I have a resulting 2D matrix containing values:
a = [[1, 2, 3, 4, 5, 6, 9, 8],
[2, 5, 8, 9, 7, 8, 9, 7],
[4, 5, 6, 4, 7, 8, 9, 8]
I then take the "average" of these.. I.e. (1+2....., n)/NFFT/2+1. This is only going to give one value, for each of the bins?
OK, so now onto Scaling. Can we therefore scale to the Frequency? If so, do we scale by Fs * window_size
EDIT:
I don't think my last question made that much sense..
So far I have a 2D vector containing 146x128
..
block1 = [1, 2, 3, ...... n] (where `n` is 128)
block2 = [0, 5, 6, ...... n]
.....
block146 = [1, 4, 5, .... n]
Now, if I calculate the average magnitude of each of these blocks, i.e. sum(sqrt(re*re+im*im))/128
this will give me the averages for each block (a vector containing 146 samples).
For each of my samples in the blocks
do I therefore divide each of these by the average magnitude for that particular block? I.e.
block1 = [0, 1, 2, ...... n]
average = (0+1+2)/128 = 0.58 (for example)
result = [0/0.58, 1/0.58, 2/0.58.... n]
Therefor my resulting vector will still be a 2D vector.
Does this make sense?
EDIT:
1 - What type of signal are you trying to process? I.e: What type of application is it? Is it a unidimensional signal, like Audio sampled at 44.1Khz? Or is it a multidimensional signal like an image or video?
I'm trying to process a Audio signal.
What I have done so far is:
Taken a 1D input vector, this vector is then split into blocks (of size 256
with an overlap/hop of 128
this gives me a total of 146 blocks with each containing 128 samples. These samples have been multiplied by the Hamming Window and then a DFT has been passed through.
These blocks contain the resulting values after Windowing and DFT has taken place.
2 - Where are you performing these calculations? Matlab, Octave, python, or some other?
I am performing these calculations using C++
3 - When you say 2D vector (maybe could be less confusing to call this a Matrix), What does the numbers inside mean? what are the columns and what are the rows?
The columns and rows contain the DFT results of the Hanning windowing of each block at 256/2 = 128
Essentially, what I'm attempting to do is create a spectrogram. Currently, I can plot the spectrogram with magnitudes of the values, i.e. (re*re+im*im)
however, instead of plotting the magnitudes, I want to plot the overall power spectrum (PSD) so I can determine where the power is mainly within the signal.
If I do this using matplotlib
I can get the following results:
So in theory, instead of showing the magnitudes
which is what I'm currently doing, I want to show the PSD (like in the graph). I looked at the matplotlib source code and this is what was given:
result, windowVals = apply_window(result, window, axis=0,
return_window=True)
result = np.fft.fft(result, n=pad_to, axis=0)[:numFreqs, :]
result = np.conjugate(result) * result
# for PSD
result /= (np.abs(windowVals)**2).sum()
result[1:-1] *= scaling_factor
# where scaling_factor is either 1 or 2
To me, the above code looks like they are taking the squared sum of the WindowValues (Which is Hamming) and dividing each value from the result
(which is the resulting DFT) and finally multiplying by the scaling factor.
Does this look correct?