数字图像处理DSP实验二第三章.docx
- 文档编号:5477447
- 上传时间:2022-12-16
- 格式:DOCX
- 页数:21
- 大小:171.38KB
数字图像处理DSP实验二第三章.docx
《数字图像处理DSP实验二第三章.docx》由会员分享,可在线阅读,更多相关《数字图像处理DSP实验二第三章.docx(21页珍藏版)》请在冰豆网上搜索。
数字图像处理DSP实验二第三章
第三章
M3.2求解并画出当N=10时,习题3.18中的序列的DTFT的实部和虚部以及幅度和相位谱
解析:
根据3.18题中计算得到的DTFT,分子分母化为exp(-j*w*n)的多项式之和,将分子分母每一项的系数分别组成矩阵num,dem,再利用freqz函数进行运算,可得的如下结果
【图形】
可以看出该频率响应的虚部实际为0,即其应为零相位谱,故使用zerophase函数更加准确的画出的其频率响应如下所示:
【图形】
解析:
本题采用定义式求其傅里叶变换,再变换为相应形式求得需要的系数向量,相应代码如下:
【代码】
%ReadinthedesirednumberoffrequencysamplesN
k=input('Numberoffrequencypoints=');
%ReadinthelengthN
N=input('N=');
nn=0:
(1/N):
1;
mm=(1-(1/N)):
(-1/N):
0
num=[nn,mm];
den=[zeros(1,N),1];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=freqz(num,den,w);
%Plotthefrequencyresponse
subplot(2,2,1)
plot(w/pi,real(h));grid
title('Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,2)
plot(w/pi,imag(h));grid
title('Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,3)
plot(w/pi,abs(h));grid
title('MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(h));grid
title('PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
同理,使用zerophase函数得出准确的频率响应如下
【代码】
k=input('Numberoffrequencypoints=');
%Readinthenumeratoranddenominatorcoefficients
N=input('N=');
nn=0:
(1/N):
1;
mm=(1-(1/N)):
(-1/N):
0;
num=[nn,mm];
den=[zeros(1,N),1,zeros(1,N)];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=zerophase(num,den,w);%求频率响应
%Plotthefrequencyresponse
plot(w/pi,real(h));grid
xlabel('\omega/\pi');ylabel('Amplitude')
【图形】
【图形】
同理,使用zerophase函数得出准确的频率响应如下
以上(a)~(d)四道题的代码如下,(输入时按照括号里相应的题号提示进行矩阵输入)
使用freqz函数的代码如下
%Program3_2
%Discrete-TimeFourierTransformComputation
%
%Readinthedesirednumberoffrequencysamples
k=input('Numberoffrequencypoints=');
%Readinthenumeratoranddenominatorcoefficients
num=input('Numeratorcoefficients=(a)[1zeros(1,41)-1](b)[1-1zeros(1,8)-1](c)[-0.5zeros(1,9)1zeros(1,9)-0.5](d)[-0.25000-11zeros(1,4)0.5zeros(1,4)1-1000-0.25]');
den=input('Denominatorcoefficients=(a)[zeros(1,20)10-1]...(b)[1,-1]...(c)[zeros(1,8)-5,10,-5](d)[zeros(1,9)-12-1]');
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=freqz(num,den,w);%求频率响应
%Plotthefrequencyresponse
subplot(2,2,1)
plot(w/pi,real(h));grid
title('Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,2)
plot(w/pi,imag(h));grid
title('Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,3)
plot(w/pi,abs(h));grid
title('MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(h));grid
title('PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
使用zerophase函数的代码如下:
k=input('Numberoffrequencypoints=');
%Readinthenumeratoranddenominatorcoefficients
num=input('Numeratorcoefficients=(a)[1zeros(1,41)-1]...(d)[-0.25000-11zeros(1,4)0.5zeros(1,4)1-1000-0.25]');
den=input('Denominatorcoefficients=(a)[zeros(1,20)10-1zeros(1,23)]...(d)[zeros(1,9)-12-1zeros(1,9)]');
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=zerophase(num,den,w);%求频率响应
%Plotthefrequencyresponse
plot(w/pi,real(h));grid
xlabel('\omega/\pi');ylabel('Amplitude')
解析:
本题发现其傅里叶变换并不好求,于是换种思路,采用定义式求其傅里叶变换,再变换为相应形式求得需要的系数向量,相应代码如下:
【代码】
%ReadinthedesirednumberoffrequencysamplesN
k=input('Numberoffrequencypoints=');
%ReadinthelengthN
N=input('N=');
nn=-N:
1:
N;
num=cos(pi/20.*nn);
den=[zeros(1,N),1];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=freqz(num,den,w);
%Plotthefrequencyresponse
subplot(2,2,1)
plot(w/pi,real(h));grid
title('Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,2)
plot(w/pi,imag(h));grid
title('Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,3)
plot(w/pi,abs(h));grid
title('MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(h));grid
title('PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
同理,使用zerophase函数得出准确的频率响应如下
【代码】
k=input('Numberoffrequencypoints=');
%Readinthenumeratoranddenominatorcoefficients
N=input('N=');
nn=-N:
1:
N;
num=cos(pi/20.*nn);
den=[zeros(1,N),1,zeros(1,N)];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=zerophase(num,den,w);%求频率响应
%Plotthefrequencyresponse
plot(w/pi,real(h));grid
xlabel('\omega/\pi');ylabel('Amplitude')
【图形】
M3.3画出如下DTFT的实部和虚部以及幅度和相位谱
解析:
下面两题的表示形式已经满足了freqz函数对分子分母表达形式的要求,直接可提取分子的系数num=0.1323.*[1,0.1444,-0.4519,0.1444,1],分母的系数den=[1,0.1386,0.8258,0.1393,0.4153]进行运算
【代码】
%ReadinthedesirednumberoffrequencysamplesN
k=input('Numberoffrequencypoints=');
num=0.1323.*[1,0.1444,-0.4519,0.1444,1];
den=[1,0.1386,0.8258,0.1393,0.4153];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=freqz(num,den,w);
%Plotthefrequencyresponse
subplot(2,2,1)
plot(w/pi,real(h));grid
title('Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,2)
plot(w/pi,imag(h));grid
title('Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,3)
plot(w/pi,abs(h));grid
title('MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(h));grid
title('PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
解析:
同理,直接可提取分子的系数num=0.3192.*[1,0.1885,-0.1885,-1],分母的系数den=[1,0.7856,1.4654,-0.2346],利用freqz函数进行运算
【代码】
%ReadinthedesirednumberoffrequencysamplesN
k=input('Numberoffrequencypoints=');
num=0.3192.*[1,0.1885,-0.1885,-1];
den=[1,0.7856,1.4654,-0.2346];
%Computethefrequencyresponse
w=0:
pi/(k-1):
pi;
h=freqz(num,den,w);
%Plotthefrequencyresponse
subplot(2,2,1)
plot(w/pi,real(h));grid
title('Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,2)
plot(w/pi,imag(h));grid
title('Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(2,2,3)
plot(w/pi,abs(h));grid
title('MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(h));grid
title('PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
M3.5验证表3.2列出的复序列的DTFT的对称关系
1.反转
【codes】
%property1
N=8;
a=0.5;
n=0:
N-1;
x=exp(j*a*n);
[Xw]=freqz(x,1,512);
y=exp(j*a*fliplr(n));
m=0:
511;w=-pi*m/512;
[Y1w]=freqz(y,1,w);
Y=exp(j*w*(N-1)).*Y1;
subplot(2,2,1)
plot(w/pi,abs(X));grid
title('x[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,2)
plot(w/pi,angle(X));grid
title('x[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
subplot(2,2,3)
plot(w/pi,abs(Y));grid
title('y[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(Y));grid
title('y[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
2.共轭反转对称
【codes】
%property2
N=8;
a=0.5;
n=0:
N-1;
x=exp(j*a*n);
[Xw]=freqz(x,1,512);
y=exp(-j*a*fliplr(n));
[Y1w]=freqz(y,1,512);
Y=conj(exp(j*w*(N-1)).*Y1);
subplot(2,2,1)
plot(w/pi,abs(X));grid
title('x[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,2)
plot(w/pi,angle(X));grid
title('x[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
subplot(2,2,3)
plot(w/pi,abs(Y));grid
title('y[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(Y));grid
title('y[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
性质三
【代码】
%property3
N=8;
a=0.5;
n=0:
N-1;
x=exp(j*a*n);
y=real(x);
[X0w]=freqz(x,1,512);%X0为x的DTFT
[Yw]=freqz(y,1,512);%Y为y的DTFT
m=0:
511;
w0=-pi*m/512;
[X1w]=freqz(x,1,w0);%X1为x(exp(-j*w))的DTFT
M=conj(X1);
X=0.5*(X0+M');
subplot(2,2,1)
plot(w/pi,abs(X));grid
title('x[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,2)
plot(w/pi,angle(X));grid
title('x[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
subplot(2,2,3)
plot(w/pi,abs(Y));grid
title('y[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(Y));grid
title('y[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
性质四
【代码】
%property3
N=8;
a=0.5;
n=0:
N-1;
x=exp(j*a*n);
y=j*imag(x);
[X0w]=freqz(x,1,512);%X0为x的DTFT
[Yw]=freqz(y,1,512);%Y为y的DTFT
m=0:
511;
w0=-pi*m/512;
[X1w]=freqz(x,1,w0);%X1为x(exp(-j*w))的DTFT
M=conj(X1);
X=0.5*(X0-M');
subplot(2,2,1)
plot(w/pi,abs(X));grid
title('x[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,2)
plot(w/pi,angle(X));grid
title('x[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
subplot(2,2,3)
plot(w/pi,abs(Y));grid
title('y[n]MagnitudeSpectrum')
xlabel('\omega/\pi');ylabel('Magnitude')
subplot(2,2,4)
plot(w/pi,angle(Y));grid
title('y[n]PhaseSpectrum')
xlabel('\omega/\pi');ylabel('Phase,radians')
【图形】
结论:
由图可以看出,性质三、四得到的傅里叶变换也满足对称关系,由于计算存在精度误差,此处未比较相位谱,采用比较的是实部和虚部,尽管图形不一定完全吻合,但也能基本看出其的大小关系。
性质五、六
【代码】
%property3
N=8;
a=0.5;
n=0:
N-1;
x=exp(j*a*n);
y=exp(j*a*fliplr(n));
xcs=0.5*[zeros(1,N-1)x]+0.5*[yzeros(1,N-1)];
xca=0.5*[zeros(1,N-1)x]-0.5*[yzeros(1,N-1)];
[Y1w]=freqz(xcs,1,512);
[Y2w]=freqz(xca,1,512);
[Xw]=freqz(x,1,512);%X为x的DTFT
Y1=Y1.*exp(j*w*(N-1));%Y1为xcs的DTFT
Y2=Y2.*exp(j*w*(N-1));%Y2为xca的DTFT
subplot(3,2,1)
plot(w/pi,real(X));grid
title('XRealpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(3,2,2)
plot(w/pi,imag(X));grid
title('XImaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(3,2,3)
plot(w/pi,real(Y1));grid
title('Y1Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(3,2,4)
plot(w/pi,imag(Y1));grid
title('Y1Imaginarypart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(3,2,5)
plot(w/pi,real(Y2));grid
title('Y2Realpart')
xlabel('\omega/\pi');ylabel('Amplitude')
subplot(3,2,
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字图像 处理 DSP 实验 第三
![提示](https://static.bdocx.com/images/bang_tan.gif)