【通信】Matlab实现多同步压缩变换

1 内容介绍

时频分析作为非平稳信号处理的一种重要方法,它提供了时域与频域的联合分布信息,清楚的刻画了信号频率随时间的变化关系,能给出每个时刻下信号的瞬时频率及其幅值。本文把时频分析作为研究对象,在对传统时频分析方法做出深入分析的基础上,重点研究同步压缩连续小波变换的时频分析方法,并将其应用到解决多分量非平稳信号的时频聚集性、分离以及瞬时频率估计的问题上。

2 部分代码

clear

load('vib_data1.mat')

fs = 12000; N = 1200;      % sampling frequency and points

time = (1:N)/fs;              % time sequence

fre = (fs/2)/(N/2):(fs/2)/(N/2):(fs/2);    % frequency sequence

data=data(1:N);

[MTs tfr]=MSST_Y_new(data,95,50);

gamma = 10^(-2);

sigma = 0.025;

[~,~,~,~,~,Ts1,Ts2,Ts3,Ts4] = sstn(data,gamma,sigma);

.................................................................

figure

suptitle('Fig. 19');

subplot(211)

plot(time,data);

ylabel('Amp / g');

xlabel('Time / s');

axis([0 0.1 -4 4]);

subplot(212)

ft=abs(fft(data))/600;

plot(fre,ft(1:end/2));

ylabel('Amp / g');

xlabel('Fre / Hz');

axis([0 6000 0 0.3]);

...........................................................................

x1=0.06; x2=0.065;

y1=2.6;   y2=3.7;

dd=0.1;

figure

suptitle('Fig. 20');

subplot(511);

imagesc(time,fre/1000,abs(tfr));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy;axis ([0 0.1 y1-dd y2+dd]);

rectangle('Position',[x1 y1 x2-x1 y2-y1],'EdgeColor','red','Linewidth',1);

axes('position',[0.92,0.805,0.07,0.15]); 

imagesc(time,fre/1000,abs(tfr));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy

axis off;

%set(ha,'xlim',[x1 x2],'ylim',[y1 y2]);

xlim([x1,x2]);ylim([y1,y2]);

subplot(512);

imagesc(time,fre/1000,abs(Ts1));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy;axis ([0 0.1 y1-dd y2+dd]);

rectangle('Position',[x1 y1 x2-x1 y2-y1],'EdgeColor','red','Linewidth',1);

axes('position',[0.92,0.61,0.07,0.15]); 

imagesc(time,fre/1000,abs(Ts1));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy

axis off;

%set(ha,'xlim',[x1 x2],'ylim',[y1 y2]);

xlim([x1,x2]);ylim([y1,y2]);

subplot(513);

imagesc(time,fre/1000,abs(Ts2));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy;axis ([0 0.1 y1-dd y2+dd]);

rectangle('Position',[x1 y1 x2-x1 y2-y1],'EdgeColor','red','Linewidth',1);

axes('position',[0.92,0.44,0.07,0.15]); 

imagesc(time,fre/1000,abs(Ts2));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy

axis off;

%set(ha,'xlim',[x1 x2],'ylim',[y1 y2]);

xlim([x1,x2]);ylim([y1,y2]);

subplot(514);

imagesc(time,fre/1000,abs(Ts4));

%xlabel('Time / s');

ylabel('Fre / kHz');

rectangle('Position',[x1 y1 x2-x1 y2-y1],'EdgeColor','red','Linewidth',1);

axis xy;axis ([0 0.1 y1-dd y2+dd]);

axes('position',[0.92,0.26,0.07,0.15]); 

imagesc(time,fre/1000,abs(Ts4));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis off;

%set(ha,'xlim',[x1 x2],'ylim',[y1 y2]);

xlim([x1,x2]);ylim([y1,y2]);

subplot(515);

imagesc(time,fre/1000,abs(MTs));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy;axis ([0 0.1 y1-dd y2+dd]);

rectangle('Position',[x1 y1 x2-x1 y2-y1],'EdgeColor','red','Linewidth',1);

axes('position',[0.92,0.09,0.07,0.15]); 

imagesc(time,fre/1000,abs(MTs));

%xlabel('Time / s');

ylabel('Fre / kHz');

axis xy

axis off;

%set(ha,'xlim',[x1 x2],'ylim',[y1 y2]);

xlim([x1,x2]);ylim([y1,y2]);

..........................................................................

MTs1=MTs;

MTs1(1:310,:)=0;

MTs1(360:end,:)=0;

MTs2=MTs;

MTs2(1:260,:)=0;

MTs2(310:end,:)=0;

[Cs2] = brevridge_mult(abs(MTs2), (1:2)/13.9, 2, 1, 5);

[Cs1] = brevridge_mult(abs(MTs1), (1:2)/13.9, 2, 1, 5);

IF1=[Cs1(1,1:1130) Cs1(2,1131:end)];

IF2=[Cs2(2,1:176) Cs2(1,177:end)];

Cs(1,:)=IF1;Cs(2,:)=IF2;

n=1200;

ds=3;

for k=1:2

for j=1:n

MTs_sig(k,j)=sum(real(MTs(max(1,Cs(k,j)-ds):min(round(n/2),Cs(k,j)+ds),j)));

end

end

figure

suptitle('Fig. 21');

subplot(411)

plot(time,fre(IF1)/1000,'b-');hold on;plot(time,fre(IF2)/1000,'b-');

axis ([0 0.1 y1-dd y2+dd]);

ylabel('Fre / kHz');

subplot(412)

plot(time,MTs_sig(1,:),'b-');

ylabel('Amp / g');

axis ([0 0.1 -2.5 2.5]);

subplot(413)

plot(time,MTs_sig(2,:),'b-');

ylabel('Amp / g');

axis ([0 0.1 -1.5 1.5]);

subplot(414)

plot(time,data,'k-','linewidth',1);hold on;

plot(time,sum(MTs_sig),'r-');

ylabel('Amp / g');

axis ([0 0.1 -2.5 2.5]);

legend('Measured signal','S1 + S2');

.......................................................

figure

suptitle('Fig. 22');

subplot(211)

plot(time,real(MTs_sig(1,:))/max(real(MTs_sig(1,:))),'b-','Linewidth',2);hold on;

plot(time,detrend(fre(Cs(1,:)))/max(detrend(fre(Cs(1,:))))/1.25,'r-','Linewidth',1.5);

axis ([0 0.1 -1 1]);

xlabel('Time / s');

subplot(212)

plot(time,real(MTs_sig(2,:))/max(real(MTs_sig(2,:))),'b-','Linewidth',2);hold on;

plot(time,detrend(fre(Cs(2,:)))/max(detrend(fre(Cs(2,:))))/1.25,'r-','Linewidth',1.5);

axis ([0 0.1 -1 1]);

legend('Mono-component mode','Detected IF');

xlabel('Time / s');​

3 运行结果

4 参考文献

[1]王智慧. 小波时频同步压缩变换方法研究及应用[D]. 西安电子科技大学, 2020.

博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机、雷达通信、无线传感器等多种领域的Matlab仿真，相关matlab代码问题可私信交流。