信号小波分析平面图谱
对信号做cwt想得到时间频率图类似伪彩色图一样的平面图谱，不知道怎么操作?问题能说的再详细一点吗 谢谢 请参考以下说明与程序：
此程序实现的是将50HZ和100HZ正弦波分离的程序。并对小波谱分析的实现给出了启发性程序。由于最近太忙，无法解释其中的理论，紧把程序和图形贴出来，供大家讨论和批评指正。
程序如下：
%% 本程序实现小波用作谱分析
clc;clear;
%% 1.正弦波定义
f1=50; % 频率1
f2=100; % 频率2
fs=2*(f1+f2); % 采样频率
Ts=1/fs; % 采样间隔
N=120; % 采样点数
n=1:N;
y=sin(2*pi*f1*n*Ts)+sin(2*pi*f2*n*Ts); % 正弦波混合
figure(1)
plot(y);
title('两个正弦信号')
figure(2)
stem(abs(fft(y)));
title('两信号频谱')
%% 2.小波滤波器谱分析
h=wfilters('db30','l'); % 低通
g=wfilters('db30','h'); % 高通
h=; % 补零（圆周卷积，且增大分辨率变于观察）
g=; % 补零（圆周卷积，且增大分辨率变于观察）
figure(3);
stem(abs(fft(h)));
title('低通滤波器图')
figure(4);
stem(abs(fft(g)));
title('高通滤波器图')
%% 3.MALLET分解算法(圆周卷积的快速傅里叶变换实现)
sig1=ifft(fft(y).*fft(h)); % 低通(低频分量)
sig2=ifft(fft(y).*fft(g)); % 高通(高频分量)
figure(5); % 信号图
subplot(2,1,1)
plot(real(sig1));
title('分解信号1')
subplot(2,1,2)
plot(real(sig2));
title('分解信号2')
figure(6); % 频谱图
subplot(2,1,1)
stem(abs(fft(sig1)));
title('分解信号1频谱')
subplot(2,1,2)
stem(abs(fft(sig2)));
title('分解信号2频谱')
%% 4.MALLET重构算法
sig1=dyaddown(sig1); % 2抽取
sig2=dyaddown(sig2); % 2抽取
sig1=dyadup(sig1); % 2插值
sig2=dyadup(sig2); % 2插值
sig1=sig1(1,); % 去掉最后一个零
sig2=sig2(1,); % 去掉最后一个零
hr=h(end:-1:1); % 重构低通
gr=g(end:-1:1); % 重构高通
hr=circshift(hr',1)'; % 位置调整圆周右移一位
gr=circshift(gr',1)'; % 位置调整圆周右移一位
sig1=ifft(fft(hr).*fft(sig1)); % 低频
sig2=ifft(fft(gr).*fft(sig2)); % 高频
sig=sig1+sig2; % 源信号
%% 5.比较
figure(7);
subplot(2,1,1)
plot(real(sig1));
title('重构低频信号');
subplot(2,1,2)
plot(real(sig2));
title('重构高频信号');
figure(8);
subplot(2,1,1)
stem(abs(fft(sig1)));
title('重构低频信号频谱');
subplot(2,1,2)
stem(abs(fft(sig2)));
title('重构高频信号频谱');
figure(9)
plot(real(sig),'r','linewidth',2);
hold on;
plot(y);
legend('重构信号','原始信号')
title('重构信号与原始信号比较')
顺便说一下，这是分解两个信号的。要分解多个信号，可以采用多层，但注意2抽取的危险。还有一个方法，先为低频2个，高频1个。用一层正反变换。再把低频两个分开，再用一层正反变换。每次都用一层分解和重构。这有些像svm分两类问题，其实也可分三类。
首先，最重要的一点确定采样间隔。50hz，100hz分开它？首先明确子代滤波概念。50hz要用低通，100hz要用高通。因此最佳分割点75hz。因为小波二分频带，所以最高频率150hz。二分后，。但采样频率fs=2*f=300hz。
其次，什么样的小波。首先最主要的支撑长度，为什么？那就如同说，FIR滤波器为什么高阶效果好一样。所以我选了DB30，支撑长度60。再次，最好线性相位。最后，怎么看滤波器H(N),G(N)频谱。把原滤波器补零。补多少零，分析的信号长度减去滤波其长度。
最后，进行完MALLET分解后，也就是经过低通和高通滤波后。我们的要求达到了。可以通过小波系数的福利叶变换，看频谱。
这和福利叶的滤波器一样的。但小波变换始终是二分的，而福利叶变换是点频的。
关于多层特别小心，因为2抽取会使频谱混叠。所以最好采用些不抽取的算法。
还有小波滤波器是有过渡带的，在75HZ附近的信号频谱是危险的，它会混到两频带中。所以，任何小波分析前，要对信号频谱份量进行分析，任何盲目的分解都是不可行的。 cwt函数功能：实现一维连续小波变换的函数。
cwt函数语法格式： COEFS=cwt(S, SCALES, 'wname')
COEFS=cwt(S, SCALES, 'wname', 'plot')
COEFS=cwt(S, SCALES, 'wname', 'PLOTMODE')
COEFS=cwt(S, SCALES, 'wname', 'PLOTMODE', XLIM)
使用说明：cwt为一维小波变换的函数。
格式 COEFS=cwt(S, SCALES, 'wname') 采用'wname'小波，在正、实尺度SCALES下计算向量一维小波系数。
格式 COEFS=cwt(S, SCALES, 'wname', 'plot') 除了计算小波系数外，还加以图形显示。
格式 COEFS=cwt(S, SCALES, 'wname', 'PLOTMODE') 计算并画出连续小波变换的系数，并使用PLOTMODE对图形着色。
格式 COEFS=cwt(S, SCALES, 'wname', 'plot') 相当于 格式 COEFS=cwt(S, SCALES, 'wname', 'PLOTMODE') 中的语法 COEFS=cwt(S, SCALES, 'wname', 'absglb')
格式 COEFS=cwt(S, SCALES, 'wname', 'PLOTMODE', XLIM) 能够计算并画出连续小波变换的系数。系数使用PLOTMODE和XLIM进行着色。其中：XLIM=，并且有如下关系：1<=x1<=x2<=length(S)。
MODE值 含义
'lvl' scale-by-scale着色模式
'glb' 考虑所有尺度的着色模式
'abslvl'或'lvlabs' 使用系数绝对值的scale-by-scale着色模式
'absglb'或'glbabs' 使用系数绝对值并考虑所有尺度的着色模式
例子：cwt函数的应用举例
t=linspace(-1:1:512);
s=1-abs(t);
c=cwt(s,1:32,'cgau4');
c=cwt(s,);
c=cwt(s,,'db2');
c=cwt(s,1:64,'sym4','abslvl',);
猫头鹰先生 发表于 2016-6-22 12:52
请参考以下说明与程序：
此程序实现的是将50HZ和100HZ正弦波分离的程序。并对小波谱分析的实现给出了启发 ...
学习了真不错
页:
[1]