DSP研究性学习报告Word文档格式.docx
《DSP研究性学习报告Word文档格式.docx》由会员分享,可在线阅读,更多相关《DSP研究性学习报告Word文档格式.docx(17页珍藏版)》请在冰豆网上搜索。
![DSP研究性学习报告Word文档格式.docx](https://file1.bdocx.com/fileroot1/2023-1/6/f74831ab-9013-41e6-a245-040fed0d6bb3/f74831ab-9013-41e6-a245-040fed0d6bb31.gif)
【仿真程序】
N=30;
L=512;
f1=100;
f2=120;
fs=600;
T=1/fs;
t=(0:
N-1)*T;
x=cos(2*pi*f1*t)+cos(2*pi*f2*t);
y=fft(x,L);
mag=abs(y);
f=(0:
length(y)-1)*fs/length(y);
plot(f(1:
L/2),mag(1:
L/2));
N=600;
N=20;
A=1;
B=1用哈明窗
%数据的长度
%DFT的点数
%抽样频率
%抽样间隔
ws=2*pi*fs;
wh=(hamming(N));
X=x*wh;
X=fftshift(fft(x,L));
w=(-ws/2+(0:
L-1)*ws/L)/(2*pi);
plot(w,abs(X))
使用矩形窗
wh=(boxcar(N));
A=1,B=0.2
矩形窗
%Ê
ý
¾
Ý
µ
Ä
³
¤
¶
È
%DFTµ
ã
Ê
%³
é
Ñ
ù
Æ
Â
¼
ä
¸
ô
x=cos(2*pi*f1*t)+0.2*cos(2*pi*f2*t);
用哈明窗
2.分析矩形窗、Hann窗、Hamming窗、Blackman窗、Kaiser窗的频谱,并进行比较。
【题目分析】
矩形窗属于时间变量的零次幂窗。
矩形窗使用最多,习惯上不加窗就是使信号通过了矩形窗。
这种窗的优点是主瓣比较集中,缺点是旁瓣较高,并有负旁瓣,导致变换中带进了高频干扰和泄漏,甚至出现负谱现象。
Hanning窗
汉宁窗又称升余弦窗,汉宁窗可以看作是3个矩形时间窗的频谱之和,或者说是3个sinc(t)Hamming窗型函数之和,而括号中的两项相对于第一个谱窗向左、右各移动了π/T,从而使旁瓣互相抵消,消去高频干扰和漏能。
可以看出,汉宁窗主瓣加宽并降低,旁瓣则显著减小,从减小泄漏观点出发,汉宁窗优于矩形窗.但汉宁窗主瓣加宽,相当于分析带宽加宽,频率分辨力下降
Hamming窗
海明窗也是余弦窗的一种,又称改进的升余弦窗。
海明窗与汉宁窗都是余弦窗,只是加权系数不同。
海明窗加权的系数能使旁瓣达到更凯泽(Kaiser)窗
小。
分析表明,海明窗的第一旁瓣衰减为一42dB.海明窗的频谱也是由3个矩形时窗的频谱合成,但其旁瓣衰减速度为20dB/(10oct),这比汉宁窗衰减速度慢。
海明窗与汉宁窗都是很有用的窗函数。
Kaiser窗是一种可以任意调整窗谱主瓣宽度与旁瓣衰减比例的函数。
以瞬时无功功率理论为基础,提出一种谐波电流的实时检测新方法,用Kaiser滑动时窗截取谐波电流信号,通过对窗函数参数的选定,提高了谐波电流检测的准确性。
1.矩形窗:
利用w=boxcar(n)的形式得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,它的元素由窗函数的值组成。
‘w=boxcar(n)’等价于‘w=ones(1,n)’.
2.三角窗:
利用w=triang(n)的形式得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,它的元素由窗函数的值组成。
w=triang(N-2)等价于bartlett(N)。
3.汉宁窗:
利用w=hanning(n)得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,包含了窗函数的n个系数。
4.海明窗:
利用w=hamming(n)得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,包含了窗函数的n个系数。
它和汉宁窗的主瓣宽度相同,但是它的旁瓣进一步被压低。
5.布拉克曼窗:
利用w=blackman(n)得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,包含了窗函数的n个系数。
它的主瓣宽度是矩形窗主瓣宽度的3倍,为12*pi/N,但是它的最大旁瓣值比主瓣值低57dB。
6.切比雪夫窗:
它是等波纹的,利用函数w=(N,R)方式设计出N阶的切比雪夫2窗函数,函数的主瓣值比旁瓣值高RdB,且旁瓣是等波纹的。
7.巴特里特窗:
利用w=bartlett(n)的形式得到窗函数,其中n为窗函数的长度,而返回值w为一个n阶的向量,包含了窗函数的n个系数。
8.凯塞窗:
利用w=kaiser(n,beta)的形式得到窗函数。
1.矩形窗(RectangleWindow)调用格式:
w=boxcar(n),根据长度n产生一个矩形窗w。
K=300;
t=0:
K-1;
w=boxcar(N);
Y=abs(fftshift(fft(w,K)));
plot(t,Y);
2.三角窗(TriangularWindow)调用格式:
w=triang(n),根据长度n产生一个三角窗w。
w=triang(N);
3.汉宁窗(HanningWindow)调用格式:
w=hanning(n),根据长度n产生一个汉宁窗w。
w=hanning(N);
4.海明窗(HammingWindow)调用格式:
w=hamming(n),根据长度n产生一个海明窗w。
w=hamming(N);
5.布拉克曼窗(BlackmanWindow)调用格式:
w=blackman(n),根据长度n产生一个布拉克曼窗w。
w=blackman(N);
6.恺撒窗(KaiserWindow)调用格式:
w=kaiser(n,beta),根据长度n和影响窗函数旁瓣的β参数产生一个恺撒窗w。
w=kaiser(N,50);
【研讨题目】
3.利用DFT分析x(t)=e-|t|u(t)的频谱。
【题目分析】首先利用CTFT公式计算其模拟频谱的理论值;
然后对其进行等间隔理想采样,得到序列,利用DTFT公式计算采样序列的数字连续频谱理论值,通过与模拟频谱的理论值对比,理解混叠误差形成的原因及减小误差的措施.接下来是对序列进行加窗处理,得到有限长加窗序列,再次利用DTFT公式计算加窗后序列的数字连续频谱,并与加窗前的数字连续频谱进行对比,理解截断误差形成的原因及减小误差的措施;
最后是对加窗序列进行DFT运算,得到加窗后序列的DFT值,它是对数字连续频谱进行等间隔采样的采样值.
【仿真结果】
该程序过程清晰、容易理解,程序运行结果如图2所示。
图2(a)是信号的时域波形,图2(b)是对应的幅度谱图。
由于在归一化频率为0.5的地方还有较大幅度,所以对信号进行理想采样后存在较大的混叠失真,表现在图2(d)中归一化频率为-0.5~0.5范围内的波形与图2(b)的波形明显不同。
图2(e)是理想采样序列加矩形窗得到的图形,对应的幅度谱线如图2(f)中实线所示。
该谱线与图2(d)相比有明显的不同(如波动现象),这是加窗截断的结果。
窗口长度越短截断效应会越明显。
对加窗序列进行DFT运算,只要DFT点数大于等于窗口长度,算出的DFT值就是对加窗序列的连续频谱在一个周期内进行的等间隔采样的采样值。
在图2(f)中表现为DFT幅值在加窗序列连续幅谱的谱线上,是连续频谱曲线上的有限个数据点,即所谓栅栏效应。
图2(f)中用虚线画出了连续信号的幅谱理论值,与DFT的幅值对比,存在一定误差,但只要采样频率足够高,时窗长度足够长,FFT点数足够大,得到的DFT值越逼近实际频谱。
【自主学习内容】
clc;
closeall;
clear;
%CTFT程序,以x(t)=exp(-t)t>
=0为例
%利用数值运算计算并绘制连续信号波形
L=4,%定义信号波形显示时间长度
fs=4,T=1/fs;
%定义采样频率和采样周期
t_num=linspace(0,L,100);
%取若干时点,点数决定作图精度
xt_num=exp(-1*t_num);
%计算信号在各时点的数值
subplot(3,2,1);
plot(t_num,xt_num),%绘信号波形
xlabel('
时间(秒)'
),ylabel('
x(t)'
),%加标签
grid,title('
(a)信号时域波形'
),%加网格和标题
%利用符号运算和数值运算计算连续信号幅度谱的理论值
symstW%定义时间和角频率符号对象
xt=exp(-1*t)*heaviside(t),%连续信号解析式
XW=fourier(xt,t,W),%用完整调用格式计算其傅氏变换
%在0两边取若干归一化频点,点数决定作图精度
w1=[linspace(-0.5,0,50),linspace(0,1.5,150)];
XW_num=subs(XW,W,w1*2*pi*fs);
%利用置换函数求频谱数值解
mag1=abs(XW_num);
%计算各频点频谱的幅值
subplot(3,2,2);
plot(w1,mag1),%绘制归一化频率幅值谱线
频率(*2*pi*fs)rad/s'
幅度'
(b)连续信号幅频理论值'
%DTFT程序,以x(n)=exp(-nT)n>
%利用数值运算计算并绘制离散信号图形
N=L*fs+1;
n_num=0:
N-1;
%生成信号波形采样点
xn_num=exp(-1*n_num*T);
%计算信号理想采样后的序列值
subplot(3,2,3);
stem(n_num,xn_num,'
b.'
),%绘序列图形
n'
x(n)'
(c)理想采样图形'
%利用符号运算和数值运算计算离散信号幅度谱的理论值
symsnzw%定义符号对象
xn=exp(-n*T),%定义离散信号
Xz=ztrans(xn,n,z),%用完整调用格式计算其Z变换
%利用复合函数计算序列傅里叶变换的解析解
Z=exp(j*w);
Hejw=compose(Xz,Z,z,w);
Hejw_num=subs(Hejw,w,w1*2*pi);
%求频谱数值解
mag2=abs(Hejw_num);
%计算各频点频谱的幅度
subplot(3,2,4);
plot(w1,mag2*T),%绘制频谱幅度曲线
频率(*2*pi)rad'
(d)离散信号幅频理论值'
%序列加窗图示及频谱幅值绘制
%利用数值运算计算并绘制加窗后序列xw(n)的图形
M=8;
win=(window(@rectwin,M))'
;
%定义窗点和窗型
xwn=xn_num.*[win,zeros(1,N-M)];
%给离散信号加窗
subplot(3,2,5);
stem(n_num,xwn,'
),%加窗序列图示
xw(n)'
(e)加窗序列图形'
%利用符号运算和数值运算计算加窗序列的频谱幅值
%先求加窗序列的Z变换,注意表达式长度限制问题
Xwz=0;
forn=0:
(M-1);
Xwz=Xwz+xwn(n+1)*z^(-n);
end
%利用复合函数计算加窗序列傅里叶变换的解析解
Zw=exp(j*w);
HejwM=compose(Xwz,Zw,z,w);
HejwM_num=subs(HejwM,w,w1*2*pi);
mag3=abs(HejwM_num);
subplot(3,2,6);
plot(w1,mag3*T),%绘频谱幅度曲线
%利用DFT计算加窗序列xw(n)的离散谱幅值
Ndft=16,Xk=fft(xwn,Ndft);
%定义DFT点数和DFT运算
Xk0=fftshift(Xk)*T;
%将DFT值0对称和幅值加权处理
ifmod(Ndft,2)==0;
N1=Ndft;
elseN1=Ndft-1;
end;
k=[0:
(Ndft-1)]-N1/2;
wk=k/Ndft;
%0对称取值并归一化
holdon;
stem(wk,abs(Xk0),'
r.'
),%绘制DFT图形
legend('
幅谱'
'
DFT'
0),%加响应图例,位置自动最佳
归一化频率'
(f)加窗序列幅谱及其DFT幅值'
),
plot(w1,mag1,'
k:
'
),%与连续信号幅谱的理论值比较
【研讨题目】
4.利用DFT分析计算x(t)=cos2pit的频谱,设w0=pi/2
【题目分析】因为f=w/2pi=1Hz,所以抽样频率需要大于2Hz.
【序列频谱计算的基本方法】
【仿真结果】
【结果分析】
【自主学习内容】已知一连续信号为
其中
。
试利用DFT分析其频谱。
解:
信号
的最高频率
,抽样频率
,取抽样频率
;
最低的频率分辨率为
,最少的信号样点数为
的MATLAB程序如下:
%programexa_1_1.m,利用矩形窗计算有限长余弦信号频谱
length(y)-1)'
*fs/length(y);
频率(Hz)'
)
ylabel('
幅度谱'
程序运行结果如下图所示。
由图可见,频谱图显示出两个较为明显的峰值(对应
)。
结论:
当截取信号
样点时,频率分辨率
,刚好能够分辨出
和
两个频谱分量,但频谱泄漏较严重。
若取
(达不到最低的频率分辨率)。
%programexa_1_2.m,利用矩形窗计算有限长余弦信号频谱
由图可见,频谱图显示不出对应
的两个明显峰值。
,达不到最低的频率分辨率
,故分辨不出
两个频谱分量,且频谱泄漏更为严重。
若取频率分辨率
,则对应的信号样点数为
%programexa_1_3.m,利用矩形窗计算有限长余弦信号频谱
y=fft(x);
length(y)/2),mag(1:
length(y)/2));
由图可见,频谱图显示出两个特别明显的峰值(对应
,高分辨率的频谱图具有较高的质量,频谱分析时必须保证获取足够的信号数据长度。
L=8,%定义信号波形显示时间长度
xt_num=cos(2*pi*t_num);
xt=cos(2*pi*t)*heaviside(t),%连续信号解析式
xn_num=cos(2*pi*n_num*T);
%