分享最近自己用matlab编写的连续小波变换程序

xinglong-liu

分享最近自己编写的一个连续小波变换程序，希望对大家有用。

close all
clear all;
clc;
dt=1e-2;
t=0:dt:4-dt;
s=sin(12*pi*t);
s(1:end/2) =s(1:end/2)+sin(8*pi*t(1:end/2));
s(end/2+1:end) =s(end/2+1:end)+sin(24*pi*t(end/2+1:end));

N0=length(t)                %should be a even number
figure(2)
subplot(211)
plot(t,s)

fc=2;
fb=2;
T=8;
M0=T/dt;
wc=fc*2*pi;

fmax=20                   %input the highest frequency
fmin=2                      %input the lowest frequency
nscale=64*4;             %input the scale number  
amax=fc/fmin;
amin=fc/fmax;
amax=log2(amax);
amin=log2(amin);
scale=linspace(amin,amax,nscale)';
a=2.^scale;
freq=fc./a;

Mmax=floor(M0*max(a))
Mmax=(2^ceil(log2(2*Mmax+N0))-N0)/2
s=[zeros(1,Mmax) s zeros(1,Mmax)];
N=N0+2*Mmax
S=fft(s);

for m=1:length(a)
    aa=a(m);
    M=floor(M0*aa);
    if ~mod(M,2)
        M=M+1;
    end
    if M<64
        disp('the scale is too small')
    end
    t1=[-(M-1)/2:1:(M-1)/2]*dt;
    t1=t1/aa;
    morlet=(pi*fb)^(-1/2)*exp(i*2*pi*fc*t1).*exp(-t1.^2/fb);
    morlet_R=morlet(1:(M-1)/2);
    morlet_L=morlet((M+1)/2:M);
    morlet=[morlet_L zeros(1,N-M) morlet_R];
    MORLET=fft(morlet);
    Coff=MORLET.*S;
    cwt_morlet(m,:)=ifft(Coff);
end
cwt_coff=cwt_morlet(:,Mmax+1:Mmax+N0);
figure(2)
subplot(212)
set(gcf,'Color','w');
pcolor(t,freq,abs(cwt_coff).^2);
colormap jet;shading interp;

wbh4158

太感谢了，最近正研究这个

wht0211

回复 xinglong-liu 的帖子
不错

   

stao136

非常不错，很有价值，

secondye

楼主，你的贡献很大啊，谢谢

xiaofeilonghe

刚踏进信号这方面的专业，哎，看着前辈们的成果，很是羡慕

shinhwaforever

受益了,谢谢楼主分享

c2019

受益了,谢谢楼主分享

c2019

有些关键步骤的物理意义还是不太懂，楼主可以在加些注释吗，

dong1104

非常不错

水野遥

已复制保留，呵呵，慢慢学习。