Estimate PSD from CWT
5 views (last 30 days)
Show older comments
Good afternoon! I want to estimate the power spectral density (psd) using the cwt coefficients. I have read that a good estimator is the global spectrum wavelet (gws). But I have the following problem: the sample rate is 1000Hz, so the original psd goes from 0 to 500hz. When I convert the scales of my interest (from 1 to 16), with the db6 wavelet, to frequencies, these are 45 to 700 Hz. How do I compare the traditional psd with gws? Any algorithm to help me? thank you
0 Comments
Answers (1)
Drishti
on 10 Oct 2024
Hi Gisela,
I understand that you are trying to estimate Power Spectral Density (PSD) and wants to compare the traditional PSD with Global Wavelet Spectrum (GWS) using ‘db6’ wavelet.
To estimate PSD, you can utilize the ‘pwelch’ function that follows Welch’s method.
Refer to the code snippet below for estimating PSD:
% Compute traditional PSD using Welch's method
[pxx, f] = pwelch(signal, [], [], [], fs);
Also, for the generation of GWS, you can utilize the ‘wavedec’ function with ‘db6’ wavelet.
Refer to the implementation below for better understanding:
% Perform DWT using db6
maxLevel = 16; % Maximum level of decomposition
[coeffs, levels] = wavedec(signal, maxLevel, 'db6');
% Calculate the corresponding frequencies for each scale
center_frequency = 0.666 * (fs / 2);
% Calculate frequencies for each scale
frequencies = center_frequency ./ (2.^(0:maxLevel-1));
% Estimate PSD from DWT coefficients (GWS)
gws = zeros(size(frequencies));
for i = 1:maxLevel
% Reconstruct signal from detail coefficients at each level
reconstructed = wrcoef('d', coeffs, levels, 'db6', i);
% Compute power of reconstructed signal
gws(i) = sum(reconstructed.^2);
end
After calculating the values of both PSD and GWS, you can compare them by simply plotting the values.
For more information, you can refer to the MATLAB Documentation of ‘pwelch’ and ‘wavedec’ functions.
I hope this helps!
0 Comments
See Also
Categories
Find more on Parametric Spectral Estimation in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!