Plot Wavelet FFT in Hz

4 次查看（过去 30 天）
Alexander Voznesensky
评论： Robert U 2017-8-10
Hi! I'm trying to do fft of wavelet function. But I'm confused: how to define x-axis f in Hz? Here is my code:
scales=1:100;
wname='db20';
N=10;
[phi, psi, tgrid] = wavefun (wname, N);
f=???
figure;
hold on
plot(f,abs(fft(phi)));
plot(f,abs(fft(psi)));
before, in my practice, i used:
f=0:Fs/(N-1):Fs;
What should I do now?

请先登录，再回答此问题。

回答（2 个）

Robert U
Robert U 2017-8-8
编辑：Robert U 2017-8-8
Hello Alexander,
if you want to display wavelet analysis in terms of frequency one common way is to use the center frequency of the used wavelet.
Using a web search engine you could have searched for "wavelet equivalent frequency" which is providing several useful explanations as: https://de.mathworks.com/videos/understanding-wavelets-part-1-what-are-wavelets-121279.html
Kind regards,
Robert
  1 个评论
Alexander Voznesensky
Alexander Voznesensky 2017-8-8
编辑：Alexander Voznesensky 2017-8-8
Thank you. I tried:
[freq,xval,recfreq] = centfrq(wname, N, 'plot');
And I know the points number:
N=length(xval);
But what about Fs? How to define f vector?

请先登录，再进行评论。

Robert U
Robert U 2017-8-9
Hi Alexander,
if you want to plot the Fourier spectrum of the wavelet you can apply the usual fft algorithm (described here: https://de.mathworks.com/help/matlab/math/basic-spectral-analysis.html ).
%%Create wavelet
wname='db20';
N=12;
[phi, psi, tgrid] = wavefun (wname, N);
%%Calculate FFT of wavelet output
SampleFrq = 2^N; % calculate sample frequency of signals
SigL = length(psi); % phi & psi have same length
NFFT = 2^nextpow2(SigL); % standard fft calculation
SigMagPsi = fft(psi,NFFT)/SigL;
SigMagPsi = 2*abs(SigMagPsi(1:NFFT/2+1));
SigMagPhi = fft(phi,NFFT)/SigL;
SigMagPhi = 2*abs(SigMagPhi(1:NFFT/2+1));
f = SampleFrq / 2 * linspace(0,1,NFFT/2+1);
%%Plot Figure
figure;
hold on
plot(f,SigMagPhi);
plot(f,SigMagPsi);
ah = gca;
ah.XScale = 'log';
ah.YLabel.String = 'Magnitude';
ah.XLabel.String = 'Frequency [Hz]';
legend(['Phi';'Psi'])
figure;
centfrq(wname,N,'plot')
  2 个评论
Alexander Voznesensky
Alexander Voznesensky 2017-8-9
编辑：Alexander Voznesensky 2017-8-9
Thank you very much! But, the main question is still here: why
SampleFrq = 2^N?
I have some stuff about it:
Fs=10^3; % частота дискретизации, Гц
T=1/Fs; % период дискретизации, сек
N=1e4; % количество отсчетов
t=(0:N-1)*T; % массив отсчетов времени
f=0:Fs/(N-1):Fs; % массив частот
f1=10;
f2=25;
f3=50;
f4=100;
x1=cos(2*pi*f1.*t)+cos(2*pi*f2.*t)+cos(2*pi*f3.*t)+cos(2*pi*f4.*t); % стационарный сигнал
figure;
subplot(2,1,1);
plot(t,x1);
xlabel('t, сек');
ylabel('x');
title('Стационарный сигнал. Временная область');
grid;
subplot(2,1,2);
plot(f,abs(fft(x1)));
xlabel('f, Гц');
ylabel('X');
title('Стационарный сигнал. Частотная область');
grid;
x2_1=cos(2*pi*f4.*t(1:N/4));
x2_2=cos(2*pi*f3.*t(N/4+1:N/2));
x2_3=cos(2*pi*f2.*t(N/2+1:3*N/4));
x2_4=cos(2*pi*f1.*t(3*N/4+1:end));
x2=[x2_1 x2_2 x2_3 x2_4]; % нестационарный сигнал
figure;
subplot(2,1,1);
plot(t,x2);
xlabel('t, сек');
ylabel('x');
title('Нестационарный сигнал. Временная область');
grid;
subplot(2,1,2);
plot(f,abs(fft(x2)));
xlabel('f, Гц');
ylabel('X');
title('Нестационарный сигнал. Частотная область');
grid;
So my Fs is 1000 Hz. What in this case? Is Fs the same for the Fourier spectrum of the wavelet? So:
f=0:Fs/(length(psi)-1):Fs;
Robert U
Robert U 2017-8-10
Hello Alexander:
The basic idea is to understand the generated wavelet as a time-dependent signal that is sampled with a certain (fixed) sampling frequency even though it could be a function depending on any physical (scalar) value (e.g. distance).
Sampling frequency is the step size of your points grid XVAL (tgrid). Assuming that they are equidistantly distributed you can check the equivalent sampling frequency by
Fs = 1 / mean( diff( tgrid ));
You will see that the sampling frequency increases with the number of iterations. You can calculate your sampling frequency either with the formula above or you convince yourself that the scaling is according to
Fs = 2^N;
This feature was documented on earlier versions ( http://radio.feld.cvut.cz/matlab/toolbox/wavelet/wavefun.html ) and is still for wavefun2 ( https://de.mathworks.com/help/wavelet/ref/wavefun2.html).
[PHI,PSI,XVAL] = wavefun('wname',ITER) returns the scaling and wavelet functions on the 2^ITER points grid XVAL.
Kind regards,
Robert

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Continuous Wavelet Transforms 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by