matlab 指数信号的重构,小波变换重构信号求信噪比

clc

clear all

%% 注释

% 分离多频率信号(将低频信号分离出来)并进行重构

% 多分辨率分析：用于分解信号低频部分。在处于低频时，频率分辨率高，而处于高频时，有较低的频率分辨率

%主要调用函数：[C,L]=wavedec(x,n,wname),x为需要处理的信号(加噪信号),n为小波分解层数,wname为小波名

%                  C为小波分解向量(实值矢量)，L为正整数的bookkeeping vector

%              A = appcoef(___,N)返回第N层近似系数，如果[C,L]是一维信号的M层小波分解，则0<=N<=M

%              D = detcoef(C,L,N)从小波分解结构提取第N层的细节系数

%              X = waverec(C,L,'wname')基于多小波分解结构[C,L]和小波名重建信号

%% 产生原始信号

N=10000;

fs=10000000;

n=0:N-1;

t=0:(1/fs):(N-1)*(1/fs);

f=100000;

s=zeros(1,N);%生成1乘N全零矩阵，用于存储数据

s=s+1.001*sin(2*pi*f*t)+0.0004*sin(6*pi*f*t)+0.0002*sin(10*pi*f*t);

noise=0.0005*randn(1,N);%生成均值为0，方差为0.0005的一个1*N的随机数

s=s+noise;

%% 画图，原始信号

figure;

subplot(1,2,1),plot(s(1:1000));

xlabel('采样点数');%对坐标轴x进行标注

ylabel('样点数值');title('10M采集(原始信号)');grid on;%对坐标轴y进行标注，设置图形标题，打开网格

%% 加窗画功率谱图

w = blackman(length(s(1:1000)));%布莱克曼窗，常用来检测两个频率相近幅度不同的信号，根据长度length(s(1:1000))产生一个布莱克曼窗w

subplot(1,2,2),periodogram(s(1:1000),w,numel(s(1:1000)),fs,'power');%返回功率谱，numel(A)返回数组A中元素个数

xlabel('频率(MHz)');

ylabel('振幅(dB)');title('功率谱');grid on;

%% 求信噪比和信纳比

[Sxx, F] = periodogram(s(1:1000),w,numel(s(1:1000)),fs,'power');

rbw = enbw(w,fs);

% THD = thd(Sxx,F,rbw,'power')

% SFDR = sfdr(Sxx,F,rbw,'power')

SNR = snr(Sxx,F,rbw,'power');

SINAD = sinad(Sxx,F,rbw,'power');

%% 第一次小波分解，分离高频与低频

[c,l] = wavedec(s,4,'sym4');

a4=appcoef(c,l,'sym4',4);%第四层分解低频系数

% a3=appcoef(c,l,'sym8',3);

% a2=appcoef(c,l,'sym8',2);

% a1=appcoef(c,l,'sym8',1);

d4=detcoef(c,l,4);

d3=detcoef(c,l,3);

d2=detcoef(c,l,2);%第二层分解高频系数

d1=detcoef(c,l,1);

%% 第一次对分离出的低频信号进行重构

dd4=zeros(1,length(d4));

dd3=zeros(1,length(d3));

dd2=zeros(1,length(d2));

dd1=zeros(1,length(d1));

c1=[a4 dd4 dd3 dd2 dd1];

s1=waverec(c1,l,'sym4');

% figure;

% plot(t(1:1000),s1(1:1000));

% title('重构信号');

%% 绘制原始信号和分离后信号的频谱图

% X = fft(s1);                  %计算信号傅里叶变换

% P2 = abs(X/N);               %计算双侧频谱P2

% P1 = P2(1:N/2+1);            %基于P2和偶数信号长度N计算单侧频谱P1

% P1(2:end-1) = 2*P1(2:end-1);

% f = fs*(0:(N/2))/N;          %定义频域f并绘制单侧幅值频谱

% figure;

% plot(f(1:1000),P1(1:1000));

%% 画图，第一次小波重构

figure;

subplot(1,2,1),plot(s1(1:1000));

xlabel('采样点数');%对坐标轴x进行标注

ylabel('样点数值');title('第一次小波重构信号');grid on;%对坐标轴y进行标注，设置图形标题，打开网格

%% 加窗画功率谱图

w = blackman(length(s1(1:1000)));%布莱克曼窗，常用来检测两个频率相近幅度不同的信号，根据长度length(s(1:1000))产生一个布莱克曼窗w

subplot(1,2,2),periodogram(s1(1:1000),w,numel(s1(1:1000)),fs,'power');%返回功率谱，numel(A)返回数组A中元素个数

xlabel('频率(MHz)');

ylabel('振幅(dB)');title('第一次小波重构信号功率谱');grid on;

%% 求信噪比和信纳比

[Sxx, F] = periodogram(s1(1:1000),w,numel(s1(1:1000)),fs,'power');

rbw = enbw(w,fs);

% THD = thd(Sxx,F,rbw,'power')

% SFDR = sfdr(Sxx,F,rbw,'power')

SNR1 = snr(Sxx,F,rbw,'power');

SINAD1 = sinad(Sxx,F,rbw,'power');