本帖最后由 1393107100 于 2019-5-1 11:10 編輯

clear,clc

close all

%%%%%%　小波變換　%%%%%%%%%%%%%%%%

fs=1024;

t=1/fs:1/fs:1;

f1=100;f2=200;f3=300;

s=cos(2*pi*f1*t.*(t>=0&t<0.3))+2*cos(2*pi*f2*t.*(t>=0.3&t<0.8))+3*sin(2*pi*f3*t.*(t>=0.8&t<=1));

%s=cos(2*pi*f1*t);

figure(1);

plot(t,s,'b')

title('原始單頻信號')

Y = awgn(s,6);%加高斯白噪聲

subplot(412);

plot(t,Y,'r')

title('加信噪比為6的白噪聲')

%%%%%%%%%%%%%%%%%小波時頻圖繪制 %%%%%%%%%%%%%%%%%%%%%

wavename='cmor4-4';? ? %%選用帶寬參數和中心頻率均為4 的復morlet小波

totalscal=256;? ?? ?? ? %尺度序列的長度，即scal的長度

fc=centfrq(wavename);%% 小波的中心頻率

%%%%%%%%% 繪制時間域小波波形 %%%%%%%%%%%%%%%%%%%

cparam=2*fc*totalscal;%為得到合適的尺度所求出的參數

a=totalscal:-1:1;

scal=cparam./a;%得到各個尺度，以使轉換得到頻率序列為等差序列

coefs=cwt(Y,scal,wavename,1/fs);%求連續小波系數,小波系數coefs(系數是復數時要取模)

f=scal2frq(scal,wavename,1/fs);% 將尺度轉換為頻率.F = scal2frq(A,'wname',DELTA),該函數能將尺度轉換為實際頻率，其中A為尺度，wname為小波名稱，DELTA為采樣周期。

figure(2);

imagesc(t,f,abs(coefs)); %繪制色譜圖

colormap(jet)

colorbar

xlabel('時間 t/s')

ylabel('頻率 f/Hz')

title('小波時頻圖')

abs_coefs=abs(coefs).^2;%abs函數：數值的絕對值和復數的幅值,

bw=abs_coefs>=0.7*mean(mean(abs_coefs));%mean數組的均值.mean(A)如果A為矩陣，那么返回包含每列均值的行向量。

bw=bwlabel(bw);%在BW數組中，0代表黑背景，1代表白

for i=1:max(max(bw))

if length(find(bw==i))>1000

bw(find(bw==i))=-1;

end

end

bw=bwlabel(bw==-1);

figure(3)

imagesc (bw); title('脊線區域')

out=zeros(size(coefs));% 輸出：大小為輸入的矩陣大小，各元素的值為0

for i=1:max(max(bw))

out=out+local_maxima(abs_coefs,bw,i);%小波系數模局部極大值

end

out=out.*bw;

RIDGE=zeros(size(coefs));

for i=1:max(max(out))

[x,y]=find(out==i);

x=x(1:100:end);

if length(x)<3

break;

end

cs = spline(y(1:100:end),x);

%spline三次方樣條數據插值

%pp = spline(x,y) 返回一個分段多項式結構體以用于 ppval 和樣條實用工具 unmkpp

x_new=round(ppval(cs,y(1):y(end)));

x_new(find(x_new<=0))=1;

for i=1:length(x_new)

RIDGE(x_new(i),y(i))=1;

end

end

figure,imagesc (RIDGE);??title('脊線')

時頻圖.jpg

(22.45 KB, 下載次數: 2)

2019-5-1 11:10 上傳

時頻圖

脊線.jpg

(14.9 KB, 下載次數: 0)

2019-5-1 11:08 上傳

脊線圖

總結

以上是生活随笔為你收集整理的matlab小波脊线,小波脊线提取，模极大值法。运行的结果不太对，代码有些地方我也没完全看懂...的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。