Reputation: 109
Below I present my code how I calculate power spectrum for my signal:
Fs=100;
n = length(x) (in my example always n=160);
whann = window(@hann,n).';
x = whann.*x';
xdft = fft(x);
xdft = (2/length(x))*xdft(1:length(x)/2+1);
xdft = abs(xdft);
xdft=xdft.^2;
freq = 0:Fs/length(x):Fs/2;
Now, I would like to calculate area under power spectrum but only for frequency range 4-6 Hz. The 32 first elements of vector freq look like this:
freq = [0,00 0,28 0,56 0,83 1,11 1,39 1,67 1,94 2,22 2,50 2,78 3,06 3,33 3,61 3,89 4,17 4,44 4,72 5,00 5,28 5,56 5,83 6,11 6,39 6,67 6,94 7,22 7,50 7,78 8,06 8,33]
So, I can find only area between 4,17 Hz- 6,11 Hz.
Could you advice me, how to calculate area under spectrum for specific frequency range (as I mentioned above for example 4-6 Hz)?
Thanks in advance for any help
Upvotes: 1
Views: 3088
Reputation: 2750
I would proceed as follows:
idx = find(freq>=4 & freq<=6);
trapz(freq(idx),spectrum(idx))
If I understood your question right, what stated above should lead you to the results you wanna estimate.
EDIT
Since you don't have spectrum values for freq=4Hz
and freq=6Hz
, I would suggest to interpolate values like this:
int_spec = exp(interp1(log(freq),log(spec),log(4:.1:6),'linear','extrap'))
and then call
trapz(4:.1:6,int_spec)
Upvotes: 2