西安电子科技大学数字信号处理大作业
西安电子科技大学数字信号处理上机作业
数字信号处理MATLAB上机作业M 2.21.题目The square wave and the sawtooth wave are two periodic sequences as sketched in figure ing the function stem. The input data specified by the user are: desired length L of the sequence, peak value A, and the period N. For the square wave sequence an additional user-specified parameter is the duty cycle, which is the percent of the period for which the signal is positive. Using this program generate the first 100 samples of each of the above sequences with a sampling rate of 20 kHz ,a peak value of 7, a period of 13 ,and a duty cycle of 60% for the square wave.2.程序% 用户定义各项参数参数A = input('The peak value =');L = input('Length of sequence =');N = input('The period of sequence =');FT = input('The desired sampling frequency =');DC = input('The square wave duty cycle = ');% 产生所需要的信号t = 0:L-1;T = 1/FT;x = A*sawtooth(2*pi*t/N);y = A*square(2*pi*(t/N),DC);% Plotsubplot(2,1,1)stem(t,x);ylabel('幅度');xlabel('n');subplot(2,1,2)stem(t,y);ylabel('幅度');xlabel('n');3.结果4.结果分析M 2.41.题目(a)Write a matlab program to generate a sinusoidal sequence x[n]= Acos(ω0 n+Ф) and plot thesequence using the stem function. The input data specified by the user are the desired length L, amplitude A, the angular frequency ω0 , and the phase Фwhere 0<ω0 <pi and 0<=Ф<=2pi. Using this program generate the sinusoidal sequences shown in figure 2.15. (b)Generate sinusoidal sequences with the angular frequencies given in Problem 2.22.Determine the period of each sequence from the plot and verify the result theoretically. 2.程序%用户定义的参数L = input('Desired length = ');A = input('Amplitude = ');omega = input('Angular frequency = ');phi = input('Phase = ');%信号产生n = 0:L-1;x = A*cos(omega*n + phi);stem(n,x);xlabel('n');ylabel('幅度');title(['\omega_{o} = ',num2str(omega)]);3.结果(a)ω0=0ω0=0.1πω0=0.8πω0=1.2π(b)ω0=0.14πω0=0.24πω0=0.34πω0=0.68πω0=0.75π4.结果分析M 2.51.题目Generate the sequences of problem 2.21(b) to 2.21(e) using matlab.2.程序(b)n = 0 : 99;x=sin(0.6*pi*n+0.6*pi);stem(n,x);xlabel('n');ylabel('幅度');(c)n = 0 : 99;x=2*cos(1.1*pi*n-0.5*pi)+2*sin(0.7*pi*n);stem(n,x);xlabel('n');ylabel('幅度');(d)n = 0 : 99;x=3*sin(1.3*pi*n-4*cos(0.3*pi*n+0.45*pi));stem(n,x);xlabel('n');ylabel('幅度');(e)n = 0 : 99;x=5*sin(1.2*pi*n+0.65*pi)+4*sin(0.8*pi*n)-cos(0.8*pi*n);stem(n,x);xlabel('n');ylabel('幅度');(f)n = 0 : 99;x=mod(n,6);stem(n,x);xlabel('n');ylabel('幅度');3.结果(b)(c)(d)(e)(f)4.结果分析M 2.61.题目Write a matlab program to plot a continuous-time sinusoidal signal and its sampled version and verify figure 2.19. You need to use the hold function to keep both plots.2.程序%用户定义的参数fo = input('Frequency of sinusoid in Hz = ');FT = input('Samplig frequency in Hz = ');%产生信号t = 0:0.001:1;g1 = cos(2*pi*fo*t);plot(t,g1,'-')xlabel('时间t');ylabel('幅度')holdn = 0:1:FT;gs = cos(2*pi*fo*n/FT);plot(n/FT,gs,'o');hold off3.结果4.结果分析M 3.11.题目Using program 3_1 determine and plot the real and imaginary parts and the magnitude and phase spectra of the following DTFT for various values of r and θ:G(e jω)=1, 0<r<1.1−2r(cosθ)e−jω+r2e−2jω2.程序%program 3_1%discrete-time fourier transform computatition%k=input('Number of frequency points = ');num=input('Numerator coefficients= ');den=input('Denominator coefficients= ');%computer the frequency responsew=0:pi/k:pi;h=freqz(num,den,w);%plot the frequency responsesubplot(221)plot(w/pi,real(h));gridtitle('real part')xlabel('\omega/\pi');ylabel('Amplitude') subplot(222)plot(w/pi,imag(h));gridtitle('imaginary part')xlabel('\omega/\pi');ylabel('Amplitude') subplot(223)plot(w/pi,abs(h));gridtitle('magnitude spectrum')xlabel('\omega/\pi');ylabel('magnitude') subplot(224)plot(w/pi,angle(h));gridtitle('phase spectrum')xlabel('\omega/\pi');ylabel('phase,radians')3.结果(a)r=0.8 θ=π/6(b)r=0.6 θ=π/34.结果分析M 3.41.题目Using matlab verify the following general properties of the DTFT as listed in Table 3.2:(a) Linearity, (b) time-shifting, (c) frequency-shifting, (d) differentiation-in-frequency, (e) convolution, (f) modulation, and (g) Parseval’s relation. Since all data in matlab have to be finite-length vectors, the sequences to be used to verify the properties are thus restricted to be of finite length.2.程序%先定义两个信号N = input('The length of the sequence = ');k = 0:N-1;%g为正弦信号g = 2*sin(2*pi*k/(N/2));%h为余弦信号h = 3*cos(2*pi*k/(N/2));[G,w] = freqz(g,1);[H,w] = freqz(h,1);%*************************************************************************%% 线性性质alpha = 0.5;beta = 0.25;y = alpha*g+beta*h;[Y,w] = freqz(y,1);figure(1);subplot(211),plot(w/pi,abs(Y));xlabel('\omega/\pi');ylabel('|Y(e^j^\omega)|');title('线性叠加后的频率特性');grid;% 画出Y 的频率特性subplot(212),plot(w/pi,alpha*abs(G)+beta*abs(H));xlabel('\omega/\pi');ylabel('\alpha|G(e^j^\omega)|+\beta|H(e^j^\omega)|');title('线性叠加前的频率特性');grid;% 画出alpha*G+beta*H 的频率特性%*************************************************************************% % 时移性质n0 = 10;%时移10个的单位y2 = [zeros([1,n0]) g];[Y2,w] = freqz(y2,1);G0 = exp(-j*w*n0).*G;figure(2);subplot(211),plot(w/pi,abs(G0));xlabel('\omega/\pi');ylabel('|G0(e^j^\omega)|');title('G0的频率特性');grid;% 画出G0的频率特性subplot(212),plot(w/pi,abs(Y2));xlabel('\omega/\pi');ylabel('|Y2(e^j^\omega)|');title('Y2的频率特性');grid;% 画出Y2 的频率特性%*************************************************************************% % 频移特性w0 = pi/2; % 频移pi/2r=256; %the value of w0 in terms of number of samplesk = 0:N-1;y3 = g.*exp(j*w0*k);[Y3,w] = freqz(y3,1);% 对采样的512个点分别进行减少pi/2,从而生成G(exp(w-w0))k = 0:511;w = -w0+pi*k/512;G1 = freqz(g,1,w);figure(3);subplot(211),plot(w/pi,abs(Y3));xlabel('\omega/\pi');ylabel('|Y3(e^j^\omega)|');title('Y3的频率特性');grid;% 画出Y3的频率特性subplot(212),plot(w/pi,abs(G1));xlabel('\omega/\pi');ylabel('|G1(e^j^\omega)|');title('G1的频率特性');grid;% 画出G1 的频率特性%*************************************************************************% % 频域微分k = 0:N-1;y4 = k.*g;[Y4,w] = freqz(y4,1);%在频域进行微分y0 = ((-1).^k).*g;G2 = [G(2:512)' sum(y0)]';delG = (G2-G)*512/pi;figure(4);subplot(211),plot(w/pi,abs(Y4));xlabel('\omega/\pi');ylabel('|Y4(e^j^\omega)|');title('Y4的频率特性');grid;% 画出Y4的频率特性subplot(212),plot(w/pi,abs(delG));xlabel('\omega/\pi');ylabel('|delG(e^j^\omega)|');title('delG的频率特性');grid;% 画出delG的频率特性%*************************************************************************% % 相乘性质y5 = conv(g,h);%时域卷积[Y5,w] = freqz(y5,1);figure(5);subplot(211),plot(w/pi,abs(Y5));xlabel('\omega/\pi');ylabel('|Y5(e^j^\omega)|');title('Y5的频率特性');grid;% 画出Y5的频率特性subplot(212),plot(w/pi,abs(G.*H));%频域乘积xlabel('\omega/\pi');ylabel('|G.*H(e^j^\omega)|');title('G.*H的频率特性');grid;% 画出G.*H的频率特性%*************************************************************************% % 帕斯瓦尔定理y6 = g.*h;%对于freqz函数,在0到2pi直接取样[Y6,w] = freqz(y6,1,512,'whole');[G0,w] = freqz(g,1,512,'whole');[H0,w] = freqz(h,1,512,'whole');% Evaluate the sample value at w = pi/2% and verify with Y6 at pi/2H1 = [fliplr(H0(1:129)') fliplr(H0(130:512)')]';val = 1/(512)*sum(G0.*H1);% Compare val with Y6(129) i.e sample at pi/2 % Can extend this to other points similarly% Parsevals theoremval1 = sum(g.*conj(h));val2 = sum(G0.*conj(H0))/512;% Comapre val1 with val23.结果(a)(b)(c)(d)(e)4.结果分析M 3.81.题目Using matlab compute the N-point DFTs of the length-N sequences of Problem 3.12 for N=3, 5, 7, and 10. Compare your results with that obtained by evaluating the DTFTs computed in Problem 3.12 at ω= 2pik/N, k=0, 1,……N-1.2.程序%用户定义N的长度N = input('The value of N = ');k = -N:N;y1 = ones([1,2*N+1]);w = 0:2*pi/255:2*pi;Y1 = freqz(y1, 1, w);%对y1做傅里叶变换Y1dft = fft(y1);k = 0:1:2*N;plot(w/pi,abs(Y1),k*2/(2*N+1),abs(Y1dft),'o');grid;xlabel('归一化频率');ylabel('幅度');(a)clf;N = input('The value of N = ');k = -N:N;y1 = ones([1,2*N+1]);w = 0:2*pi/255:2*pi;Y1 = freqz(y1, 1, w);Y1dft = fft(y1);k = 0:1:2*N;plot(w/pi,abs(Y1),k*2/(2*N+1),abs(Y1dft),'o');xlabel('Normalized frequency');ylabel('Amplitude');(b)%用户定义N的长度N = input('The value of N = ');k = -N:N;y1 = ones([1,2*N+1]);y2 = y1 - abs(k)/N;w = 0:2*pi/255:2*pi;Y2 = freqz(y2, 1, w);%对y1做傅里叶变换Y2dft = fft(y2);k = 0:1:2*N;plot(w/pi,abs(Y2),k*2/(2*N+1),abs(Y2dft),'o');grid;xlabel('归一化频率');ylabel('幅度');(c)%用户定义N的长度N = input('The value of N = ');k = -N:N;y3 =cos(pi*k/(2*N));w = 0:2*pi/255:2*pi;Y3 = freqz(y3, 1, w);%对y1做傅里叶变换Y3dft = fft(y3);k = 0:1:2*N;plot(w/pi,abs(Y3),k*2/(2*N+1),abs(Y3dft),'o');grid;xlabel('归一化频率');ylabel('幅度');3.结果(a)N=3N=5 N=7N=10 (b)N=3N=5 N=7N=10 (c)N=3N=5 N=7N=104.结果分析M 3.191.题目Using Program 3_10 determine the z-transform as a ratio of two polynomials in z-1 from each of the partial-fraction expansions listed below:(a)X1(z)=−2+104+z−1−82+z−1,|z|>0.5,(b)X2(z)=3.5−21−0.5z−1−3+z−11−0.25z−2,|z|>0.5,(c)X3(z)=5(3+2z−1)2−43+2z−1+31+0.81z−2,|z|>0.9,(d)X4(z)=4+105+2z−1+z−16+5z−1+z−2,|z|>0.5.2.程序% Program 3_10% Partical-Fraction Expansion to rational z-Transform %r = input('Type in the residues = ');p = input('Type in the poles = ');k = input('Type in the constants = ');[num, den] = residuez(r,p,k);disp('Numberator polynominal coefficients');disp(num) disp('Denominator polynomial coefficients'); disp(den)4.结果分析M 4.61.题目Plot the magnitude and phase responses of the causal IIR digital transfer functionH(z)=0.0534(1+z−1)(1−1.0166z−1+z−2) (1−0.683z−1)(1−1.4461z−1+0.7957z−2).What type of filter does this transfer function represent? Determine the difference equation representation of the above transfer function.2.程序b=[0.0534 -0.00088644 -0.00088644 0.0534];a=[1 -2.1291 1.7833863 -0.5434631];figure(1)freqz(b,a);figure(2)[H,w]=freqz(b,a);plot(w/pi,abs(H)),grid;xlabel('Normalized Frequency (\times\pi rad/sample)'),ylabel('Magnitude');幅度化成真值之后:4.结果分析H(z)=0.0534−0.00088644z−1−0.00088644z−2+0.0534z−31−2.1291z−1+1.7833863z−2−0.5434631z−3M 4.71.题目Plot the magnitude and phase responses of the causal IIR digital transfer functionH(z)=(1−z−1)4(1−1.499z−1+0.8482z−2)(1−1.5548z−1+0.6493z−2).2.程序b=[1 -4 6 -4 1];a=[1 -3.0538 3.8227 -2.2837 0.5472]; figure(1)freqz(b,a);figure(2)[H,w]=freqz(b,a);plot(w/pi,abs(H)),grid;xlabel('Normalized Frequency (\times\pi rad/sample)'), ylabel('Magnitude');3.结果4.结果分析。
数字信号处理西电
数字信号处理上机第一次实验实验一:设给定模拟信号()1000t a x t e -=,的单位是ms 。
(1) 利用MATLAB 绘制出其时域波形和频谱图(傅里叶变换),估计其等效带宽(忽略谱分量降低到峰值的3%以下的频谱)。
(2) 用两个不同的采样频率对给定的进行采样。
○1。
○2 。
比较两种采样率下的信号频谱,并解释。
实验一MATLAB 程序:(1)○1 clc; fs=5000;ts=1/fs;N=1000;t=(-N:N)*ts;s=exp(-abs(t));plot(t,s,'linewidth',1.5)xlabel('时间')ylabel('幅度')set(gca,'fontweight','b','fontsize',12)SPL=N*100;figuresp=fftshift(fft(s,SPL));sp=sp/max(sp)*100;freqb=-fs/2:fs/SPL:fs/2-fs/SPL;plot(freqb,abs(sp))xlabel('频率')ylabel('频谱幅度')set(gca,'fontweight','b','fontsize',12)yy=abs(abs(sp)-3);[aa,freqind]=min(yy);(freqind-SPL/2)*fs/SPLt ()a x t ()()15000s a f x t x n =以样本秒采样得到。
()()11j x n X e ω画出及其频谱()()11000s a f x t x n =以样本得到。
()()11j x n X e ω画出及其频谱○2 clc;fs=1000;ts=1/fs;N=1000;t=(-N:N)*ts;s=exp(-abs(t));plot(t,s,'linewidth',1.5)xlabel('时间')ylabel('幅度')set(gca,'fontweight','b','fontsize',12) SPL=N*100;figuresp=fftshift(fft(s,SPL));sp=sp/max(sp)*100;freqb=-fs/2:fs/SPL:fs/2-fs/SPL;plot(freqb,abs(sp))xlabel('频率')ylabel('频谱幅度')set(gca,'fontweight','b','fontsize',12)yy=abs(abs(sp)-3);[aa,freqind]=min(yy);(freqind-SPL/2)*fs/SPL实验三:设,,编写MATLAB 程序,计算:(1) 5点圆周卷积;(2) 6点圆周卷积;(3) 线性卷积;(4) 画出的,和时间轴对齐。
西安电子科技大学数字信号处理大作业
数字信号处理大作业班级:021231学号:姓名:指导老师:吕雁一写出奈奎斯特采样率和和信号稀疏采样的学习报告和体会1、采样定理在进行A/D信号的转换过程中,当采样频率fs.max大于信号中最高频率fmax的2倍时(fs.max>2fmax),采样之后的数字信号完整地保留了原始信号中的信息,一般实际应用中保证采样频率为信号最高频率的5~10倍;采样定理又称奈奎斯特定理。
(1)在时域频带为F的连续信号 f(t)可用一系列离散的采样值f(t1),f(t1±Δt),f(t1±2Δt),...来表示,只要这些采样点的时间间隔Δt≤1/2F,便可根据各采样值完全恢复原始信号。
(2)在频域当时间信号函数f(t)的最高频率分量为fmax时,f(t)的值可由一系列采样间隔小于或等于1/2fo的采样值来确定,即采样点的重复频率fs ≥2fmax。
2、奈奎斯特采样频率(1)概述奈奎斯特采样定理:要使连续信号采样后能够不失真还原,采样频率必须大于信号最高频率的两倍(即奈奎斯特频率)。
奈奎斯特频率(Nyquist frequency)是离散信号系统采样频率的一半,因哈里·奈奎斯特(Harry Nyquist)或奈奎斯特-香农采样定理得名。
采样定理指出,只要离散系统的奈奎斯特频率高于被采样信号的最高频率或带宽,就可以真实的还原被测信号。
反之,会因为频谱混叠而不能真实还原被测信号。
采样定理指出,只要离散系统的奈奎斯特频率高于采样信号的最高频率或带宽,就可以避免混叠现象。
从理论上说,即使奈奎斯特频率恰好大于信号带宽,也足以通过信号的采样重建原信号。
但是,重建信号的过程需要以一个低通滤波器或者带通滤波器将在奈奎斯特频率之上的高频分量全部滤除,同时还要保证原信号中频率在奈奎斯特频率以下的分量不发生畸变,而这是不可能实现的。
在实际应用中,为了保证抗混叠滤波器的性能,接近奈奎斯特频率的分量在采样和信号重建的过程中可能会发生畸变。
电子科技大学 数字信号处理 第二次编程作业
第二次编程作业加载信号(EEGdata.txt);分析信号的幅度谱,确定低通FIR数字滤波器的指标;分别利用各种窗函数(Rectangular, Hamming, Kaiser)设计此FIR滤波器;对信号加随机噪声,并用设计的FIR滤波器对含噪声信号进行滤波。
要求:1)给出程序,画出滤波器幅度谱及其增益图(dB);分析设计的滤波器是否达到指标。
2)利用设计的三种滤波器对加载的信号分别进行滤波,对比分析滤波结果。
程序:data = load('EEGdata.txt');n=1:length(data);xn=0.1*rand(length(data),1)-0.05;x=xn+data;[h0,w0]=freqz(data);[h1,w1]=freqz(x);figure(1)plot(w0/pi,abs(h0),w1/pi,abs(h1),'r');legend('Ô-dataÐźÅƵÆ×ͼ','¼ÓÔëÐźÅƵÆ×ͼ');wp = 0.15*pi;ws = 0.2*pi;DB = ws-wp;N = ceil(2*3.11*pi/DB+1);wc = (wp+ws)/2;wc = wc/pi;win = rectwin(N+1);b = fir1(N, wc, 'low', win);figure(2)freqz(b,1);wp = 0.15*pi;ws = 0.2*pi;DB = ws-wp;N = ceil(2*3.11*pi/DB+1);wc = (wp+ws)/2;wc = wc/pi;win = hamming(N+1);b1 = fir1(N, wc, 'low', win);figure(3)freqz(b1,1);fpts = [0.15, 0.2];mag = [1, 0];dev = [0.01, 0.01];[N, Wn, beta, ftype] = kaiserord(fpts, mag, dev)Kwin = kaiser(N+1, beta);b2 = fir1(N, Wn, Kwin);[h,omega] = freqz(b2,1,512);figure(4)freqz(b2,1);y=filter(b, 1, x);figure(5)subplot(2,1,1);plot(n,data,n,x,'g',n,y,'r');title('¾ØÐδ°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ');lege nd('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy');[H,W]=freqz(y);subplot(2,1,2);plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('¾ØÐδ°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ') ;legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy');y1=filter(b1, 1, x);figure(6)subplot(2,1,1);plot(n,data,n,x,'g',n,y1,'r');title('Hamming´°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ'); legend('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy1');[H,W]=freqz(y1);subplot(2,1,2);plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('Hamming´°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy1');y2=filtfilt(b2, 1, x); figure(7) subplot(2,1,1);plot(n,data,n,x,'g',n,y2,'r');title('Kaiser´°Â˲¨·ù¶ÈÆ×Ч¹û¶Ô±Èͼ');l egend('Ô-ÐźÅdata','¼ÓÔëÐźÅx','Â˲¨ºóÐźÅy2'); [H,W]=freqz(y2); subplot(2,1,2);plot(w0/pi,abs(h0),W/pi,abs(H),'r');title('Kaiser´°Â˲¨ÆµÆ×Ч¹û¶Ô±Èͼ');legend('Ô-ÐźÅdata','Â˲¨ºóÐźÅy2');figure(8)plot(n,data,n,y,'c',n,y1,'r',n,y2,'g');title('ÈýÖÖ´°º¯ÊýÂ˲¨Ð§¹û¶Ô±Èͼ');legend('Ô-dataÐźÅ','¾ØÐδ°Â˲¨','Hamming´°Â˲¨','Kaiser´°Â˲¨');0.10.20.30.40.50.60.70.80.9102468101214161800.10.20.30.40.50.60.70.80.91-3000-2000-1000Normalized Frequency (⨯π rad/sample)P h a s e (d e g r e e s )00.10.20.30.40.50.60.70.80.91-100-5050Normalized Frequency (⨯π rad/sample)M a g n i t u d e (d B )00.10.20.30.40.50.60.70.80.91-3000-2000-1000Normalized Frequency (⨯π rad/sample)P h a s e (d e g r e e s )0.10.20.30.40.50.60.70.80.91-150-100-50050Normalized Frequency (⨯π rad/sample)M a g n i t u d e (d B )0.10.20.30.40.50.60.70.80.91-2000-1500-1000-500Normalized Frequency (⨯π rad/sample)P h a s e (d e g r e e s )0.10.20.30.40.50.60.70.80.91-150-100-50050Normalized Frequency (⨯π rad/sample)M a g n i t u d e (d B )-0.2-0.100.10.2矩形窗滤波幅度谱效果对比图0.10.20.30.40.50.60.70.80.91051015矩形窗滤波频谱效果对比图-0.2-0.100.10.2Hamming 窗滤波幅度谱效果对比图0.10.20.30.40.50.60.70.80.91051015Hamming 窗滤波频谱效果对比图-0.2-0.100.10.2Kaiser 窗滤波幅度谱效果对比图0.10.20.30.40.50.60.70.80.9105101520Kaiser 窗滤波频谱效果对比图结果分析:滤波效果皆不尽如人意,原因是加载的噪声信号是随机信号,各频率皆有,滤波器通带部分的噪声无法滤掉。
西安电子科技大学数字信号处理上机报告
数字信号处理大作业院系:电子工程学院学号:02115043姓名:邱道森实验一:信号、系统及系统响应一、实验目的(1) 熟悉连续信号经理想采样前后的频谱变化关系,加深对时域采样定理的理解。
(2) 熟悉时域离散系统的时域特性。
(3) 利用卷积方法观察分析系统的时域特性。
(4) 掌握序列傅里叶变换的计算机实现方法,利用序列的傅里叶变换对连续信号、离散信号及系统响应进行频域分析。
二、实验原理采样是连续信号数字处理的第一个关键环节。
对连续信号()a x t 进行理想采样的过程可用(1.1)式表示:()()()ˆa a xt x t p t =⋅ 其中()t xa ˆ为()a x t 的理想采样,()p t 为周期冲激脉冲,即 ()()n p t t nT δ∞=-∞=-∑()t xa ˆ的傅里叶变换()j a X Ω为 ()()s 1ˆj j j a a m X ΩX ΩkΩT ∞=-∞=-∑进行傅里叶变换,()()()j ˆj e d Ωt a a n X Ωx t t nT t δ∞∞--∞=-∞⎡⎤=-⎢⎥⎣⎦∑⎰ ()()j e d Ωtan x t t nT t δ∞∞--∞=-∞=-∑⎰()j e ΩnTan x nT ∞-=-∞=∑式中的()a x nT 就是采样后得到的序列()x n , 即()()a x n x nT =()x n 的傅里叶变换为()()j j e enn X x n ωω∞-=-∞=∑比较可知()()j ˆj e aΩTX ΩX ωω==为了在数字计算机上观察分析各种序列的频域特性,通常对()j e X ω在[]0,2π上进行M 点采样来观察分析。
对长度为N 的有限长序列()x n ,有()()1j j 0eekk N nn X x n ωω--==∑其中2π,0,1,,1k k k M Mω==⋅⋅⋅-一个时域离散线性时不变系统的输入/输出关系为()()()()()m y n x n h n x m h n m ∞=-∞=*=-∑上述卷积运算也可以转到频域实现()()()j j j e e e Y X H ωωω= (1.9)三、实验内容及步骤(1) 认真复习采样理论、 离散信号与系统、 线性卷积、 序列的傅里叶变换及性质等有关内容, 阅读本实验原理与方法。
西电DSP大作业任务报告
DSP实验课程序设计报告学院:电子工程学院学号:1202121013姓名:赵海霞指导教师:苏涛DSP实验课大作业设计一实验目的在DSP上实现线性调频信号的脉冲压缩、动目标显示(MTI)和动目标检测(MTD),并将结果与MATLAB上的结果进行误差仿真。
二实验内容2.1 MATLAB仿真设定带宽、脉宽、采样率、脉冲重复频率,用MATLAB产生16个脉冲的LFM,每个脉冲有4个目标(静止,低速,高速),依次做2.1.1 脉压2.1.2 相邻2脉冲做MTI,产生15个脉冲2.1.3 16个脉冲到齐后,做MTD,输出16个多普勒通道2.2 DSP实现将MATLAB产生的信号,在visual dsp中做脉压,MTI、MTD,并将结果与MATLAB作比较。
三实验原理3.1 线性调频线性调频脉冲压缩体制的发射信号其载频在脉冲宽度内按线性规律变化即用对载频进行调制(线性调频)的方法展宽发射信号的频谱,在大时宽的前提下扩展了信号的带宽。
若线性调频信号中心频率为f,脉宽为τ,带宽为B,幅度为A,μ为调频斜率,则其表达式如下:]212cos[)()(20t t f t rect A t x μπτ+••=;)(为矩形函数rect 在相参雷达中,线性调频信号可以用复数形式表示,即)]212(exp[)()(20t t f j t rect A t x μπτ+••= 在脉冲宽度内,信号的角频率由220μτπ-f 变化到220μτπ+f 。
3.2 脉冲压缩原理脉冲雷达信号发射时,脉冲宽度τ决定着雷达的发射能量,发射能量越大, 作用距离越远;在传统的脉冲雷达信号中,脉冲宽度同时还决定着信号的频率宽度B ,即带宽与时宽是一种近似倒数的关系。
脉冲越宽,频域带宽越窄,距离分辨率越低。
脉冲压缩的主要目的是为了解决信号的作用距离和信号的距离分辨率之间的矛盾。
为了提高信号的作用距离,我们就需要提高信号的发射功率,因此,必须提高发射信号的脉冲宽度,而为了提高信号的距离分辨率,又要求降低信号的脉冲宽度。
西安电子科大数字信号处理课后实验答案
一系统响应及稳定性的实验报告一. 实验目的:(1)掌握 求系统响应的方法。
(2)掌握时域离散系统的时域特性。
(3)分析、观察及检验系统的稳定性。
二. 实验原理与方法:1.在时域中,描写系统特性的方法是差分方程和单位脉冲响应,在频域可以用系统函数描述系统特性。
已知输入信号可以由差分方程、单位脉冲响应或系统函数求出系统对于该输入信号的响应。
在计算机上可用filter 函数求差分方程的解, conv 函数计算输入信号和系统的单位脉冲响应的线性卷积,求出系统的响应。
2.系统的时域特性指的是系统的线性时不变性质、因果性和稳定性。
3.系统的稳定性是指对任意有界的输入信号,系统都能得到有界的系统响应。
或者系统的单位脉冲响应满足绝对可和的条件。
系统的稳定性由其差分方程的系数决定。
实际中检查系统是否稳定,不可能检查系统对所有有界的输入信号,输出是否都是有界输出,或者检查系统的单位脉冲响应满足绝对可和的条件。
可行的方法是在系统的输入端加入单位阶跃序列,如果系统的输出趋近一个常数(包括零),就可以断定系统是稳定的[19]。
系统的稳态输出是指当∞→n 时,系统的输出。
如果系统稳定,信号加入系统后,系统输出的开始一段称为暂态效应,随n 的加大,幅度趋于稳定,达到稳态输出。
三.实验内容及步骤:1.编制程序,包括产生输入信号、单位脉冲响应序列的子程序,用filter 函数或conv 函数求解系统输出响应的主程序。
程序中要有绘制信号波形的功能。
2.给定一个低通滤波器的差分方程为)1(9.0)1(05.0)(05.0)(-+-+=n y n x n x n y a) 分别求出系统对)()(81n R n x =和)()(2n u n x =的响应序列,并画出其波形。
b) 求出系统的单位冲响应,画出其波形。
3.给定系统的单位脉冲响应为)()(101n R n h =)3()2(5.2)1(5.2)()(2-+-+-+=n n n n n h δδδδ 用线性卷积法分别求系统h 1(n)和h 2(n)对)()(81n R n x =的输出响应,并画出波形。
西电数字信号处理大作业
西电数字信号处理⼤作业第⼆章2.25 已知线性时不变系统的差分⽅程为若系统的输⼊序列x(x)={1,2,3,4,2,1}编写利⽤递推法计算系统零状态响应的MATLAB程序,并计算出结果。
代码及运⾏结果:>> A=[1,-0.5];>> B=[1,0,2];>> n=0:5;>> xn=[1,2,3,4,2,1];>> zx=[0,0,0];zy=0;>> zi=filtic(B,A,zy,zx);>> yn=filter(B,A,xn,zi);>> figure(1)>> stem(n,yn,'.');>> grid on;2.28图所⽰系统是由四个⼦系统T1、T2、T3和T4组成的,分别⽤单位脉冲响应或差分⽅程描述为T1:其他T2:其他T3:T4:编写计算整个系统的单位脉冲响应h(n),0≤n≤99的MATLAB程序,并计算结果。
代码及结果如下:>> a=0.25;b=0.5;c=0.25;>> ys=0;>> xn=[1,zeros(1,99)];>> B=[a,b,c];>> A=1;>> xi=filtic(B,A,ys);>> yn1=filter(B,A,xn,xi);>> h1=[1,1/2,1/4,1/8,1/16,1/32]; >> h2=[1,1,1,1,1,1];>> h3=conv(h1,h2);>> h31=[h3,zeros(1,89)]; >> yn2=yn1+h31;>> D=[1,1];C=[1,-0.9,0.81]; >> xi2=filtic(D,C,yn2,xi); >> xi2=filtic(D,C,ys);>> yn=filter(D,C,yn2,xi); >> n=0:99;>> figure(1)>> stem(n,yn,'.');>> title('单位脉冲响应'); >> xlabel('n');ylabel('yn');2.30 利⽤MATLAB画出受⾼斯噪声⼲扰的正弦信号的波形,表⽰为其中v(n)是均值为零、⽅差为1的⾼斯噪声。
西安电子科技大学-数字信号处理-试卷C答案
Answer to “Digital Signal Processing of 2005”Problem 1(a) even part: };5.0,1,7,7,5,7,7,1,5.0{---=e X odd part: };5.0,1,3,1,0,1,3,1,5.0{----=o X(b) };20,16,11,94,36,40,31,16,12,0{-----=y (c) MATLAB Programn=-4:2;x=[1 -2 4 6 -5 8 10]; [x11,n11]=sigshift(x,n,2); [x12,n12]=sigshift(x,n,-1); [x13,n13]=sigfold(x,n); [x13,n13]=sigshift(x13,n13,-2); [x12,n12]=sigmult(x,n,x12,n12); [y,n]=sigadd(2*x11,n11,x12,n12); [y,n]=sigadd(y,n,-1*x13,n13)Problem 2(a)w j w j w j w j jw jw e e e e e e X 65424210124)(-----++++++=,()j X e ωis periodic in ω with period 2π(b) MATLAB Program :clear; close all;n = 0:6; x = [4,2,1,0,1,2,4]; w = [0:1:1000]*pi/1000;X = x*exp(-j*n'*w); magX = abs(X); phaX = angle(X); % Magnitude Response Plotsubplot(2,1,1); plot(w/pi,magX);grid;xlabel('frequency in pi units'); ylabel('|X|'); title('Magnitude Response'); % Phase response plotsubplot(2,1,2); plot(w/pi,phaX*180/pi);grid;xlabel('frequency in pi units'); ylabel('Degrees'); title('Phase Response'); axis([0,1,-180,180])(c) Because the given sequence x (n)={4,2,1,0,1,2,4} (n=0,1,2,3,4,5,6) is symmetric about 132N α-==,the phase response ()j H e ω< satisfied the condition :()3j H e ωαωω<=-=- so the phase response is a linear function in ω.(d) 150,350Hz Hz Ω=-;(e) The difference of amplitude and magnitude response:Firstly, the amplitude response is a real function, and it may be both positive and negative. The magnitude response is always positive.Secondly, the phase response associated with the magnitude response is a discontinuous function. While the associated with the amplitude is a continuous linear function.Problem 3(a) )9.09.01/()1()(211------=z z z z HZero:0 and 1; Pole:-0.6 and 1.5; (b)1116151()212110.61 1.5H z z z--=⨯+⨯+-, 165()((0.6)(1.5))()2121n nh n u n =-+ (c) ROC : 0.6 1.5Z <<,()()()163531215212n nh n u n u n ⎛⎫⎛⎫=⨯--⨯--- ⎪ ⎪⎝⎭⎝⎭ Problem 4(a) y(n)={50,44,34,52};(b) y(n)={5,16,34,52,45,28,0}; (c) N=6;(d) MATLAB Program :Function y=circonv(x1,x2,N) If (length(x1)>N)error(“N must not be smaller than the length of sequence ”) elsex1=[x1,zeros(1,N-length(x1))]; endif(length(x2)>N)error(“N must not be smaller than the length of sequence ”)elsex2=[x2,zeros(1,N-length(x2))]; endy1=dft(x1,N).*dft(x2,N); y=idft(y,N);(e) DTFT is discrete in time domain, but continuous in frequency domain. The DFT is discrete both in time and frequency domain.The FFT is a very efficient method for calculating DFT.Problem 5(a) Direct form II uses the little delay and it can decrease the space of the compute. (b)The advantage of the linear-phase form:1. For frequency-selective filters, linear-phase structure is generally desirable to have a phase-responsethat is a linear function of frequency.2. This structure requires 50% fewer multiplications than the direct form. (c) Block diagrams are shown as under:1z -1z -1z -1z -1z -1z -()x n )n()x n1-1-Problem 6(a) we use Hamming or Blackman window to design the bandpass filter because it can provide us attenuationexceed 60dB .(b) According to Blackman window :first, Determine transition width =p s W W - ;second, Determine the type of the window according to s A ;third, Compute M using the formula MW W p s π2=- ; fourth, Compute ideal LPF2sp W W Wc +=;fifth, design the window needed, multiply point by point; sixth, determines p A R ,(c) MATLAB Program :%% Specifications about Blackman window: ws1 = 0.2*pi; % lower stopband edge wp1 = 0.3*pi; % lower passband edge wp2 = 0.6*pi; % upper passband edge ws2 = 0.7*pi; % upper stopband edge Rp = 0.5; % passband ripple As = 60; % stopband attenuation %tr_width = min((wp1-ws1),(ws2-wp2));M = ceil(6.6*pi/tr_width); M = 2*floor(M/2)+1, % choose odd M n = 0:M-1;w_ham = (hamming(M))';wc1 = (ws1+wp1)/2; wc2 = (ws2+wp2)/2; hd = ideal_lp(wc2,M)-ideal_lp(wc1,M); h = hd .* w_ham;[db,mag,pha,grd,w] = freqz_m(h,1); delta_w = pi/500;Asd = floor(-max(db([1:floor(ws1/delta_w)+1]))), % Actual AttnRpd = -min(db(ceil(wp1/delta_w)+1:floor(wp2/delta_w)+1)), % Actual passband ripple (5) %%% Filter Response Plotssubplot(2,2,1); stem(n,hd); title('Ideal Impulse Response: Bandpass'); axis([-1,M,min(hd)-0.1,max(hd)+0.1]); xlabel('n'); ylabel('hd(n)') set(gca,'XTickMode','manual','XTick',[0;M-1],'fontsize',10) subplot(2,2,2); stem(n,w_ham); title('Hamming Window'); axis([-1,M,-0.1,1.1]); xlabel('n'); ylabel('w_ham(n)')set(gca,'XTickMode','manual','XTick',[0;M-1],'fontsize',10) set(gca,'YTickMode','manual','YTick',[0;1],'fontsize',10)subplot(2,2,3); stem(n,h); title('Actual Impulse Response: Bandpass'); axis([-1,M,min(hd)-0.1,max(hd)+0.1]); xlabel('n'); ylabel('h(n)') set(gca,'XTickMode','manual','XTick',[0;M-1],'fontsize',10)subplot(2,2,4); plot(w/pi,db); title('Magnitude Response in dB');axis([0,1,-As-30,5]); xlabel('frequency in pi units'); ylabel('Decibels') set(gca,'XTickMode','manual','XTick',[0;0.3;0.4;0.5;0.6;1])set(gca,'XTickLabelMode','manual','XTickLabels',['0';'0.3';'0.4';'0.5';'0.6';'1'],... 'fontsize',10)set(gca,'TickMode','manual','YTick',[-50;0])set(gca,'YTickLabelMode','manual','YTickLabels',['-50';'0']);gridProblem 7Firstly, we use the given specifications ofs p s p A R ,,,ωω to design an analog lowpass IIR filter.Secondly, we change the analog lowpass IIR filter into the analog highpass IIR filter. Thirdly, we change the analog highpass IIR filter into the digital highpass IIR filter.。
电子科技大学20春《数字信号处理》在线作业3.doc
1.题目及选项如下图所示A.AB.BC.CD.D【参考答案】: A2.下列关于FFT的说法中错误的是()。
A.FFT是一种新的变换B.FFT是DFT的快速算法C.FFT基本上可以分成按时间抽取法和按频率抽取法两类D.基2 FFT要求序列的点数为(其中为整数)【参考答案】: A3.下列关于FIR滤波器的说法中正确的是()。
A.FIR滤波器难以设计成线性相位特性B.FIR滤波器的单位冲激响应长度是无限的C.FIR滤波器的单位冲激响应长度是有限的D.对于相同的幅频特性要求,用FIR滤波器实现要比用IIR滤波器实现阶数低【参考答案】: C4.题目及选项如下图所示A.AB.BC.CD.D【参考答案】: A5.题目及选项如下图所示A.AB.BC.CD.D 【参考答案】: A6.题目及选项如下图所示A.AB.BC.CD.D 【参考答案】: A7.题目及选项如下图所示A.AB.BC.CD.D 【参考答案】: A8.题目及选项如下图所示A.AB.BC.CD.D【参考答案】: A9.在IIR数字滤波器的设计中,需将模拟参考滤波器从S平面到Z平面作单值映射,应采用()。
A.脉冲响应不变法B.双线性变换法C.窗函数法D.频率采样法【参考答案】: B10.已知一连续时间信号的最高截止频率是4kHz,则要从该连续时间信号的均匀离散采样值无失真地恢复原信号,则采样频率为()。
A.2kHzB.4kHzC.6kHzD.8kHz【参考答案】: D11.下列关于窗函数设计法的说法中错误的是()。
A.窗函数的长度增加,则主瓣宽度减小,旁瓣宽度减小B.窗函数的旁瓣相对幅度取决于窗函数的形状,与窗函数的长度无关C.为减小旁瓣相对幅度而改变窗函数的形状,通常主瓣的宽度会增加D.对于长度固定的窗,只要选择合适的窗函数就可以使主瓣宽度足够窄,旁瓣幅度足够小【参考答案】: B12.已知某序列的Z变换收敛域为全Z平面,则该序列为()。
A.有限长序列B.双边序列C.左边序列D.右边序列【参考答案】: A13.题目及选项如下图所示A.AB.BC.CD.D【参考答案】: A14.使用FFT算法计算两有限长序列的线性卷积,则需要调用()次FFT算法。
数字信号处理实验报告(西电)
数字信号处理实验报告班级:****姓名:郭**学号:*****联系方式:*****西安电子科技大学电子工程学院绪论数字信号处理起源于十八世纪的数学,随着信息科学和计算机技术的迅速发展,数字信号处理的理论与应用得到迅速发展,形成一门极其重要的学科。
当今数字信号处理的理论和方法已经得到长足的发展,成为数字化时代的重要支撑,其在各个学科和技术领域中的应用具有悠久的历史,已经渗透到我们生活和工作的各个方面。
数字信号处理相对于模拟信号处理具有许多优点,比如灵活性好,数字信号处理系统的性能取决于系统参数,这些参数很容易修改,并且数字系统可以分时复用,用一套数字系统可以分是处理多路信号;高精度和高稳定性,数字系统的运算字符有足够高的精度,同时数字系统不会随使用环境的变化而变化,尤其使用了超大规模集成的DSP 芯片,简化了设备,更提高了系统稳定性和可靠性;便于开发和升级,由于软件可以方便传送,复制和升级,系统的性能可以得到不断地改善;功能强,数字信号处理不仅能够完成一维信号的处理,还可以试下安多维信号的处理;便于大规模集成,数字部件具有高度的规范性,对电路参数要求不严格,容易大规模集成和生产。
数字信号处理用途广泛,对其进行一系列学习与研究也是非常必要的。
本次通过对几个典型的数字信号实例分析来进一步学习和验证数字信号理论基础。
实验一主要是产生常见的信号序列和对数字信号进行简单处理,如三点滑动平均算法、调幅广播(AM )调制高频正弦信号和线性卷积。
实验二则是通过编程算法来了解DFT 的运算原理以及了解快速傅里叶变换FFT 的方法。
实验三是应用IRR 和FIR 滤波器对实际音频信号进行处理。
实验一●实验目的加深对序列基本知识的掌握理解●实验原理与方法1.几种常见的典型序列:0()1,00,0(){()()(),()sin()j n n n n u n x n Aex n a u n a x n A n σωωϕ+≥<====+单位阶跃序列:复指数序列:实指数序列:为实数 正弦序列:2.序列运算的应用:数字信号处理中经常需要将被加性噪声污染的信号中移除噪声,假定信号 s(n)被噪声d(n)所污染,得到了一个含噪声的信号()()()x n s n d n =+。
《数字信号处理》第三版课后实验答案_西安电子科技大学出版社
程序清单及波形显示: clc;close all;clear all;%======内容1:调用filter 解差分方程,由系统对u(n)的响应判断稳定性======A=[1,-0.9];B=[0.05,0.05]; %系统差分方程系数向量B 和A x1n=[1 1 1 1 1 1 1 1 zeros(1,50)]; %产生信号x1(n)=R8(n) x2n=ones(1,128); %产生信号x2(n)=u(n)y1n=filter(B,A,x1n); %求系统对x1(n)的响应y1(n) n=0:length(y1n)-1;subplot(2,2,1);stem(n,y1n,'.'); title('(a) 系统对R8(n)的响应y1(n)');xlabel('n');ylabel('y1(n)');y2n=filter(B,A,x2n); %求系统对x2(n)的响应y2(n) n=0:length(y2n)-1;subplot(2,2,2);stem(n,y2n,'.'); title('(b) 系统对u(n)的响应y2(n)');xlabel('n');ylabel('y2(n)');hn=impz(B,A,58); %求系统单位脉冲响应h(n) n=0:length(hn)-1;subplot(2,2,3);y=hn;stem(n,hn,'.'); title('(c) 系统单位脉冲响应h(n)');xlabel('n');ylabel('h(n)');(a) 系统对R8(n)的响应y1(n)ny 1(n )(b) 系统对u(n)的响应y2(n)ny 2(n )(c) 系统单位脉冲响应h(n)nh (n )%===内容2:调用conv 函数计算卷积============================x1n=[1 1 1 1 1 1 1 1 ]; %产生信号x1(n)=R8(n) h1n=ones(1,10); h2n=[1 2.5 2.5 1 ];y21n=conv(h1n,x1n); y22n=conv(h2n,x1n); figure(2)n=0:length(h1n)-1;subplot(2,2,1);stem(n,h1n); title('(d) 系统单位脉冲响应h1n');xlabel('n');ylabel('h1(n)'); n=0:length(y21n)-1;subplot(2,2,2);stem(n,y21n); title('(e) h1(n)与R8(n)的卷积y21n');xlabel('n');ylabel('y21(n)'); n=0:length(h2n)-1;subplot(2,2,3);stem(n,h2n); title('(f) 系统单位脉冲响应h2n');xlabel('n');ylabel('h2(n)'); n=0:length(y22n)-1;subplot(2,2,4);stem(n,y22n); title('(g) h2(n)与R8(n)的卷积y22n');xlabel('n');ylabel('y22(n)');nh 1(n )ny 21(n )(f) 系统单位脉冲响应h2nnh 2(n)(g) h2(n)与R8(n)的卷积y22nny 22(n )%=========内容3:谐振器分析======================== un=ones(1,256); %产生信号u(n) n=0:255;xsin=sin(0.014*n)+sin(0.4*n); %产生正弦信号A=[1,-1.8237,0.9801];B=[1/100.49,0,-1/100.49]; %系统差分方程系数向量B 和Ay31n=filter(B,A,un); %谐振器对u(n)的响应y31(n) y32n=filter(B,A,xsin); %谐振器对u(n)的响应y31(n) figure(3)n=0:length(y31n)-1;subplot(2,1,1);stem(n,y31n,'.'); title('(h) 谐振器对u(n)的响应y31n');xlabel('n');ylabel('y31(n)'); n=0:length(y32n)-1;subplot(2,1,2);stem(n,y32n,'.'); title('(i) 谐振器对正弦信号的响应y32n');xlabel('n');ylabel('y32(n)');050100150200250300(h) 谐振器对u(n)的响应y31nny 31(n )(i) 谐振器对正弦信号的响应y32nny 32(n )程序清单及波形显示:% DTMF 双频拨号信号产生6位电话号码 %clear all;clc;tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68]; % DTMF 信号代表的16个数N=205;K=[18,20,22,24,31,34,38,42];f1=[697,770,852,941]; % 行频率向量 f2=[1209,1336,1477,1633]; % 列频率向量 TN=input('键入6位电话号码= '); % 输入6位数字TNr=0; %接收端电话号码初值为零for l=1:6;d=fix(TN/10^(6-l)) TN=TN-d*10^(6-l); for p=1:4; for q=1:4;if tm(p,q)==abs(d); break,end % 检测码相符的列号q endif tm(p,q)==abs(d); break,end % 检测码相符的行号p endn=0:1023; % 为了发声,加长序列x = sin(2*pi*n*f1(p)/8000) + sin(2*pi*n*f2(q)/8000);% 构成双频信号sound(x,8000); % 发出声音 pause(0.1)% 接收检测端的程序X=goertzel(x(1:205),K+1); % 用Goertzel 算法计算八点DFT样本val = abs(X); % 列出八点DFT向量subplot(3,2,l);stem(K,val,'.');grid;xlabel('k');ylabel('|X(k)|') % 画出DFT(k)幅度axis([10 50 0 120])limit = 80; %for s=5:8;if val(s) > limit, break, end % 查找列号endfor r=1:4;if val(r) > limit, break, end % 查找行号endTNr=TNr+tm(r,s-4)*10^(6-l);enddisp('接收端检测到的号码为:') % 显示接收到的字符disp(TNr)显示结果:键入6位电话号码= 123456d = 1d = 2d = 3d = 4d = 5d = 6接收端检测到的号码为:123456050100|X (k )||X (k )|050100k|X (k )|050100k|X (k )|050100k|X (k )|050100k|X (k )|% DTMF 双频拨号信号产生8位电话号码 %clear all;clc;tm=[1,2,3,65;4,5,6,66;7,8,9,67;42,0,35,68]; % DTMF 信号代表的16个数 N=205;K=[18,20,22,24,31,34,38,42];f1=[697,770,852,941]; % 行频率向量 f2=[1209,1336,1477,1633]; % 列频率向量 TN=input('键入8位电话号码= '); % 输入8位数字TNr=0; %接收端电话号码初值为零for l=1:8;d=fix(TN/10^(8-l)) TN=TN-d*10^(8-l); for p=1:4; for q=1:4;if tm(p,q)==abs(d); break,end % 检测码相符的列号q endif tm(p,q)==abs(d); break,end % 检测码相符的行号p endn=0:1023; % 为了发声,加长序列 x = sin(2*pi*n*f1(p)/8000) + sin(2*pi*n*f2(q)/8000);% 构成双频信号 sound(x,8000); % 发出声音 pause(0.1)% 接收检测端的程序X=goertzel(x(1:205),K+1); % 用Goertzel 算法计算八点DFT 样本val = abs(X); % 列出八点DFT 向量 subplot(4,2,l);stem(K,val,'.');grid;xlabel('k');ylabel('|X(k)|') % 画出DFT(k)幅度axis([10 50 0 120])limit = 80; % for s=5:8;if val(s) > limit, break, end % 查找列号 endfor r=1:4;if val(r) > limit, break, end % 查找行号 endTNr=TNr+tm(r,s-4)*10^(8-l); enddisp('接收端检测到的号码为:') % 显示接收到的字符 disp(TNr) 显示结果:键入8位电话号码= 12345678 d = 1 d = 2 d = 3 d = 4 d = 5 d = 6 d = 7 d = 8接收端检测到的号码为:12345678|X (k )||X (k )||X (k )||X (k )||X (k )||X (k )|k|X (k )|k|X (k )|程序清单及波形显示: % 时域采样理论验证程序Tp=64/1000; %观察时间Tp=64微秒Fs=1000;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,1); stem(n,xat,'.'); xlabel('n');ylabel('x1(n)'); title('(a) Fs=1000Hz');k=0:M-1;fk=k/Tp;subplot(3,2,2);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=1000Hz');xlabel('\omega/hz');ylabel('(H1(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))]);Fs=300;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,3); stem(n,xat,'.'); xlabel('n');ylabel('x2(n)'); title('(b)Fs=300Hz');k=0:M-1;fk=k/Tp;subplot(3,2,4);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=300Hz');xlabel('\omega/hz');ylabel('(H2(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))]);Fs=200;T=1/Fs; M=Tp*Fs;n=0:M-1;t=n*T;A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;xat=A*exp(-alph*t).*sin(omega*t);Xk=T*fft(xat,M); %M点FFT[xat)]subplot(3,2,5); stem(n,xat,'.'); xlabel('n');ylabel('x3(n)'); title('(c) Fs=200Hz');k=0:M-1;fk=k/Tp;subplot(3,2,6);plot(fk,abs(Xk));title('(a) T*FT[xa(nT)],Fs=200Hz');xlabel('\omega/hz');ylabel('(H3(ejw))');axis([0,Fs,0,1.2*max(abs(Xk))])nx 1(n )(a) Fs=1000Hz(a) T*FT[xa(nT)],Fs=1000Hzω/hz(H 1(e j w ))nx 2(n )(b) Fs=300Hz(a) T*FT[xa(nT)],Fs=300Hzω/hz(H 2(e j w ))nx 3(n )(c) Fs=200Hz(a) T*FT[xa(nT)],Fs=200Hzω/hz(H 3(e j w ))%频域采样理论验证程序 clc;clear;close all; M=27;N=32;n=0:M; xa=0:(M/2); xb= ceil(M/2)-1:-1:0; xn=[xa,xb]; %产生M 长三角波序列x(n) Xk=fft(xn,1024); %1024点FFT[x(n)], 用于近似序列x(n)的TF X32k=fft(xn,32) ;%32点FFT[x(n)] x32n=ifft(X32k); %32点IFFT[X32(k)]得到x32(n) X16k=X32k(1:2:N); %隔点抽取X32k 得到X16(K) x16n=ifft(X16k,N/2); %16点IFFT[X16(k)]得到x16(n) subplot(3,2,2);stem(n,xn,'.'); title('(b) 三角波序列x(n)');xlabel('n');ylabel('x(n)');axis([0,32,0,20]) k=0:1023;wk=2*k/1024; subplot(3,2,1);plot(wk,abs(Xk));title('(a)FT[x(n)]'); xlabel('\omega/\pi');ylabel('|X(e^j^\omega)|');axis([0,1,0,200]) k=0:N/2-1; subplot(3,2,3);stem(k,abs(X16k),'.'); title('(c) 16点频域采样');xlabel('k');ylabel('|X_1_6(k)|');axis([0,8,0,200]) n1=0:N/2-1; subplot(3,2,4);stem(n1,x16n,'.'); title('(d) 16点IDFT[X_1_6(k)]');xlabel('n');ylabel('x_1_6(n)');axis([0,32,0,20]); k=0:N-1; subplot(3,2,5);stem(k,abs(X32k),'.'); title('(e) 32点频域采样');xlabel('k');ylabel('|X_3_2(k)|');axis([0,16,0,200]);n1=0:N-1; subplot(3,2,6);stem(n1,x32n,'.');box on title('(f) 32点IDFT[X_3_2(k)]');xlabel('n');ylabel('x_3_2(n)');axis([0,32,0,20])(b) 三角波序列x(n)nx (n )0100200(a)FT[x(n)]ω/π|X (e j ω)|(c) 16点频域采样k|X 16(k)|102030(d) 16点IDFT[X 16(k)]nx 16(n)(e) 32点频域采样k|X 32(k )|(f) 32点IDFT[X 32(k)]nx 32(n )程序清单及波形显示:% 用FFT 对信号作频谱分析 clear all;close all %实验内容(1)=================================================== x1n=[ones(1,4)]; %产生序列向量x1(n)=R4(n) M=8;xa=1:(M/2); xb=(M/2):-1:1; x2n=[xa,xb]; %产生长度为8的三角波序列x2(n) x3n=[xb,xa]; X1k8=fft(x1n,8); %计算x1n 的8点DFT X1k16=fft(x1n,16); %计算x1n 的16点DFT X2k8=fft(x2n,8); %计算x1n 的8点DFT X2k16=fft(x2n,16); %计算x1n 的16点DFT X3k8=fft(x3n,8); %计算x1n 的8点DFT X3k16=fft(x3n,16); %计算x1n 的16点DFT %以下绘制幅频特性曲线 subplot(1,2,1);stem(X1k8,'.'); %绘制8点DFT 的幅频特性图 title('(1a) 8点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度'); subplot(1,2,2);stem(X1k16,'.'); %绘制16点DFT 的幅频特性图 title('(1b)16点DFT[x_1(n)]');xlabel('ω/π');ylabel('幅度');figure(2)subplot(2,2,1);stem(X2k8,'.'); %绘制8点DFT 的幅频特性图 title('(2a) 8点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,2);stem(X2k16,'.'); %绘制16点DFT 的幅频特性图 title('(2b)16点DFT[x_2(n)]');xlabel('ω/π');ylabel('幅度'); subplot(2,2,3);stem(X3k8,'.'); %绘制8点DFT 的幅频特性图 title('(3a) 8点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,4);stem(X3k16,'.'); %绘制16点DFT 的幅频特性图 title('(3b)16点DFT[x_3(n)]');xlabel('ω/π');ylabel('幅度');-2.5-2-1.5-1-0.500.511.522.5(1a) 8点DFT[x 1(n)]ω/π幅度-2.5-2-1.5-1-0.500.511.522.5(1b)16点DFT[x 1(n)]ω/π幅度2468-4-2024(2a) 8点DFT[x 2(n)]ω/π幅度5101520-20-1001020(2b)16点DFT[x 2(n)]ω/π幅度2468-4-2024(3a) 8点DFT[x 3(n)]ω/π幅度-10-50510(3b)16点DFT[x 3(n)]ω/π幅度%实验内容(2) 周期序列谱分析==================================N=8;n=0:N-1; %FFT 的变换区间N=8 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k8=fft(x4n); %计算x4n 的8点DFT X5k8=fft(x5n); %计算x5n 的8点DFT N=16;n=0:N-1; %FFT 的变换区间N=16 x4n=cos(pi*n/4);x5n=cos(pi*n/4)+cos(pi*n/8);X4k16=fft(x4n); %计算x4n 的16点DFT X5k16=fft(x5n); %计算x5n 的16点DFT figure(3)subplot(2,2,1);stem(X4k8,'.'); %绘制8点DFT 的幅频特性图 title('(4a) 8点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,3);stem(X4k16,'.'); %绘制16点DFT 的幅频特性图 title('(4b)16点DFT[x_4(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,2);stem(X5k8,'.'); %绘制8点DFT 的幅频特性图 title('(5a) 8点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');subplot(2,2,4);stem(X5k16,'.'); %绘制16点DFT 的幅频特性图title('(5b)16点DFT[x_5(n)]');xlabel('ω/π');ylabel('幅度');2468-1012x 10-15(4a) 8点DFT[x 4(n)]ω/π幅度5101520-4-2024x 10-15(4b)16点DFT[x 4(n)]ω/π幅度2468-4-2024(5a) 8点DFT[x 5(n)]ω/π幅度5101520-4-2024x 10-15(5b)16点DFT[x 5(n)]ω/π幅度%实验内容(3) 模拟周期信号谱分析=============================== figure(4) Fs=64;T=1/Fs; N=16;n=0:N-1; %FFT 的变换区间N=16 x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)16点采样 X6k16=fft(x6nT); %计算x6nT 的16点DFT X6k16=fftshift(X6k16); %将零频率移到频谱中心 Tp=N*T;F=1/Tp; %频率分辨率F k=-N/2:N/2-1;fk=k*F; %产生16点DFT 对应的采样点频率(以零频率为中心) subplot(3,1,1);stem(fk,abs(X6k16),'.');box on %绘制8点DFT 的幅频特性图 title('(6a) 16点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度'); axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k16))])N=32;n=0:N-1; %FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)32点采样X6k32=fft(x6nT); %计算x6nT的32点DFTX6k32=fftshift(X6k32); %将零频率移到频谱中心Tp=N*T;F=1/Tp; %频率分辨率Fk=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率(以零频率为中心)subplot(3,1,2);stem(fk,abs(X6k32),'.');box on %绘制8点DFT的幅频特性图title('(6b) 32点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k32))])N=64;n=0:N-1; %FFT的变换区间N=16x6nT=cos(8*pi*n*T)+cos(16*pi*n*T)+cos(20*pi*n*T); %对x6(t)64点采样X6k64=fft(x6nT); %计算x6nT的64点DFTX6k64=fftshift(X6k64); %将零频率移到频谱中心Tp=N*T;F=1/Tp; %频率分辨率Fk=-N/2:N/2-1;fk=k*F; %产生16点DFT对应的采样点频率(以零频率为中心)subplot(3,1,3);stem(fk,abs(X6k64),'.'); box on%绘制8点DFT的幅频特性图title('(6a) 64点|DFT[x_6(nT)]|');xlabel('f(Hz)');ylabel('幅度');axis([-N*F/2-1,N*F/2-1,0,1.2*max(abs(X6k64))])-30-20-1001020300510(6a) 16点|DFT[x 6(nT)]|f(Hz)幅度010(6b) 32点|DFT[x 6(nT)]|f(Hz)幅度020(6a) 64点|DFT[x 6(nT)]|f(Hz)幅度程序清单及波形显示: clc;clear all;close all;fc1=250; fm1=15; fc2=500; fm2=50; fc3=1000; fm3=100;N=800; Fs=10000;Ts=1/Fs; n=[0:N-1];t=n*Ts;x11=cos(2*pi*fc1*t); x12=cos(2*pi*fm1*t); x1=x11.*x12; subplot(3,1,1);plot(t,x11,'g');plot(t,x12,'r');plot(t,x1,'b'); legend('载波','调制波 ','已调 ');xlabel('t/s');ylabel('波形')x=cos(2*pi*fc1*t).*cos(2*pi*fm1*t)+cos(2*pi*fc2*t).*cos(2*pi*fm2*t)+cos(2*pi*fc3*t).*cos(2*pi*fm3*t); subplot(3,1,2);plot(t,x);X=fft(x)subplot(3,1,3)k=[0:(N-1)/2]stem(k*2/N,abs(X(k+1))/max(abs(X(k+1))),'.');axis([0,0.3,0,1]);xlabel('\omeg a/\pi');ylabel('幅度');wp=[0.04,0.06];ws=[0.03,0.07];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形')subplot(3,1,2);[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-80,1]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为250Hz的频率响应');subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为250H的滤波信号')wp=[0.08,0.12];ws=[0.07,0.13];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形');subplot(3,1,2)[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-90,2]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为500Hz的频率响应')subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为500H的滤波信号')wp=[0.17,0.23];ws=[0.16,0.24];rp=0.1;rs=60;[N1,wp]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N1,rp,rs,wp);y1=filter(B,A,x);figuresubplot(3,1,1);plot(t,x);xlabel('t/s');title('3路混合信号波形');subplot(3,1,2);[H1,w1]=freqz(B,A,N);plot(w1/pi,20*log10(abs(H1)));axis([0,0.5,-100,10]);xlabel('\omega/\pi');ylabel('|H(e^j\omega)|');title('中心频率为1000Hz的频率响应')subplot(3,1,3);plot(t,y1);xlabel('t/s');ylabel('y1(t)');title('中心频率为1000H的滤波信号')00.010.020.030.040.050.060.070.08-101t/s波形00.010.020.030.040.050.060.070.08-505ω/π幅度-505t/s3路混合信号波形-80-60-40-200ω/π|H (e j ω)|中心频率为250Hz 的频率响应-202t/sy 1(t )中心频率为250H 的滤波信号0.010.020.030.040.050.060.070.08-505t/s3路混合信号波形0.050.10.150.20.250.30.350.40.450.5-80-60-40-200ω/π|H (e j ω)|中心频率为500Hz 的频率响应0.010.020.030.040.050.060.070.08-202t/sy 1(t )中心频率为500H 的滤波信号0.010.020.030.040.050.060.070.08-505t/s3路混合信号波形0.050.10.150.20.250.30.350.40.450.5-100-500ω/π|H (e j ω)|中心频率为1000Hz 的频率响应0.010.020.030.040.050.060.070.08-202t/sy 1(t )中心频率为1000H 的滤波信号程序清单及波形显示: clc;clear;clear allN=1000; Fs=1000;T=1/Fs;Tp=N*T; t=0:T:(N-1)*T;fc=Fs/10;f0=fc/10; %载波频率fc=Fs/10,单频调制信号频率为f0=Fc/10; mt=cos(2*pi*f0*t); %产生单频正弦波调制信号mt ,频率为f0 ct=cos(2*pi*fc*t); %产生载波正弦波信号ct ,频率为fc xt=mt.*ct; %相乘产生单频调制信号xt nt=2*rand(1,N)-1; %产生随机噪声ntfp=150; fs=200;Rp=0.1;As=70; % 滤波器指标fb=[fp,fs];m=[0,1]; % 计算remezord 函数所需参数f,m,devdev=[10^(-As/20),(10^(Rp/20)-1)/(10^(Rp/20)+1)];[n,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez 函数所需参数hn=remez(n,fo,mo,W); % 调用remez 函数进行设计,用于滤除噪声nt 中的低频成分yt=filter(hn,1,10*nt); %滤除随机噪声中低频成分,生成高通噪声yt xt=xt+yt; %噪声加信号 fst=fft(xt,N);k=0:N-1;f=k/Tp;subplot(2,1,1);plot(t,xt);grid;xlabel('t/s');ylabel('x(t)');axis([0,Tp/5,min(xt),max(xt)]);title('(a) 信号加噪声波形')subplot(2,1,2);plot(f,abs(fst)/max(abs(fst)));grid;title('(b) 信号加噪声的频谱')axis([0,Fs/2,0,1.2]);xlabel('f/Hz');ylabel('幅度')-10-50510t/sx (t )(a) 信号加噪声波形0501001502002503003504004505000.51(b) 信号加噪声的频谱f/Hz幅度%==调用xtg 产生信号xt, xt 长度N=1000,并显示xt 及其频谱,========= N=1000;xt=xtg(N);fp=120; fs=150;Rp=0.2;As=60;Fs=1000; % 输入给定指标% (1) 用窗函数法设计滤波器wc=(fp+fs)/Fs; %理想低通滤波器截止频率(关于pi归一化)B=2*pi*(fs-fp)/Fs; %过渡带宽度指标Nb=ceil(11*pi/B); %blackman窗的长度Nhn=fir1(Nb-1,wc,blackman(Nb));Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性ywt=fftfilt(hn,xt,N); %调用函数fftfilt对xt滤波subplot(4,1,1);plot(f,Hw);xlabel('f/Hz');ylabel('幅度'); title('(a)低通滤波器幅频特性');subplot(4,1,2);plot(t,ywt); title('(b)滤除噪声后的信号波形');xlabel('t/s');ylabel('ywt');% (2) 用等波纹最佳逼近法设计滤波器fb=[fp,fs];m=[1,0]; % 确定remezord函数所需参数f,m,dev dev=[(10^(Rp/20)-1)/(10^(Rp/20)+1),10^(-As/20)];[Ne,fo,mo,W]=remezord(fb,m,dev,Fs); % 确定remez函数所需参数hn=remez(Ne,fo,mo,W); % 调用remez函数进行设计Hw=abs(fft(hn,1024)); % 求设计的滤波器频率特性yet=fftfilt(hn,xt,N); % 调用函数fftfilt对xt滤波subplot(3,1,1);plot(f,Hw);xlabel('f/Hz');ylabel('幅度'); title('(c)低通滤波器幅频特性');subplot(3,1,2);plot(t,yet);title('(d)滤除噪声后的信号波形');xlabel('t/s');ylabel('yet');。
数字信号处理课后答案西安电子
第2章 时域离散信号和系统的频域分析 解: 假设输入信号x(n)=ejω0n,系统单位脉冲响应为h(n) 则系统,输出为
上式说明当输入信号为复指数序列时, 输出序列仍是复 指数序列, 且频率相同, 但幅度和相位取决于网络传 输函数。 利用该性质解此题:
第2章 时域离散信号和系统的频域分析
第2章 时域离散信号和系统的频域分析
上式中|H(ejω)|是ω的偶函数, 相位函数是ω的奇函数
|H(ej,ω)|=|H(e-
θ(ω)=-θ(-ω), 故
jω)|,
4. 设
第2章 时域离散信号和系统的频域分析
将x(n)以4为周期进行周期延拓, 形成周期序列
,
画出x(n)和
的波形, 求出
的离散傅里叶级数
和傅里叶变换。
解: 画出x(n) 和
由z3(z-1)=0, 得极点为 z1, 2=0, 1 零极点图和收敛域如题15解图(a)所示, 点相互对消。
图中, z=1处的零极
第2章 时域离散信号和系统的频域分析 题15解图
第2章 时域离散信号和系统的频域分析 (2)
第2章 时域离散信号和系统ຫໍສະໝຸດ 频域分析零点为极点 为
极零点分布图如题15解图(b)所示。 (3) 令y(n)=R4(n), 则
(4) δ(n)
(6) 2-n[u(n)-u(n-10)]
第2章 时域离散信号和系统的频域分析 解 (1)
(2)
第2章 时域离散信号和系统的频域分析 (3)
(4) ZT[δ(n)]=10≤|z|≤∞ (5) ZT[δ(n-1)]=z-10<|z|≤∞
(6)
≤
第2章 时域离散信号和系统的频域分析
16. 已 知
数字信号处理-西安电子科技大学出版(_高西全丁美玉)第三版_课后习题答案(全)
18
第 1 章 时域离散信号和时域离散系统
(3) 这是一个延时器, 延时器是线性非时变系统, 下面证明。 令输入为
输出为
x(n-n1)
y′(n)=x(n-n1-n0) y(n-n1)=x(n-n1-n0)=y′(n) 故延时器是非时变系统。 由于
T[ax1(n)+bx2(n)]=ax1(n-n0)+bx2(n-n0) =aT[x1(n)]+bT[x2(n)]
x(m)h(n-m)
m
第 1 章 时域离散信号和时域离散系统
题7图
28
第 1 章 时域离散信号和时域离散系统
y(n)={-2,-1,-0.5, 2, 1, 4.5, 2, 1; n=-2, -1, 0, 1, 2, 3, 4, 5}
第 1 章 时域离散信号和时域离散系统
解法(二) 采用解析法。 按照题7图写出x(n)和h(n)的表达式分别为
5. 设系统分别用下面的差分方程描述, x(n)与y(n)分别表示系统输入和输 出, 判断系统是否是线性非时变的。
(1)y(n)=x(n)+2x(n-1)+3x(n-2) (2)y(n)=2x(n)+3 (3)y(n)=x(n-n0) n0 (4)y(n)=x(-n)
15
第 1 章 时域离散信号和时域离散系统
非零区间如下:
0≤m≤3 -4≤m≤n
第 1 章 时域离散信号和时域离散系统
根据非零区间, 将n分成四种情况求解: ① n<0时, y(n)=0
② 0≤n≤3时, y(n)= ③ 4≤n≤7时, y(n)= ④ n>7时, y(n)=0
1=n+1
n
1=8-m n0
西电数字信号处理大作业-浅谈奈奎斯特频率采样和压缩感知概要
浅谈奈奎斯特频率采样和压缩感知信息技术的飞速发展使得人们对信息的需求量剧增。
现实世界的模拟化和信号处理工具的数字化决定了信号采样是从模拟信源获取数字信息的必经之路。
在信号和图像处理领域,凡是涉及到计算机作为处理工具的场合,所面临的首要问题就是模拟信号的数字化问题,然后再对得到的离散的样本进行各种处理。
连续信号转化为离散的数字化信号的过程称为采样。
对模拟信号采样所得的离散数字信号能否代表并恢复成原来的连续模拟信号呢?如能恢复应具备什么样的条件呢?这个问题直接关系到是否可以用数字处理工具和数字化的方法处理模拟信号。
一奈奎斯特频率采样奈奎斯特采样定理给我们提供了如何采样的重要理论基础。
它指出,如果信号是带限的,采样速率必须达到信号带宽的两倍以上才能精确重构信号。
事实上,在音频和可视电子设备、医学图像设备、无线接收设备等设备中的所有信号采样协议都隐含了这样的限制。
奈奎斯特采样定理至出现以来一直是数字信号和图像处理领域的重要理论基础,它支撑着几乎所有的信号和图像处理过程,包括信号和图像的获取、存储、处理、传输等。
采样定理,又称香农采样定理,奈奎斯特采样定理,是信息论,特别是通讯与信号处理学科中的一个重要基本结论.E.T.Whittaker (1915年发表的统计理论),克劳德·香农与Harry Nyquist都对它作出了重要贡献。
另外,V. A. Kotelnikov 也对这个定理做了重要贡献。
采样是将一个信号(即时间或空间上的连续函数)转换成一个数值序列(即时间或空间上的离散函数)。
采样定理指出,如果信号是带限的,并且采样频率高于信号带宽的一倍,那么,原来的连续信号可以从采样样本中完全重建出来。
带限信号变换的快慢受到它的最高频率分量的限制,也就是说它的离散时刻采样表现信号细节的能力是有限的。
采样定理是指,如果信号带宽小于奈奎斯特频率(即采样频率的二分之一),那么此时这些离散的采样点能够完全表示原信号。
西电数字信号处理上机实验报告
数字信号处理上机实验报告14020710021 张吉凯第一次上机实验一:设给定模拟信号()1000t a x t e -=,t 的单位是ms 。
(1) 利用MATLAB 绘制出其时域波形和频谱图(傅里叶变换),估计其等效带宽(忽略谱分量降低到峰值的3%以下的频谱)。
(2) 用两个不同的采样频率对给定的()a x t 进行采样。
○1()()15000s a f x t x n =以样本秒采样得到。
()()11j x n X e ω画出及其频谱。
○2()()11000s a f x t x n =以样本秒采样得到。
()()11j x n X e ω画出及其频谱。
比较两种采样率下的信号频谱,并解释。
(1)MATLAB 程序:N=10; Fs=5; T s=1/Fs;n=[-N:T s:N];xn=exp(-abs(n)); w=-4*pi:0.01:4*pi; X=xn*exp(-j*(n'*w)); subplot(211) plot(n,xn);title('x_a(t)时域波形');xlabel('t/ms');ylabel('x_a(t)'); axis([-10, 10, 0, 1]); subplot(212);plot(w/pi,abs(X)); title('x_a(t)频谱图');xlabel('\omega/\pi');ylabel('X_a(e^(j\omega))');ind = find(X >=0.03*max(X))*0.01; eband = (max(ind) -min(ind)); fprintf('等效带宽为%fKHZ\n',eband); 运行结果:等效带宽为12.110000KHZ(2)MATLAB程序:N=10;omega=-3*pi:0.01:3*pi;%Fs=5000Fs=5;T s=1/Fs;n=-N:T s:N;xn=exp(-abs(n));X=xn*exp(-j*(n'*omega));subplot(2,2,1);stem(n,xn);grid on;axis([-10, 10, 0, 1.25]); title('时域波形(f_s=5000)');xlabel('n');ylabel('x_1(n)');subplot(2,2,2);plot(omega/pi,abs(X));title('频谱图(f_s=5000)');xlabel('\omega/\pi');ylabel('X_1(f)');grid on;%Fs=1000Fs=1;T s=1/Fs;n=-N:T s:N;xn=exp(-abs(n));X=xn*exp(-j*(n'*omega));subplot(2,2,3);stem(n,xn);grid on;axis([-10, 10, 0, 1.25]); title('时域波形(f_s=1000)');xlabel('n');ylabel('x_2(n)');grid on;subplot(2,2,4);plot(omega/pi,abs(X));title('频谱图(f_s=1000)');xlabel('\omega/\pi');ylabel('X_2(f)');grid on;运行结果:实验二:给定一指数型衰减信号()()0cos 2at x t e f t π-=,采样率1s f T=,T 为采样周期。
西电DSP实验报告
西安电子科技大学数字信号处理(DSP)课程实验报告实验名称 VISUAL DSP++的使用入门电子工程学院 1402071 班同作者崔健孟智超杨伟祺实验日期 2017 年 5 月 16 日实验一:VISUAL DSP++的使用入门一、实验目的:熟悉VISUAL DSP++的开发环境二、实验内容:练习一:启动Visual DSP++,建立一个用C源代码的工程(Project),同时用调试器来评估用C语言所编写代码的性能;练习二:创立一个新的工程,修改源码来调用一个汇编(asm)程序,重新编译工程,用调试器来评估用汇编语言所写程序的性能;练习三:利用调试器的绘图(plot)功能来图形显示一个卷积算法中的多个数据的波形;练习四:利用调试器的性能统计功能(Statistical profile)来检查练习三中卷积算法的效率。
利用所收集到的性能统计数据就能看出算法中最耗时的地方。
三、实验步骤及实验结果:练习一:1)新建工程进入 Visual DSP++,显示Visual DSP++的集成开发和调试环境窗口,选择菜单File 中Open 打开文件:…\unit_1\dot_product_c \dotprodc.dpj。
2)编译 dotprodc工程在菜单 Project中选择 Build Project来对工程进行编译。
在本例子中,编译器会检测到一个未定义的错误,显示为:“.\dotprod_main.c”,line 115:error #20:identifier“itn”is undefined itn i;将该错误改正后,保存并重新编译,没有错误出现,编译成功。
3)运行VsualDSP++调试器在编译完成后,环境将自动进入调试状态,对于初次进入debugger,将显示对象选择对话框,在其中指定对象和处理器信息。
4)运行dotprod.c从 Debug菜单中选择 Run项,程序将被执行,其输出结果在 Output window中显示。
数字信号处理大作业语音信号提取程序加噪声程序
% 数字信号处理大作业语音信号提取程序加噪声程序
% 功能:从文件名为''语音文件中提取出数据,并加入声音最大音量五分之一的噪声% 参数:name-原始信号;Y-加噪声后的信号;
% 大作业内容:请自己录音,编写FIR滤波程序,消减噪声的影响。
% 作者:西安电子科技大学雷达信号处理国家重点实验室李军
% 时间:2015年
close all;clear all;clc;
namex=1e5;
% 读取原始声音文件
[name,fs0,bit]=wavread('',[100,namex]); %读取源文件
sound(name,fs0,bit);
namex=length(name);
noise=(max(name(:,1))/5)*randn(namex,2);
Y=name+noise; % 加入噪声
wavwrite(Y,fs0,bit,'mingyun_'); % 写入噪声文件
sound(Y,fs0,bit);
% 作图分析声音信号的时域特性和频域特性
figure(1)
subplot(2,1,1);
plot(abs(name));
title('(a)原声音信号的时域特性');
subplot(2,1,2);
plot(abs(Y));
title('(b)含有噪声声音信号的时域特性');。
西电-数字信号处理大作业
数字信号处理上机大作业实验一:信号、系统及系统响应(1) 简述实验目的及实验原理。
1.实验目的●熟悉连续信号经理想采样前后的频谱变化关系,加深对时域采样定理的理解。
●熟悉时域离散系统的时域特性。
●利用卷积方法观察分析系统的时域特性。
●掌握序列傅里叶变换的计算机实现方法,利用序列的傅里叶变换对连续信号、离散信号及系统响应进行频域分析。
2.实验原理与方法●时域采样。
● LTI系统的输入输出关系。
(2)按实验步骤附上实验过程中的信号序列、系统单位脉冲响应及系统响应序列的时域和幅频特性曲线,并对所得结果进行分析和解释。
Matlab源程序如下:A=1;T1=1/1000;T2=1/300;T3=1/200;a=25*pi;w0=30*pi;n=0:99;x1=A*exp(-a*n*T1).*sin(w0*n*T1);x2=A*exp(-a*n*T2).*sin(w0*n*T2);x3=A*exp(-a*n*T3).*sin(w0*n*T3);m=linspace(-pi,pi,10000);X1=x1*exp(-j*n'*m);%n'与m构造矩阵,xi向量与矩阵每一列相乘对应元素相加,构成DTFT后的矩阵X2=x2*exp(-j*n'*m);X3=x3*exp(-j*n'*m);figure(1);subplot(3,2,1)plot(m/pi,abs(X1));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为1000Hz时的幅度谱');subplot(3,2,3)plot(m/pi,abs(X2));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为300Hz时的幅度谱');subplot(3,2,5)plot(m/pi,abs(X3));xlabel('\omega/π');ylabel('|H(e^j^\omega)|');title('采样频率为200Hz时的幅度谱');subplot(3,2,2)plot(n,abs(x1));xlabel('n');ylabel('x1(t)');title('采样频率为1000Hz时的时域波形');subplot(3,2,4)plot(n,abs(x2));xlabel('n');ylabel('x2(t)');title('采样频率为300Hz时的时域波形');subplot(3,2,6)plot(n,abs(x3));xlabel('n');ylabel('x3(t)');title('采样频率为200Hz时的时域波形');波形图如下:-1-0.8-0.6-0.4-0.200.20.40.60.81ω/π|H (e j ω)|采样频率为1000Hz 时的幅度谱ω/π|H (e j ω)|采样频率为300Hz 时的幅度谱ω/π|H (e j ω)|采样频率为200Hz 时的幅度谱102030405060708090100nx 1(t )采样频率为1000Hz 时的时域波形nx 2(t )采样频率为300Hz 时的时域波形nx 3(t )采样频率为200Hz 时的时域波形② 时域离散信号、 系统和系统响应分析。
西电数字信号处理大作业
实验一、信号的采样clc,clear;dt=0.001;tf=6;t=0:dt:tf;xa=sqrt(t)+cos(t);T=0.5;n=0:tf/T;x=sqrt(n*T)+cos(n*T);figure(1)subplot(2,1,1)plot(t,xa),grid on ;title('original image')subplot(2,1,2)stem(n*T,x),grid on ,title('digital image')实验二、信号与系统的时域分析差分方程为)()2()1()(21n bx n y a n y a n y +----=,其中8.01-=a ,64.02=a ,866.0=b 。
系统单位脉冲响应)(n ha1=-0.8;a2=0.64;b=0.866;ys=0;xn=[1,zeros(1,49)];B=1;A=[1,a1,a2];xi=filtic(B,A,ys);yn=filter(B,A,xn,xi);n=0:length(yn)-1;subplot(1,1,1);stem(n,yn,'.')title('(a)');xlabel('n');ylabel('y(n)')输入x(n)=cos(n)T=0.1;z=cos(n*T);zn=conv(yn,z); figure(2);n1=1:99;stem(n1,zn,'.')实验三、系统的频域和Z域分析程序代码(画出dtft的幅度和频率谱)clc,clear;n=0:1:7;x=(0.9*exp(j*pi/3)).^n;w=0:pi/200:pi;X=x*exp(-j).^(n'*w);realX=real(X);imagX=imag(X);angX=angle(X);magX=abs(X);subplot(2,2,1);plot(w/pi,magX);grid xlabel('frequency in pi unit');title('magnitude part');subplot(2,2,2);plot(w/pi,realX);grid xlabel('frequency in pi unit');title('real part');subplot(2,2,3);plot(w/pi,imagX);grid xlabel('frequency in pi unit');title('imaginary part');subplot(2,2,4);plot(w/pi,angX);grid xlabel('frequency in pi unit');title('angel part');clc,clear;a=[1,-0.5,0.06];b=[1,1,0];m=0:length(b)-1;l=0:length(a)-1;w=0:pi/500:pi;num=b*exp(-j*m'*w);den=a*exp(-j*l'*w);H=num./den;magH=abs(H);angH=angle(H);H1=freqz(b,a,w);magH1=abs(H1);angH1=angle(H1);subplot(2,2,2);plot(w/pi,angH);grid;xlabel('w(frequency in pi units)');ylabel('Ïàλrad/w');subplot(2,2,1);plot(w/pi,magH);grid;xlabel('w(frequency in pi units)');ylabel('·ù¶È|H|');subplot(2,2,3);plot(w/pi,magH1);grid;xlabel('w(frequency in pi units)');ylabel('·ù¶È|H1|');subplot(2,2,4);plot(w/pi,angH);grid;xlabel('w(frequency in pi units)');ylabel('Ïàλrad/w');axis([0,1,-0.8,0]); figure(2);zplane(b,a);实验四、信号的频谱分析程序代码clc,clear;n=0:7;k=0:7;N=8;w=n*(2*pi)/8;x=(0.9*exp(j*pi/3)).^n;X1=[x zeros(1,8)];X2=[X1 zeros(1,16)];XK=x*exp(-j*k'*w);k1=0:15;n1=0:15;w1=n1*(2*pi)/16;XK1=X1*exp(-j*k1'*w1);k2=0:31;n2=0:31;w2=n2*(2*pi)/16;XK2=X2*exp(-j*k2'*w2);w3=0:pi/200:2*pi;X=x*exp(-j*n'*w3);magX=abs(X);angX=angle(X);magXK=abs(XK);angXK=angle(XK);magXK1=abs(XK1);angXK1=angle(XK1);magXK2=abs(XK2);angXK2=angle(XK2);subplot(4,2,1);plot(w3/pi,magX);xlabel('w/pi');ylabel('·ù¶È|X|');grid on;subplot(4,2,2);plot(w3/pi,angX);xlabel('w/pi');ylabel('Ïàλrad/pi'); subplot(4,2,3);stem(n,magXK);xlabel('K');ylabel('·ù¶È|XK|');subplot(4,2,4);stem(n,magXK);xlabel('K');ylabel('Ïàλrad/pi'); subplot(4,2,5);stem(n1,magXK1);xlabel('K1');ylabel('·ù¶È|XK1|'); subplot(4,2,6);stem(n1,magXK1);xlabel('K1');ylabel('Ïàλrad/pi'); subplot(4,2,7);stem(n2,magXK2);xlabel('K2');ylabel('·ù¶È|XK2|'); subplot(4,2,8);stem(n2,magXK2);xlabel('K2');ylabel('Ïàλrad/pi');实验五、IIR数字滤波器设计IIR汉宁窗低通高通低通巴特沃斯通带截止频率wp=0.2pi 通带最大衰减R=1dB阻带截止频率wp=0.35pi 阻带最大衰减R=10dBclc,clear;Wp=0.2;Ws=0.35;Rp=1;Rs=100;[N,Wc]=buttord(Wp,Ws,Rp,Rs);[Bz,Az]=butter(N,Wc)w=0:0.1:pi;[H,w1]=freqz(Bz,Az,w);;ang=angle(H);H=20*log10(abs(H))subplot(4,2,1); plot(w/pi,H) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('µÍͨÂ˲¨Æ÷')subplot(4,2,2);plot(w/pi,ang);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')[Bz1,Az1]=butter(N,Wc,'high')w=0:0.1:pi;[H1,w2]=freqz(Bz1,Az1,w);ang1=angle(H1);H1=20*log10(abs(H1))subplot(4,2,3); plot(w/pi,H1) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('¸ßͨÂ˲¨Æ÷')subplot(4,2,4);plot(w/pi,ang1);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')Wp1=[0.2 0.8];Ws1=[0.35 0.65];[N2,Wc1]=buttord(Wp1,Ws1,Rp,Rs);[Bz2,Az2]=butter(N2,Wc1,'stop')w=0:0.1:pi;[H2,w3]=freqz(Bz2,Az2,w);ang2=angle(H2);H2=20*log10(abs(H2))subplot(4,2,5); plot(w/pi,H2) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('´ø×èÂ˲¨Æ÷')subplot(4,2,6);plot(w/pi,ang2);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')Wp1=[0.2 0.8];Ws1=[0.35 0.65];[N2,Wc1]=buttord(Wp1,Ws1,Rp,Rs);[Bz3,Az3]=butter(N2,Wc1)w=0:0.1:pi;[H3,w4]=freqz(Bz3,Az3,w);ang3=angle(H3);H3=20*log10(abs(H3))subplot(4,2,7); plot(w/pi,H3) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('´øͨÂ˲¨Æ÷')subplot(4,2,8);plot(w/pi,ang3);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')切比雪夫1型通带截止频率wp=0.7pi 通带最大衰减R=1dB阻带截止频率wp=0.5pi 阻带最大衰减R=40dBclc,clear;Wp=0.7;Ws=0.5;Rp=1;Rs=40;[N,Wpo]=cheb1ord(Wp,Ws,Rp,Rs);[Bz,Az]=cheby1(N,Rp,Wpo)w=0:0.1:pi;[H,w1]=freqz(Bz,Az,w);ang=angle(H);H=20*log10(abs(H))subplot(4,2,1); plot(w/pi,H) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('µÍͨÂ˲¨Æ÷')subplot(4,2,2);plot(w/pi,ang);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')[Bz1,Az1]=cheby1(N,Rp,Wpo,'high');w=0:0.1:pi;[H1,w2]=freqz(Bz1,Az1,w);ang1=angle(H1);H1=20*log10(abs(H1))subplot(4,2,3); plot(w/pi,H1) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('¸ßͨÂ˲¨Æ÷')subplot(4,2,4);plot(w/pi,ang1);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')Wp1=[0.2 0.8];Ws1=[0.35 0.65];[N2,Wpo1]=cheb1ord(Wp1,Ws1,Rp,Rs);[Bz2,Az2]=cheby1(N2,Rp,Wpo1,'stop')w=0:0.1:pi;[H2,w3]=freqz(Bz2,Az2,w);ang2=angle(H2);H2=20*log10(abs(H2))subplot(4,2,5); plot(w/pi,H2) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('´ø×èÂ˲¨Æ÷')subplot(4,2,6);plot(w/pi,ang2);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')Wp1=[0.2 0.8];Ws1=[0.35 0.65];[N2,Wpo1]=cheb1ord(Wp1,Ws1,Rp,Rs);[Bz3,Az3]=cheby1(N2,Rp,Wpo1)w=0:0.1:pi;[H3,w4]=freqz(Bz3,Az3,w);ang3=angle(H3);H3=20*log10(abs(H3))subplot(4,2,7); plot(w/pi,H3) ;gridon ;xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|/dB');title('´øͨÂ˲¨Æ÷')subplot(4,2,8);plot(w/pi,ang3);gridon ;xlabel('\omega/\pi');ylabel('Phase/dB')实验六、FIR数字滤波器设计FIR汉宁窗低通高通低通% 采用Hamming窗设计一个带阻FIR滤波器阻带:0~0.5pi,阻带最小衰减Rs=40dB;通带:0.5~pi,通带最大衰减:Rp=1dB。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
数字信号处理大作业班级:021231学号:姓名:指导老师:吕雁一写出奈奎斯特采样率和和信号稀疏采样的学习报告和体会1、采样定理在进行A/D信号的转换过程中,当采样频率fs.max大于信号中最高频率fmax的2倍时(fs.max>2fmax),采样之后的数字信号完整地保留了原始信号中的信息,一般实际应用中保证采样频率为信号最高频率的5~10倍;采样定理又称奈奎斯特定理。
(1)在时域频带为F的连续信号 f(t)可用一系列离散的采样值f(t1),f(t1±Δt),f(t1±2Δt),...来表示,只要这些采样点的时间间隔Δt≤1/2F,便可根据各采样值完全恢复原始信号。
(2)在频域当时间信号函数f(t)的最高频率分量为fmax时,f(t)的值可由一系列采样间隔小于或等于1/2fo的采样值来确定,即采样点的重复频率fs ≥2fmax。
2、奈奎斯特采样频率(1)概述奈奎斯特采样定理:要使连续信号采样后能够不失真还原,采样频率必须大于信号最高频率的两倍(即奈奎斯特频率)。
奈奎斯特频率(Nyquist frequency)是离散信号系统采样频率的一半,因哈里·奈奎斯特(Harry Nyquist)或奈奎斯特-香农采样定理得名。
采样定理指出,只要离散系统的奈奎斯特频率高于被采样信号的最高频率或带宽,就可以真实的还原被测信号。
反之,会因为频谱混叠而不能真实还原被测信号。
采样定理指出,只要离散系统的奈奎斯特频率高于采样信号的最高频率或带宽,就可以避免混叠现象。
从理论上说,即使奈奎斯特频率恰好大于信号带宽,也足以通过信号的采样重建原信号。
但是,重建信号的过程需要以一个低通滤波器或者带通滤波器将在奈奎斯特频率之上的高频分量全部滤除,同时还要保证原信号中频率在奈奎斯特频率以下的分量不发生畸变,而这是不可能实现的。
在实际应用中,为了保证抗混叠滤波器的性能,接近奈奎斯特频率的分量在采样和信号重建的过程中可能会发生畸变。
因此信号带宽通常会略小于奈奎斯特频率,具体的情况要看所使用的滤波器的性能。
需要注意的是,奈奎斯特频率必须严格大于信号包含的最高频率。
如果信号中包含的最高频率恰好为奈奎斯特频率,那么在这个频率分量上的采样会因为相位模糊而有无穷多种该频率的正弦波对应于离散采样,因此不足以重建为原来的连续时间信号。
(2)奈奎斯特频率的应用除了奈奎斯特频率之外,还有一个指标非常重要,这个指标就是测量装置的带宽。
严格讲,带宽包含上限和下限两个数值,但是,由于许多宽频带的测量设备,比如说变频功率分析仪,其带宽的频率上限远远大于频率下限,或者频率下限为零,因此,一般以频率上限作为该仪器的带宽。
一般而言,带宽指-3db带宽。
-3db带宽并不表明高于带宽上限频率的信号不能通过测量仪器。
举例而言,某功率分析仪的带宽上限为100kHz,那么,100kHz的正弦波通过测量仪器的AD转换器之前的电路时,幅值衰减为原信号幅值的70.7%,功率衰减为原信号的50%。
此外,对于非正弦波形,其含有的谐波频率高于信号频率(基波频率)。
因此,不能简单的认为,100kHz带宽的仪器可以用于测量100kHz的正弦波,更不能认为100kHz带宽的仪器可以用于测量100kHz的方波或畸变波形。
要让采样过程符合奈奎斯特采样定理,测量仪器的带宽应该小于奈奎斯特频率。
若测量仪器的电路固有带宽高于奈奎斯特频率,应该在AD转换器之间加上截至频率小于奈奎斯特频率的防混叠滤波器。
对于后者,防混叠滤波器的截至频率就是仪器的带宽。
3、稀疏采样目前,Candes,Romberg,Tao和Donoho等人提出了一种全新的理论一压缩感知理论(Compressed Sensing)。
该理论是一种崭新的信号采样、信号编码和信号解码理论。
采样速率不再像Nyquist速率一样,与信号的带宽密切相关,而是与信息在信号中的结构和位置息息相关。
编码过程是围绕观测器即观测矩阵展开的,而解码过程是一个优化计算过程。
该理论已经被证明能够用较低采样速率准确的进行信号采样,并且能够以很高的概率重构原始信号。
目前国内已经有科研单位的学者对其展开研究。
如我们学校课题组基于该理论提出采用超低速率采样检测超宽带回波信号。
其CS理论如图:稀疏采样,也被称为压缩感知、压缩传感或压缩采样,是一种利用稀疏的或可压缩的信号进行信号重构的技术。
或者可以说是信号在采样的同时被压缩,从而在很大程度上降低了采样率。
稀疏采样跳过了采集N个样本这一步骤,直接获得压缩的信号的表示。
其理论利用到了许多自然信号在特定的基上具有紧凑的表示。
即这些信号是“稀疏”的或“可压缩”的。
由于这一特性,稀疏采样理论的信号编解码框架和传统的压缩过程大不一样,主要包括信号的稀疏表示、编码测量和重构算法等三个方面。
简单地说,压缩感知理论指出:只要信号是可压缩的或在某个变换域是稀疏的,那么就可以用一个与变换基不相关的观测矩阵将变换所得高维信号投影到一个低维空间上,然后通过求解一个优化问题就可以从这些少量的投影中以高概率重构出原信号,可以证明这样的投影包含了重构信号的足够信息。
在该理论框架下,采样速率不再取决于信号的带宽,而在很大程度上取决于两个基本准则:稀疏性和非相干性,或者稀疏性和等距约束性。
显然,在压缩感知理论中,图像/信号的采样和压缩同时以低速率进行,使传感器的采样和计算成本大大降低,而信号的恢复过程是一个优化计算的过程.因此,该理论指出了将模拟信号直接采样压缩为数字形式的有效途径。
从理论上讲任何信号都具有可压缩性,只要能找到其相应的稀疏表示空间,就可以有效地进行压缩采样。
当前,压缩感知理论主要涉及三个核心问题:(1) 具有稀疏表示能力的过完备字典设计;(2) 满足非相干性或等距约束性准则的测量矩阵设计;(3) 快速鲁棒的信号重建算法设计。
压缩感知理论必将给信号采样方法带来一次新的革命。
这一理论的引人之处还在于它对应用科学的许多领域具有重要的影响,如统计学、信息论、编码等。
目前,学者们已经在模拟-信息采样、合成孔径雷达成像、遥感成像、核磁共振成像、深空探测成像、无线传感器网络、信源编码、人脸识别、语音识别、探地雷达成像等诸多领域对压缩感知展开了广泛的应用研究。
Rice大学已经成功设计出了一种基于压缩感知的新型单像素相机,在实践中为取代传统相机迈出了实质性的一步。
(1)压缩感知理论框架传统的信号采集、编解码过程如图所示:编码端先对信号进行采样,再对所有采样值进行变换,并将其中重要系数的幅度和位置进行编码,最后将编码值进行存储或传输:信号的解码过程仅仅是编码的逆过程,接收的信号经解压缩、反变换后得到恢复信号。
采用这种传统的编解码方法,由于信号的采样速率不得低于信号带宽的2倍,使得硬件系统面临着很大的采样速率的压力。
此外在压缩编码过程中,大量变换计算得到的小系数被丢弃,造成了数据计算和内存资源的浪费。
传统编解码理论的框图压缩感知理论对信号的采样、压缩编码发生在同一个步骤,利用信号的稀疏性,以远低于Nyquist采样率的速率对信号进行非自适应的测量编码。
测量值并非信号本身,而是从高维到低维的投影值,从数学角度看,每个测量值是传统理论下的每个样本信号的组合函数,即一个测量值已经包含了所有样本信号的少量信息。
解码过程不是编码的简单逆过程,而是在盲源分离中的求逆思想下。
利用信号稀疏分解中已有的重构方法在概率意义上实现信号的精确重构或者一定误差下的近似重构。
解码所需测量值的数目远小于传统理论下的样本数。
压缩感知理论的编解码框图(2) 压缩感知的基本理论及核心问题 假设有一信号)(N R f f ∈,长度为N ,基向量为),...,2,1(N i i =ψ,对信号进行变换:αψψ==∑=f a f i Ni i 或1 显然f 是信号在时域的表示,α是信号在ψ域的表示。
信号是否具有稀疏性或者近似稀疏性是运用压缩感知理论的关键问题,若式中的α只有K 个是非零值)(K N >>者仅经排序后按指数级衰减并趋近于零,可认为信号是稀疏的。
信号的可稀疏表示是压缩感知的先验条件。
在已知信号是可压缩的前提下,压缩感知过程可分为两步:设计一个与变换基不相关的)(N M N M <<⨯维测量矩阵对信号进行观测,得到M 维的测量向量。
由M 维的测量向量重构信号。
(3) 信号的稀疏表示稀疏的数学定义:信号X 在正交基ψ下的变换系数向量为X T ψ=Θ,假如对于20<<p 和0>R ,这些系数满足:Rp p i i p ≤≡Θ∑/1)||(||||θ则说明系数向量Θ在某种意义下是稀疏的。
给出另一种定义:如果变换系数 >ψ=<i i X ,θ的支撑域}0;{≠i i θ的势小于等于K ,则可以说信号X 是K 项稀疏。
如何找到信号最佳的稀疏域?这是压缩感知理论应用的基础和前提,只有选择合适的基表示信号才能保证信号的稀疏度,从而保证信号的恢复精度。
在研究信号的稀疏表示时,可以通过变换系数衰减速度来衡量变换基的稀疏表示能力。
Candes 和Tao 研究表明,满足具有幂次(power-law)速度衰减的信号,可利用压缩感知理论得到恢复。
最近几年,对稀疏表示研究的另一个热点是信号在冗余字典下的稀疏分解.这是一种全新的信号表示理论:用超完备的冗余函数库取代基函数,称之为冗余字典,字典中的元素被称为原子.字典的选择应尽可能好地符合被逼近信号的结构,其构成可以没有任何限制.从从冗余字典中找到具有最佳线性组合的K 项原子来表示一个信号,称作信号的稀疏逼近或高度非线性逼近。
目前信号在冗余字典下的稀疏表示的研究集中在两个方面:(1)如何构造一个适合某一类信号的冗余字典;(2)如何设计快速有效的稀疏分解算法。
这两个问题也一直是该领域研究的热点,学者们对此已做了一些探索,其中以非相干字典为基础的一系列理论证明得到了进一步改进.西安电子科技大学的石光明教授也对稀疏表示问题进行了认真研究,并基于多组正交基级联而成的冗余字典提出一种新的稀疏分解方法。
二 关于布莱克曼窗函数的研究 1. 布莱克曼窗 布莱克曼窗的时域形式可表示为:它的频域特性为:其中)(ωR W 为矩形窗函数的幅度频率特性。
一个理想数字滤波器的频率响应为Hd(ej ω),对应的时域序列为滤波器的单位脉冲响应hd(n),是无限长非因果的。
设计FIR-DF 就是要设计一个数字系统,去逼近理想数字滤波器的频率响应为Hd(ej ω)。
窗函数法就是对无限长的hd(n)加窗(用窗函数与之相乘,从而使之变成有限长的)下面主要介绍课程设计中用到的布莱克曼窗。
其频谱函数为:其幅度函数为这样其幅度函数由五部分组成。
它们都是移位不同,且幅度也不同的Wrg(w)函数,使旁瓣再进一步抵消。