I'm trying to calculate the power spectral density (PSD) to analyze signal power. To do so, I wrote a MATLAB code with the great help of others and below is my code.
function [psd, freq] = mypsd(data, fs) N = length(data); f_delta = fs/floor(N/2)/2; % Modify N to even number freq = (0:f_delta:fs/2)'; xdft = fft(data); xdft = abs(xdft/N); xft = xdft(1:floor(N/2)+1); xft(2:end-1) = 2*xft(2:end-1); psdx = xft.^2 / (2*f_delta); psd = 10*log10(psdx); end
To confirm the code, I compared the result with the one from the
my_result = mypsd(data, fs); mat_result = pwelch(data, rectwin(length(data)), 0, length(data), fs);
In this case above, two results were same.
After that, I tried to segment the data with overlapping 50% of each segment like below, where nfft is 256, n_overlap is nfft/2 = 128.
my_psd_total = zeros(nfft/2+1, 1); for iter = 1:num_segment data_win = data((iter-1)*(nfft - n_overlap)+1:(iter-1)*(nfft - n_overlap)+nfft); [psd_win, freq_win] = mypsd(data_win.*hann(nfft), fs); my_psd_total = my_psd_total + psd_win; end my_psd_total = my_psd_total / num_segment
However, the result is different with the result of
pwelch(data, hann(nfft), n_overlap, nfft, fs)
As I know, elementwise multiplying of window with segment in time domain and average the result is the proper method of using window functions like hanning or hamming et al.
Is there anyone can explain why my result is different with built-in pwelch function?