随机信号功率的分析.docx
《随机信号功率的分析.docx》由会员分享，可在线阅读，更多相关《随机信号功率的分析.docx（13页珍藏版）》请在冰点文库上搜索。
随机信号功率的分析
本科学生实验报告
学号114090315姓名李开斌
学院物电学院专业、班级11电子
实验课程名称随机信号功率的分析
教师及职称李宏宁
开课学期2013至2014学年下学期
填报时间2014年4月17日
云南师范大学教务处编印
实验序号
06
实验名称
随机信号功率的频谱分析
实验时间
2014年4月16日
实验室
数字信号处理实验室（313）
一．实验预习
1．实验目的
深刻理解随机信号的特型，掌握随机信号功率频谱估计的基本原理，灵活运用各种随机信号功率频谱估计的基本方法。
2．实验原理、实验流程或装置示意图
实验原理：
功率谱估计是随机信号处理中的一个重要的研究和应用领域.功率谱估计基本上可以非参数估计的经典方法和参数估计的近代方法.典型功率谱估计是基于FFT算法的非参数估计,对足够长的记录数据效果较好。
在工程实际中，经典功率谱估计法获得广泛应用的是修正期图发。
该方法采取数据加窗处理再求平均的办法。
通过求各段功率谱平均，最后得到功率谱计P（m），
即：
式中：
为窗口函数ω[k]的方差。
K表示有重叠的分数段。
由于采用分段加窗求功率谱平均，有效地减少了方差和偏差，提高了估计质量，使修正周期图法在经典法中得到普遍应用。
但在估计过程存在两个与实际不符的假设，即
（1）利用有限的N个观察数据进行自相关估计，隐含着在已知N个数据之外的全部数据均为零的假设。
（2）假定数据是由N个观察数据以N为周期的周期性延拓。
同时在计算过程中采用加窗处理，使得估计的方差和功率泄露较大，频率分辨率较低，不适用于短系列的谱分析和对微弱信号的检测。
近代谱估计是建立在随机信号参数模型的基础上，通过信号参数模型或预测误差滤波器（一步预测器）参数的估计，实现功率谱估计。
由于既不需要加窗，又不需要对相关函数的估计进行如经典法那样的假设，从而减少公里泄露，提高了频谱分辨率。
常用的参数模型有自回归（AR）模型、滑动平均（MA）模型、自回归滑动平均（ARMA）模型。
其中AR模型是基本模型，求解AR模型的参数主要有L—D算法和Burg算法。
[例1.6.1]已知某随机信号由两个余弦信号和噪声构成
x（t）=cos（10*2πt）+cos（20*2πt）+s（t）
s（t）是均值为0，方差为1的白噪声。
使用周期图法和修正周期图法（Welch法）估计信号的功率谱。
【解】
周期图发是对随机系列的N点抽样数据直接经行DFT计算，然后求结果的模值得平方。
clc,clear,closeall
Fs=1000;
t=0:
1/Fs:
1;
figure
（1）;
x=sin（2*pi*50*t）+2*sin（2*pi*120*t）+randn（size（t））;
plot（t,x）;
运行结果如图1.6.1（a）所示。
粗略地估计该信号的功率谱，作N=1024点FFT：
figure
（2）;
Pxx=abs（fft（x,1024））.*2/1001;
subplot（2,1,1）;
plot（[0:
1023]*Fs/1024,10*log10（Pxx））;
axis（[0,1023,-30,30]）;
gridon
计算结果如图1.6.1（b）所示。
随着抽样数据（时窗宽度），逐渐接近无偏估计。
这种方法的问题是方差大且方差不随抽样点数增加而减小。
减小方差的办法是把抽样数据xn分段，每段做周期图，然后相加取平均值。
这种不重叠的周期图求平均的方法由Bart-lett于1953年引入，故又称为Bartlett法。
Pxx1=（abs（fft（x（1:
256）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（513:
768）））.^2）/（256*3）;
subplot（2,1,2）;
plot（[0:
255]*Fs/256,10*log10（Pxx1））;
axis（[0,1023,-30,30]）;
gridon
计算结果如图1.6.1（b）所示。
可见将数据分成3段，平均的方差是各段256点计算的1/3；随着分段的增加，方差减小。
但分段数不能太多，因为分段数增加会使得每断数据点数减少，造成频率分辨率下降。
为了减少偏差，提高频率分辨率，1970年Welch对Bartlett法进行了改进，提出了重叠平均周期图法（Welch法）。
一般选取重叠点数为分段点数的一半可以明显减小估计的方差。
若选取128点重叠，重新计算可得
figure（3）;
Pxx2=（abs（fft（x（1:
256）））.^2+abs（fft（x（129:
384）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（385:
640）））.^2+abs（fft（x（641:
896）））.^2）/（256*6）;
subplot（2,1,1）;
plot（[0:
255]*Fs/256,10*log10（Pxx2））;
axis（[0,1023,-30,30]）;
gridon
计算结果如图1.6.1（c）所示。
进一步提高周期图估计质量的方法是给每段数据加非矩形窗，然后再求各段的功率谱，即“修正周期图法”。
非矩形窗由于通过降低旁瓣的幅度从而减小旁瓣造成的“谱泄露”。
若使用256点汉宁窗重新分析得：
w=hanning（256）;
Pxx3=（abs（fft（x（1:
256）））.^2+abs（fft（x（129:
384）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（385:
640）））.^2+abs（fft（x（641:
896）））.^2）/（norm（w）^2*6）;
subplot（2,1,2）;
plot（[0:
255]*Fs/256,10*log10（Pxx3））;
axis（[0,1023,-30,30]）;
gridon
计算结果如图1.6.1（c）所示。
可见图中噪声的电平变得更加平坦，但谱峰变宽（只是由于非矩形窗的主瓣比矩形窗的主瓣宽而造成的）。
（a）时域信号波形
（b）分别用周期图法和分段的周期图法得到的功率谱
（c）分别用分段重叠的周期图法和分段重叠并使用汉明窗得到的功率谱
图1.6.1经典功率谱估计方法比较
在MATLAB信号处理工具箱中，提供了随机信号要功率谱估计的各段函数。
（1）periodogram函数可以实现周期图法的功率谱估计，起吊用格式为[Pxx,F]=PERIODOGRAM（x,WINDOW,NFFT,Fs）
其中：
x为进行功率谱估计的输入有限长序列；
WINDOW用于制定采用的窗函数，默认值为矩形窗（boxcar），窗函数的长度等于输入序列x的长度；
NFFT为DFT的点数，一般取大于输入序列x的长度，默认值为256；
Fs是绘制功率谱曲线的抽样频率，默认值为1；
Pxx为功率谱估计值；
F为Pxx值所对应的频率点。
（2）Welch-Barlett平均周期图法可以利用PSD函数实现，其调用格式为
[Pxx,F]=PSD（x,NFFT,Fs,WINDOW,NOVERLAP）
其中：
参数x，NFFT，Fs用法同periodogram函数：
WINDOW用于指定采用的窗函数，默认值为hanning窗；
NOVERLAP指定分段重叠的样函数。
如果使用boxcar窗且NOVERLAP=0，则可得到Barlett法的平均周期图。
如果NOVERLAP=length（x）/2，则可得到重叠50%的Welch法平均周期图。
（3）Welch法还可以利用pwelch函数实现，其调用格式为
[Pxx,F]=PWELCH（x,WINDOW,NOVERLAP,NFFT,Fs）
pwelch函数PSD函数进行功率谱估计的方法是一样的，只是参数设置略有不同。
（4）levinson函数可以实现L-D递推算法求解AR模型参数及白噪声系列的方差，其调用格式为
A=LEVINSON（R,ORDER）
其中：
ORDER为AR模型的阶数；
R为观察系列的自相关函数；
返回值A即为白噪声系列的方差和AR模型参数。
（5）aryule函数也可以实现L—D递推算法，他与levinson函数的区别仅在于输入参数不同，其调用格式为
[A,E,K]=ARBURG（x,ORDER）
其中：
x为观测序列；
E为预测误差；
K为反射系数；
（6）Burg算法可以利用函数arburg实现，其调用格式为
[A,E,K]=ARBURG（x,ORDER）
其中的参数与函数aryule相同。
3．实验设备及材料
装有MATLAB软件的计算机一台
4．实验方法步骤及注意事项
实验方法步骤：
（1）打开MATLAB软件
（2）根据题目要求编写程序
（3）运行程序
（4）分析实验结果
（5）关闭计算机
注意事项：
（1）对于实验仪器要轻拿轻放，遵守实验的规则。
（2）程序运行前要检查程序是否正确。
二．实验内容
1.某随机信号由两余弦信号与噪声构成x（t）=cos（20*pi*t）+cos（40*pi*t）+s（t）式中：
s（t）是均值为0、方差为1的高斯白噪声。
（1）生成信号x（t）=cos（20*pi*t）+cos（40*pi*t）+s（t），绘出其时域波形。
（2）试分别用周期图发、bartlett发和welch发分析该序列的功率频谱估计，并对上述各种的估计效果经行比较。
（3）选用不同的窗函数，利用修正周期图发估计信号的功率谱有何结论？
【解】
clc,clear,closeall
Fs=1000;
t=0:
1/Fs:
1;
figure
（1）;
x=sin（2*pi*50*t）+2*sin（2*pi*120*t）+randn（size（t））;
plot（t,x）;
figure
（2）;
Pxx=abs（fft（x,1024））.*2/1001;
subplot（2,1,1）;
plot（[0:
1023]*Fs/1024,10*log10（Pxx））;
axis（[0,1023,-30,30]）;
gridon
Pxx1=（abs（fft（x（1:
256）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（513:
768）））.^2）/（256*3）;
subplot（2,1,2）;
plot（[0:
255]*Fs/256,10*log10（Pxx1））;
axis（[0,1023,-30,30]）;
gridon
figure（3）;
Pxx2=（abs（fft（x（1:
256）））.^2+abs（fft（x（129:
384）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（385:
640）））.^2+abs（fft（x（641:
896）））.^2）/（256*6）;
subplot（2,1,1）;
plot（[0:
255]*Fs/256,10*log10（Pxx2））;
axis（[0,1023,-30,30]）;
gridon
w=boxcar（256）;
Pxx3=（abs（fft（x（1:
256）））.^2+abs（fft（x（129:
384）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（385:
640）））.^2+abs（fft（x（641:
896）））.^2）/（norm（w）^2*6）;
subplot（2,1,2）;
plot（[0:
255]*Fs/256,10*log10（Pxx3））;
axis（[0,1023,-30,30]）;
gridon
w=hamming（256）;
Pxx3=（abs（fft（x（1:
256）））.^2+abs（fft（x（129:
384）））.^2+abs（fft（x（257:
512）））.^2+abs（fft（x（385:
640）））.^2+abs（fft（x（641:
896）））.^2）/（norm（w）^2*6）;
subplot（2,1,2）;
plot（[0:
255]*Fs/256,10*log10（Pxx3））;
axis（[0,1023,-30,30]）;
gridon
2．使用MATLAB函数loadmtlb调用MATLAB中的一段语音信号，利用AR模型估计其功率谱，并讨论不同阶数的AR模型岁随机信号功率谱估计的影响。
2．对实验现象、实验结果的分析及其结论
本试验中信号的产生、提取均用MATLAB编程实现。
让我们学习到MATLAB在信号分析中的运用同时也让我们认识到了MATLAB的实用性与重要性。
在对随机信号的分析过程中，我们了解到了随机信号分析理论在实践中的运用加深了对随机信号各个方面知识的理解。
在进行微弱信号的提取过程中了解到多种提取信号的方法同时学会了多重自相关提取强噪声下弱信号的方法。
教师评语及评分：
签名：
年月日
升级会员 