源程序代码
- 格式:docx
- 大小:16.66 KB
- 文档页数:7
clear allts=0.001; %信号抽样时间间隔t=0:ts:10-ts; %时间向量fs=1/ts; %抽样频率df=fs/length(t); %fft的频率分辨率msg=randint(100,1,[-3,3],123); %生成消息序列,随机数种子为123 msg1=msg*ones(1,fs/10); %扩展成取样信号形式msg2=reshape(msg1.',1,length(t));Pm=fft(msg2)/fs; %求消息信号的频谱f=-fs/2:df:fs/2-df;subplot(2,1,1)plot(t,fftshift(abs(Pm)))title('消息信号频谱')int_msg(1)=0; %消息信号积分for ii=1:length(t)-1int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;endkf=50;fc=250; %载波频率Sfm=cos(2*pi*fc*t+2*pi*kf*int_msg); %调频信号Pfm=fft(Sfm)/fs; % FM信号频谱subplot(2,1,2)plot(f,fftshift(abs(Pfm))) % 画出已调信号频谱title('FM信号频谱')Pc=sum(abs(Sfm).^2)/length(Sfm) %已调信号功率Ps=sum(abs(msg2).^2)/length(msg2) %消息信号功率fm=50;betaf=kf*max(msg)/fm % 调制指数W=2*(betaf+1)*fm % 调制信号带宽clear allts=0.001; %信号抽样时间间隔t=0:ts:10-ts; %时间向量fs=1/ts; %抽样频率df=fs/length(t); %fft的频率分辨率msg=sin(2*pi*[0:0.01:0.99]); %生成消息序列,随机数种子为123msg1=msg'*ones(1,fs/10); %扩展成取样信号形式msg2=reshape(msg1.',1,length(t));Pm=fft(msg2)/fs; %求消息信号的频谱f=-fs/2:df:fs/2-df;subplot(2,1,1)plot(t,fftshift(abs(Pm)))title('消息信号频谱')int_msg(1)=0; %消息信号积分for ii=1:length(t)-1int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;endkf=50;fc=250; %载波频率Sfm=cos(2*pi*fc*t+2*pi*kf*int_msg); %调频信号Pfm=fft(Sfm)/fs; % FM信号频谱subplot(2,1,2)plot(f,fftshift(abs(Pfm))) % 画出已调信号频谱title('FM信号频谱')Pc=sum(abs(Sfm).^2)/length(Sfm) %已调信号功率Ps=sum(abs(msg2).^2)/length(msg2) %消息信号功率fm=50;betaf=kf*max(msg)/fm % 调制指数W=2*(betaf+1)*fm % 调制信号带宽clear allts=0.001; %信号抽样时间间隔t=0:ts:5-ts; %时间向量fs=1/ts; %抽样频率df=fs/length(t); %fft的频率分辨率msg=randint(10,1,[-3,3],123); %生成消息序列,随机数种子为123 msg1=msg*ones(1,fs/2); %扩展成取样信号形式msg2=reshape(msg1.',1,length(t));Pm=fft(msg2)/fs; %求消息信号的频谱f=-fs/2:df:fs/2-df;subplot(3,1,1)plot(t,msg2) %画出消息信号title('消息信号')int_msg(1)=0; %消息信号积分for ii=1:length(t)-1int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;endkf=50;fc=300; %载波频率Sfm=cos(2*pi*fc*t+2*pi*kf*int_msg); %调频信号phase=angle(hilbert(Sfm).*exp(-j*2*pi*fc*t)); %FM调制信号相位phi=unwrap(phase);dem=(1/(2*pi*kf)*diff(phi)/ts); %求相位微分,得到消息信号dem(length(t))=0;subplot(3,1,2)plot(t,dem);title('无噪声的解调信号')y1=awgn(Sfm,20,'measured'); %调制信号通过AWGN信道y1(find(y1>1))=1; %调制信号限幅y1(find(y1<-1))=-1;phase1=angle(hilbert(y1).*exp(-j*2*pi*fc*t)); %信号解调phi1=unwrap(phase1);dem1=(1/(2*pi*kf)*diff(phi1)/ts);dem1(length(t))=0;subplot(3,1,3)plot(t,dem1);title('信噪比为20dB时的解调信号')clear allts=0.001; %信号抽样时间间隔t=0:ts:5-ts; %时间向量fs=1/ts; %抽样频率df=fs/length(t); %fft的频率分辨率msg=sin(2*pi*[0:0.01:0.99]); %生成消息序列,随机数种子为123 msg1=msg'*ones(1,fs/20); %扩展成取样信号形式msg2=reshape(msg1.',1,length(t));Pm=fft(msg2)/fs; %求消息信号的频谱f=-fs/2:df:fs/2-df;subplot(3,1,1)plot(t,msg2) %画出消息信号title('消息信号')int_msg(1)=0; %消息信号积分for ii=1:length(t)-1int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;endkf=50;fc=300; %载波频率Sfm=cos(2*pi*fc*t+2*pi*kf*int_msg); %调频信号phase=angle(hilbert(Sfm).*exp(-j*2*pi*fc*t)); %FM调制信号相位phi=unwrap(phase);dem=(1/(2*pi*kf)*diff(phi)/ts); %求相位微分,得到消息信号dem(length(t))=0;subplot(3,1,2)plot(t,dem);title('无噪声的解调信号')y1=awgn(Sfm,20,'measured'); %调制信号通过AWGN信道y1(find(y1>1))=1; %调制信号限幅y1(find(y1<-1))=-1;phase1=angle(hilbert(y1).*exp(-j*2*pi*fc*t)); %信号解调phi1=unwrap(phase1);dem1=(1/(2*pi*kf)*diff(phi1)/ts);dem1(length(t))=0;subplot(3,1,3)plot(t,dem1);title('信噪比为20dB时的解调信号')function [v,phi]=env_phas(x,ts,f0)if nargout==2 %nargout为输出变数的个数z=loweq(x,ts,f0); %产生调制信号的正交分量phi=angle(z); %angle是对一个复数求相角的函数endv=abs(hilbert(x)); %abs用来求复数hilbert(x)的模function [M,m,df]=fftseq(m,ts,df)fs=1/ts;if nargin==2 n1=0; %nargin为输入参量的个数else n1=fs/df;endn2=length(m);n=2^(max(nextpow2(n1),nextpow2(n2))); %nextpow2(n)取n最接近的较大2次幂M=fft(m,n); %M为信号m的傅里叶变换,n为快速傅里叶变换的点数,及基n-FFT变换m=[m,zeros(1,n-n2)]; %构建新的m信号df=fs/n;%重新定义频率分辨率function x1=loweq(x,ts,f0)t=[0:ts:ts*(length(x)-1)];z=hilbert(x); %希尔伯特变换对的利用---通过实部来求虚部x1=z.*exp(-j*2*pi*f0*t); %产生信号z的正交分量,t0=0.2; %信号的持续时间,用来定义时间向量ts=0.001; %抽样间隔fs=1/ts; %抽样频率fc=300; %载波频率,fc可以任意改变t=[-t0/2:ts:t0/2]; %时间向量kf=100; %偏差常数df=0.25; %所需的频率分辨率,用在求傅里叶变换时,它表示FFT的最小频率间隔m=sin(100*t); %调制信号,m(t)可以任意更改int_m(1)=0; %求信号m(t)的积分for i=1:length(t)-1int_m(i+1)=int_m(i)+m(i)*ts;end[M,m,df1]=fftseq(m,ts,df); %对调制信号m(t)求傅里叶变换M=M/fs; %缩放,便于在频谱图上整体观察f=[0:df1:df1*(length(m)-1)]-fs/2; %时间向量对应的频率向量u=cos(2*pi*fc*t+2*pi*kf*int_m); %调制后的信号[U,u,df1]=fftseq(u,ts,df); %对调制后的信号u求傅里叶变换U=U/fs; %缩放%通过调用子程序env_phas和loweq来实现解调功能[v,phase]=env_phas(u,ts,fc); %解调,求出u的相位phi=unwrap(phase); %校正相位角,使相位在整体上连续,便于后面对该相位角求导dem=(1/(2*pi*kf))*(diff(phi)*fs); %对校正后的相位求导%再经一些线性变换来恢复原调制信号%乘以fs是为了恢复原信号,因为前面使用了缩放subplot(3,2,1) %子图形式显示结果plot(t,m(1:length(t))) %现在的m信号是重新构建的信号,%因为在对m求傅里叶变换时m=[m,zeros(1,n-n2)] axis([-0.1 0.1 -1 1]) %定义两轴的刻度xlabel('时间t')title('原调制信号的时域图')subplot(3,2,2)plot(t,u(1:length(t)))axis([-0.1 0.1 -1 1])xlabel('时间t')title('已调信号的时域图')subplot(3,2,3)plot(f,abs(fftshift(M))) %fftshift:将FFT中的DC分量移到频谱中心axis([-600 600 0 0.04])xlabel('频率f')title('原调制信号的频谱图')subplot(3,2,4)plot(f,abs(fftshift(U)))axis([-600 600 0 0.04])xlabel('频率f')title('已调信号的频谱图')subplot(3,2,5)plot(t,m(1:length(t)))axis([-0.1 0.1 -1 1])xlabel('时间t')title('原调制信号的时域图')subplot(3,2,6)plot(t,dem(1:length(t)))axis([-0.1 0.1 -1 1])xlabel('时间t')title('解调后信号的时域波形')。