滤波与功率谱估计
滤波与功率谱估计
清华大学《数字信号处理》期末作业2013 年 1 月第一题掌握去噪的方法1.1 题目描述MATLAB 中的数据文件noisdopp 含有噪声，该数据的抽样频率未知。
调出该数据，用你学过的滤波方法和奇异值分解的方法对其去噪。
要求：1．尽可能多地去除噪声，而又不损害原信号；2．给出你去噪的原理与方法；给出说明去噪效果的方法或指标；3．形成报告时应包含上述内容及必要的图形，并附上原程序。
1.2 信号特性分析MATLAB所给noisdopp信号极其频域特征如图1.1、图1.2。
图1.1含有噪声的noisdopp信号图1.2 noisdopp 信号频域特性其中横坐标f 采用归一化频率，即未知抽样频率Fs 对应2（与滤波器设计时参数一致）。
信号基本特性是一个幅值和频率逐渐增加的正弦信号叠加噪声，噪声为均匀的近似白噪声，没有周期等特点。
因为噪声信号能量在全频带均匀分布，滤波器截止频率过低则信号损失大，过高则噪声抑制小，认为频谱中含有毛刺较多的部分即为信噪比较小的部分，滤除这部分可以达到较好的滤波效果。
先给定去噪效果的评定指标。
信号开始阶段频率较高（如图1.3，红圈为信号值），一周期内采样点4~5个，即信号归一化频率达到0.4~0.5（Fs=2），难以从频域将这部分信号同噪声分离，滤波后信号损失较大，故对前128点用信噪比考察其滤波效果，定义：22()10lg (()())k k xk SN R y k x k =-∑∑其中，()x k 为原nosidopp 信号，()y k 为滤波后信号。
SNR 越大表示滤除部分能力越小，可以反映滤波后信号对原信号的跟踪能力，对前128点主要考察SNR ，越大滤波器性能越好。
对128点以后的点，信号能量集中在低频部分，滤波后效果不仅仅需考察SNR ，同时考察滤波后信号的平滑程度。
定义平滑指标[1]22((1)())((1)())k k y k y k r x k x k +-=+-∑∑r 越小表示滤波后信号相比原信号更加平滑，滤波效果更好。
功率谱计算[解说]
![功率谱计算[解说]](https://img.taocdn.com/s3/m/3cec2a18f02d2af90242a8956bec0975f465a478.png)
功率谱计算功率谱估计在现代信号处理中是一个很重要的课题，涉及的问题很多。
在这里，结合matlab，我做一个粗略介绍。
功率谱估计可以分为经典谱估计方法与现代谱估计方法。
经典谱估计中最简单的就是周期图法，又分为直接法与间接法。
直接法先取N点数据的傅里叶变换（即频谱），然后取频谱与其共轭的乘积，就得到功率谱的估计；间接法先计算N点样本数据的自相关函数，然后取自相关函数的傅里叶变换，即得到功率谱的估计.都可以编程实现，很简单。
在matlab中，周期图法可以用函数periodogram实现。
但是周期图法估计出的功率谱不够精细，分辨率比较低。
因此需要对周期图法进行修正，可以将信号序列x(n)分为n个不相重叠的小段，分别用周期图法进行谱估计，然后将这n段数据估计的结果的平均值作为整段数据功率谱估计的结果。
还可以将信号序列x(n)重叠分段，分别计算功率谱，再计算平均值作为整段数据的功率谱估计。
这2种称为分段平均周期图法，一般后者比前者效果好。
加窗平均周期图法是对分段平均周期图法的改进，即在数据分段后，对每段数据加一个非矩形窗进行预处理，然后在按分段平均周期图法估计功率谱。
相对于分段平均周期图法，加窗平均周期图法可以减小频率泄漏，增加频峰的宽度。
welch法就是利用改进的平均周期图法估计估计随机信号的功率谱，它采用信号分段重叠，加窗，FFT等技术来计算功率谱。
与周期图法比较，welch法可以改善估计谱曲线的光滑性，大大提高谱估计的分辨率。
matlab中，welch法用函数psd实现。
调用格式如下：[Pxx,F] = PSD(X,NFFT,Fs,WINDOW,NOVERLAP)X：输入样本数据NFFT：FFT点数Fs:采样率WINDOW：窗类型NOVERLAP，重叠长度现代谱估计主要针对经典谱估计分辨率低和方差性不好提出的，可以极大的提高估计的分辨率和平滑性。
可以分为参数模型谱估计和非参数模型谱估计。
参数模型谱估计有AR模型，MA模型，ARMA模型等；非参数模型谱估计有最小方差法和MUSIC法等。
matlab中 功率谱估计的函数
在matlab中，功率谱估计是信号处理和频谱分析中常用的一种方法。
通过对信号的频谱特性进行估计，可以有效地分析信号的功率分布情况，从而为信号处理和系统设计提供重要的参考信息。
在matlab中，提供了多种功率谱估计的函数，以下将对其中几种常用的函数进行介绍和分析。
1. periodogram函数periodogram函数是matlab中用于估计信号功率谱密度的函数之一。
它基于傅里叶变换将离散时间信号转换成频域信号，然后计算频域信号的功率谱密度。
其调用格式为：[Pxx, F] = periodogram(x,window,nfft,fs)其中，x为输入的离散时间信号，window为窗函数，nfft为离散傅里叶变换的点数，fs为信号的采样频率。
periodogram函数返回的Pxx 为功率谱密度估计值，F为对应的频率。
2. pwelch函数pwelch函数也是用于估计功率谱密度的函数，它采用了Welch方法，通过对信号进行分段处理，然后对各段信号进行傅里叶变换，并对各段功率谱密度进行平均。
其调用格式为：[Pxx, F] = pwelch(x,window,noverlap,nfft,fs)其中，x为输入的离散时间信号，window为窗函数，noverlap为相邻分段的重叠点数，nfft为离散傅里叶变换的点数，fs为信号的采样频率。
pwelch函数返回的Pxx为功率谱密度估计值，F为对应的频率。
3. cpsd函数cpsd函数用于估计信号的交叉功率谱密度，即两个信号之间的频谱特性。
其调用格式为：[Pxy, F] = cpsd(x,y,window,noverlap,nfft,fs)其中，x和y为输入的两个离散时间信号，window为窗函数，noverlap为相邻分段的重叠点数，nfft为离散傅里叶变换的点数，fs为信号的采样频率。
cpsd函数返回的Pxy为交叉功率谱密度估计值，F为对应的频率。
4. mscohere函数mscohere函数用于估计信号的相干函数，即两个信号之间的相关性。
第6讲 谱估计2.传统法
作长度为L>2M－1的序列：
RN (m) w(m) Qx (m) = 0 R (m L) w(m L) N 0 ≤ m ≤ M 1 M ≤ m ≤ L M +1 L M +1 < m ≤ L 1
求功率谱估计 FFT[Qx (m)] = S x (k 2π / L) k = 0,1, , L 1 2．估计质量 ． 当M<<N-1，BT法不等于周期图法，而是对周期图 的平滑，以下讨论BT谱或者周期图平滑谱的性能。
加窗后的功率谱估计值为：
S x (ω ) =
w(m) W (e )
jω
M 1
m = ( M 1)
w ( m ) R N ( S x (ω ) = ，有 RN (m) S (e jω ) = I N (ω ) 2π
∧
∫π
π
I N (θ )W (e j (ω θ ) )dθ
4.7.1 间接法
间接法是由随机信号N个观察值 x(0),, x( N 1) ，估计 出自相关函数 RN (m) ，然后求其傅里叶变换作为功 率谱的估计。
S (e ) =
∧ jω m= M
RN ( m )e jωm ∑
M
M ≤ N 1
4.7.2 直接法
设有限长实序列
x( n) (0 ≤ n ≤ N 1) 的傅里叶变换为：
ω
∧
∧
线性系统幅度频率响应的平方是 S x (e jω )，因此
sin[ωN ]2 2 此时谱估计方差为： [S (ω)] = S x2 (e jω )1 + Var N sin ω
∧
S (ω ) ≈ S x (e ) S ω (ω )
∧
jω
功率谱和功率谱密度计算公式
功率谱和功率谱密度计算公式
功率谱（Power Spectrum）
是描述随机信号或时间序列在不同频率下功率分布情况的工具。
对于离散信号，功率谱的计算通常涉及到傅里叶变换（Fourier Transform）或者更一般的傅里叶分析方法。
假设有一个离散信号(x(n))（其中(n)表示时间或样本序号），其功率谱(P(f))可以通过以下步骤计算：
傅里叶变换：首先，对信号(x(n))进行傅里叶变换，得到其频谱(X(f))：
(X(f) = \sum_{n=-\infty}^{\infty} x(n) e^{-j2\pi fn})
计算功率谱：然后，计算频谱的模的平方，即得到功率谱(P(f))：
(P(f) = |X(f)|^2)
功率谱密度（Power Spectral Density, PSD）
是单位频率范围内的平均功率，通常用于描述连续信号的功率分布。
对于连续信号(x(t))（其中(t)表示时间），其功率谱密度(S_{xx}(f))可以通过自相关函数和傅里叶变换得到：
自相关函数：首先，计算信号(x(t))的自相关函数(R_{xx}(\tau))：
(R{xx}(\tau) = \int{-\infty}^{\infty} x(t) x(t+\tau) dt)
傅里叶变换：然后，对自相关函数(R{xx}(\tau))进行傅里叶变换，得到功率谱密度(S{xx}(f))：(S{xx}(f) = \int{-\infty}^{\infty} R_{xx}(\tau) e^{-j2\pi f\tau} d\tau)。
welch方法计算功率谱
welch方法计算功率谱Welch方法是一种常用的计算功率谱密度估计的方法。
它使用了平均窗口的思想，通过将信号分段并对这些分段信号进行窗口加权和傅里叶变换，最后取各个分段的平均值来计算功率谱密度。
Welch方法具有计算简单、可靠性高的特点，在信号处理和频谱分析领域得到广泛应用。
Welch方法的基本原理是将原始信号分成多个重叠的窗口，然后对每个窗口内的信号进行加权和傅里叶变换。
加权的目的是为了在频率分析中减少数据窗口两端的边界效应。
常用的加权函数有汉宁窗、汉明窗等。
在傅里叶变换之后，得到的每个窗口的功率谱密度为傅里叶变换结果的模长的平方。
然后，将各个窗口的功率谱密度求取平均，最终得到整个信号的功率谱密度估计。
为了进一步提高频谱分辨率，可以将各个窗口的重叠程度增加，例如50%的重叠。
通常情况下，窗口的大小和重叠程度会根据具体的应用领域和信号特性进行调整。
Welch方法的主要优点是计算简单快速，并且可以在频谱分析中减少噪声的影响。
然而，由于Welch方法是基于分段信号的平均，所以在信号的时间变化较大或者边界效应较明显的情况下，可能会引入一定的偏差。
因此，在具体应用中需要根据实际情况对窗口大小和重叠程度进行适当的选择。
Welch方法常见的应用领域包括信号处理、通信系统、音频处理等。
在信号处理领域，Welch方法可以用于频谱分析、滤波器设计、噪声分析等。
在通信系统中，Welch方法可以用于无线信号的功率谱分析，用于干扰分析和频域特征提取等。
在音频处理中，Welch方法可以用于语音信号的频谱分析和音频信号的特征提取等。
总结来说，Welch方法是一种基于信号分段和平均的功率谱密度估计方法。
它通过将信号分成多个窗口，并对每个窗口内的信号进行加权和傅里叶变换来计算功率谱密度。
Welch方法具有计算简单、可靠性高的特点，在信号处理和频谱分析领域得到广泛应用。
MATLAB技术谱估计方法
MATLAB技术谱估计方法一、介绍MATLAB（Matrix Laboratory）是一种用于数值计算和数据可视化的编程语言和环境。
它提供了各种工具和函数，用于解决科学和工程领域的问题。
在信号处理领域，MATLAB被广泛应用于谱估计方法的研究和实践。
本文将详细介绍MATLAB中常用的技术谱估计方法。
二、时频分析及其应用时频分析是一种将信号在时间和频率上进行联合分析的方法。
它可以揭示信号在时间和频率域上的变化规律，对于非平稳信号的分析非常有用。
在MATLAB中，可以利用一些函数实现时频分析，如“specgram”函数可以计算信号的谱矩阵，并将其绘制成频谱图。
三、功率谱密度估计功率谱密度估计是一种用于描述信号在频域上的能量分布的方法。
MATLAB提供了多种功率谱密度估计方法，如传统的周期图方法和现代的非周期方法。
其中，最常用的有Welch方法和Yule-Walker方法。
Welch方法基于信号的分段平均，将信号划分为多个段落，然后计算每个段落的谱密度估计，最后对这些估计值进行平均。
Yule-Walker方法则基于自相关函数和线性预测，在频域上对信号进行建模。
四、自相关函数及其应用自相关函数是一种描述信号中自身相关性的方法。
在MATLAB中，可以使用“xcorr”函数计算信号的自相关函数。
自相关函数广泛应用于信号处理中的一些任务，例如信号的频率估计、信号的滤波和信号的预测。
通过自相关函数，我们可以推导出信号的自相关峰值位置，从而得到信号的周期和频率。
五、谱减法和频率追踪谱减法是一种常用的去噪方法，它利用信号的谱特性对噪声进行削减。
在MATLAB中，可以使用“pwelch”函数计算信号的功率谱密度估计，然后利用简单的减法运算去除噪声。
频率追踪是一种用于追踪信号频率变化的方法，对于非稳定信号的分析尤为重要。
MATLAB提供了多种频率追踪方法，如“pmtm”函数和“arburg”函数，可以帮助我们准确地追踪信号的频率变化。
自-信号检测与估计试题——答案(不完整版)
一、概念:1. 匹配滤波器。
概念：所谓匹配滤波器是指输出判决时刻信噪比最大的最佳线性滤波器。
应用:在数字信号检测和雷达信号的检测中具有特别重要的意义。
在输出信噪比最大准则下设计一个线性滤波器是具有实际意义的。
2. 卡尔曼滤波工作原理及其基本公式(百度百科)首先,我们先要引入一个离散控制过程的系统。
该系统可用一个线性随机微分方程（Linear StoｃhastiｃDifferｅｎce equatiｏn）来描述：X(k）=A X(k－1)+BＵ(k)+W(ｋ)再加上系统的测量值:Z(k)＝ＨＸ(k)+V(k)上两式子中，Ｘ(k)是k时刻的系统状态,U(k）是k时刻对系统的控制量。
A和B是系统参数，对于多模型系统,他们为矩阵。
Z(k)是k时刻的测量值,H是测量系统的参数，对于多测量系统,H为矩阵。
Ｗ(k)和V（k)分别表示过程和测量的噪声。
他们被假设成高斯白噪声(Whｉte ＧausｓｉaｎNoiｓe），他们的coｖaｒiancｅ分别是Q，R（这里我们假设他们不随系统状态变化而变化）。
对于满足上面的条件(线性随机微分系统,过程和测量都是高斯白噪声),卡尔曼滤波器是最优的信息处理器。
下面我们来用他们结合他们的covariａncｅs来估算系统的最优化输出（类似上一节那个温度的例子）。
首先我们要利用系统的过程模型，来预测下一状态的系统。
假设现在的系统状态是k,根据系统的模型,可以基于系统的上一状态而预测出现在状态：Ｘ(k|ｋ－1）=A X(k-１|k－1）+B U(k) ……….．(1)式(１)中，X（k|k-1)是利用上一状态预测的结果，X（k-１|k－１)是上一状态最优的结果,Ｕ(k)为现在状态的控制量，如果没有控制量,它可以为0。
到现在为止，我们的系统结果已经更新了,可是,对应于X(ｋ｜k-1)的covａｒiancｅ还没更新。
我们用P表示covａｒianｃe:P(k|ｋ-１)=A P(k-1｜ｋ－１)Ａ’+Q (2)式(2)中,P（k｜ｋ－1）是X（k|k-1)对应的covaｒｉaｎce,P(ｋ-１|k-1）是Ｘ(k-1|k-1)对应的cｏvariaｎce,A’表示A的转置矩阵,Q是系统过程的cｏvａriance。
经典功率谱估计
雷达和声呐系统
目标检测
在雷达和声呐系统中，经典功率谱估计常被用于目标检测。通过对接收到的信号进行功率 谱分析，可以判断是否存在目标以及目标的位置和速度等信息。
距离和速度测量
在雷达和声呐系统中，经典功率谱估计还可以用于距离和速度测量。通过对接收到的信号 进行功率谱分析，可以估计出目标与系统之间的距离和相对速度。
信号分类
在雷达和声呐系统中，经典功率谱估计还可以用于信号分类。通过对接收到的信号进行功 率谱分析，可以判断目标的类型，例如区分飞机、船舶或车辆等不同类型目标。
05 经典功率谱估计的改进方 法
基于小波变换的功率谱估计
1
小波变换能够将信号分解成不同频率和时间尺度 的分量，从而更好地揭示信号的内在结构和特征。
然而，这些方法通常需要较长 的数据长度和较为复杂的计算 过程，对于短数据和实时处理 的应用场景具有一定的局限性 。
研究展望
01
随着信号处理技术的发展，经典功率谱估计方法仍有进一步优化的空 间。
02
针对短数据和实时处理的应用场景，研究更为快速、准确的功率谱估 计方法具有重要的实际意义。
03
结合机器学习和人工智能技术，探索基于数据驱动的功率谱估计方法 是一个值得关注的方向。
优点
能够提供较高的频率分辨率和较低的估计误差。
原理
格莱姆-梅尔谱估计利用了信号的模型参数，通过 构造一个模型函数来描述信号的频率响应特性， 并求解该函数的极值问题得到信号的功率谱。
缺点
需要预先设定模型函数的形式和参数，且计算复 杂度较高。
03 经典功率谱估计的优缺点
优点
01
02
03
算法成熟
经典功率谱估计方法经过 多年的研究和发展，已经 相当成熟，具有较高的稳 定性和可靠性。
9现代信号处理-功率谱估计(2)
rxx (1)
rxx(2)
rxx (0) rxx (1)
rxx
(N
1)
rxx
(N
2)
0
rxx
(
N
1)
rxx(N )
rxx (1)
可以看出AR模型得到的结果与按最大熵外推rxx(N+1)得到的结果 一致，这就证明了当x(n)为高斯分布时的最大熵谱估计与AR模型 法是等价的。
上式（4.6.8)是rxx(N+1)的一次函数，由此可解得rxx(N+1)。再 用类似的方法求得rxx(N+2)， rxx(N+3)，┄，然后确定功率谱估计。
最大熵谱估计用下式计算信号功率谱：
Pxx (e j )
2 w
N
1 ake jk
k 1
（4.6.9）
二、最大似然谱估计
１、最小方差谱估计
最大似然谱估计是用一个FIR滤波器实现，该滤波器对所
关心频率的正弦信号，可以无失真地通过，而对于其它频率的
信号，让其频响尽可能地小，亦即将它们尽可能地滤除。此时,
rxx (1)
rxx(2)
rxx (0) rxx (1)
rxx ( N
1)
rxx ( N
2) 0
（4.6.7）
rxx
(
N
1)
rxx(N )
rxx (1)
上 式 是 rxx(N+1) 的 一 次 函 数 ， 可 以 解 出 rxx(N+1) 。 继 续 再 将
rxx(N+1) 代 入 Rxx(N+2) 和 det(Rxx(N+2)) 中 ， 求 det(Rxx(N+2)) 对
现代信号处理大作业
现代信号处理大型作业一．试用奇阶互补法设计两带滤波器组(高、低通互补)，进而实现四带滤波器组；并画出其频响。
滤波器设计参数为：F p ＝1.7KHz , F r ＝2.3KHz , F s ＝8KHz , A rmin ≥70dB 。
（一）、分析与通常的滤波器相比，互补滤波器具有优良的结构特性和结构特性，具有较低的噪声能量和系数敏感性，其定义如下：一组滤波器H 12(),(),.......()Z H Z H Z n 如果满足下式：He Kjw k n(),==∑110<w<2π 则称这组滤波器为幅度互补滤波器;如果满足下式：He kjw k n()=∑=121, 0<w<2π则称这组滤波器为功率互补滤波器，同时互补滤波器还应该满足：Hz A z kk n()()=∑=1其中A(z)为全通函数，适当的选择全通函数，可以使两带函数具有所需要的低通和高通特性。
（二）、设计步骤(1) 对Fp 、Fr 进行预畸);();(''FsFrtg FsFptg r p ∏=Ω∏=Ω(2) 计算'''*r p c ΩΩ=Ω，判断'c Ω是否等于1，即该互补滤波器是否为互补镜像滤波器(3) 计算相关系数⎪⎩⎪⎨⎧-==+++=+-=-=ΩΩ=--=偶数)N 为(;21奇数)N 为 (;;lg /)16/1lg(;150152;1121;1;;])110)(110[(1213090500''02'''211－min1.0min1.0i i u q k N q q q q q k k q k k k k rp Ar Ap;)2cos()1(21))12(sin()1(21)1(21'2∑∑∞=∞=+-++-=Ωm mm m m m m i u Nm q u Nm q q ππ;42⎥⎦⎤⎢⎣⎡=N N;221N N N -⎥⎦⎤⎢⎣⎡=;)/1)(1(2'2'k k v i i i Ω-Ω-=12'1212,1;12N i v i i i =Ω+=--α 22'22,1;12N i v iii =Ω+=β (4) 互补镜像滤波器的数字实现;22i ii A αα+-=;22iii B ββ+-=1221,1;1)(N i ZA Z A Z H i i i =++=∏--22212,1;1)(N i ZB Z B Z Z H i i i =++=∏--- )];()([21)(21Z H Z H Z H L +=（三）、程序与结果 1． 二带滤波器组 （1） 源程序： clear; clf;Fp=1700;Fr=2300;Fs=8000; Wp=tan(pi*Fp/Fs); Wr=tan(pi*Fr/Fs); Wc=sqrt(Wp*Wr); k=Wp/Wr;k1=sqrt(sqrt(1-k^2)); q0=0.5*(1-k1)/(1+k1);q=q0+2*q0^5+15*q0^9+150*q0^13; N=11;N2=fix(N/4); M=fix(N/2); N1=M-N2; for jj=1:M a=0;for m=0:5a=a+(-1)^m*q^(m*(m+1))*sin((2*m+1)*pi*jj/N);%N is odd, u=j end ab=0;for m=1:5b=b+(-1)^m*q^(m^2)*cos(2*m*pi*jj/N); end bW(jj)=2*q^0.25*a/(1+2*b);V(jj)=sqrt((1-k*W(jj)^2)*(1-W(jj)^2/k)); endfor i=1:N1alpha(i)=2*V(2*i-1)/(1+W(2*i-1)^2); endfor i=1:N2beta(i)=2*V(2*i)/(1+W(2*i)^2); endfor i=1:N1a(i)=(1-alpha(i)*Wc+Wc^2)/(1+alpha(i)*Wc+Wc^2); endfor i=1:N2b(i)=(1-beta(i)*Wc+Wc^2)/(1+beta(i)*Wc+Wc^2); endw=0:0.0001:0.5;LP=zeros(size(w));HP=zeros(size(w));for n=1:length(w)z=exp(j*w(n)*2*pi);H1=1;for i=1:N1H1=H1*(a(i)+z^(-2))/(1+a(i)*z^(-2)) ;endH2=1/z;for i=1:N2H2=H2*(b(i)+z^(-2))/(1+b(i)*z^(-2));endLP(n)=abs((H1+H2)/2);HP(n)=abs((H1-H2)/2);endplot(w,LP,'b',w,HP,'r');hold on;xlabel('digital frequency');ylabel('amptitude');（2）运行结果：见图1图1 二带数字滤波器组2．四带滤波器组（1）源程序：clf;Fp=1700;Fr=2300;Fs=8000;Wp=tan(pi*Fp/Fs);Wr=tan(pi*Fr/Fs);Wc=sqrt(Wp*Wr);k=Wp/Wr;k1=sqrt(sqrt(1-k^2));q0=0.5*(1-k1)/(1+k1);q=q0+2*q0^5+15*q0^9+150*q0^13;N=11;N2=fix(N/4);M=fix(N/2);N1=M-N2;for jj=1:Ma=0;for m=0:5a=a+(-1)^m*q^(m*(m+1))*sin((2*m+1)*pi*jj/N); % N is odd, u=jendb=0;for m=1:5b=b+(-1)^m*q^(m^2)*cos(2*m*pi*jj/N);endW(jj)=2*q^0.25*a/(1+2*b);V(jj)=sqrt((1-k*W(jj)^2)*(1-W(jj)^2/k));Endfor i=1:N1alpha(i)=2*V(2*i-1)/(1+W(2*i-1)^2);endfor i=1:N2beta(i)=2*V(2*i)/(1+W(2*i)^2);endfor i=1:N1a(i)=(1-alpha(i)*Wc+Wc^2)/(1+alpha(i)*Wc+Wc^2);endfor i=1:N2b(i)=(1-beta(i)*Wc+Wc^2)/(1+beta(i)*Wc+Wc^2);endw=0:0.0001:0.5;LLP=zeros(size(w));LHP=zeros(size(w));HLP=zeros(size(w));HHP=zeros(size(w));for n=1:length(w)z=exp(j*w(n)*2*pi);H1=1;for i=1:N1H1=H1*(a(i)+z^(-2))/(1+a(i)*z^(-2)) ;endH21=1;for i=1:N1H21=H21*(a(i)+z^(-4))/(1+a(i)*z^(-4)) ;H2=1/z;for i=1:N2H2=H2*(b(i)+z^(-2))/(1+b(i)*z^(-2));endH22=1/(z^2);for i=1:N2H22=H22*(b(i)+z^(-4))/(1+b(i)*z^(-4));endLP=((H1+H2)/2);HP=((H1-H2)/2);LLP(n)=abs((H21+H22)/2*LP);LHP(n)=abs((H21-H22)/2*LP);HHP(n)=abs((H21+H22)/2*HP);HLP(n)=abs((H21-H22)/2*HP);endplot(w,LLP,'b',w,LHP,'r',w,HLP,'k',w,HHP,'m')hold onxlabel('digital frequency');ylabel('amptitude');（2）运行结果：见图2图2 四带数字滤波器组二、根据《现代数字信号处理》第四章提供的数据，试用如下方法估计其功率谱，并画出不同参数情况下的功率谱曲线：1）Levison算法2）Burg算法3） ARMA 模型法 4） MUSIC 算法 1 Levinson 算法Levinson 算法用于求解Yule-Walker 方程，是一种按阶次进行递推的算法，即首先以AR （0）和AR （1）模型参数作为初始条件，计算AR （2）模型参数；然后根据这些参数计算AR （3）参数，等等，一直到计算出AR （p ）模型参数为止，需要的运算量数量级为2p ，其中p 为AR 模型的阶数。
envi对雷达影像的滤波方法
雷达影像是一种利用雷达信号获取地表信息的遥感技术。
在雷达影像处理的过程中，滤波方法是非常重要的一环，它可以有效地提高雷达影像的质量，并且在后续的特征提取和目标识别过程中起到至关重要的作用。
在众多的滤波方法中，ENVi是一款常用的遥感图像处理软件，它提供了多种滤波方法，如中值滤波、高斯滤波、维纳滤波等。
本文将针对ENVi对雷达影像的滤波方法进行详细介绍和分析。
一、中值滤波中值滤波是一种常见的非线性滤波方法，它的原理是将像素点周围的邻域像素灰度值进行排序，然后取中值作为该像素点的灰度值。
在ENVi软件中，中值滤波可以有效地去除图像中的椒盐噪声和斑点噪声，同时保持图像的边缘信息。
该方法适用于雷达影像中含有较多噪声的情况，能够有效地提高图像的质量和清晰度。
二、高斯滤波高斯滤波是一种常用的线性滤波方法，它基于高斯函数对像素点周围邻域的像素进行加权平均，以达到平滑图像的目的。
在ENVi软件中，可以通过设置滑动窗口的大小和标准差来调节高斯滤波的效果。
该方法可以有效地去除雷达影像中的高频噪声，同时能够保持图像的整体灰度分布特性，使得图像更加平滑自然。
三、维纳滤波维纳滤波是一种基于统计模型的滤波方法，它通过估计图像中信号和噪声的功率谱密度来进行滤波处理。
在ENVi软件中，维纳滤波可以根据输入图像的特性自适应调整滤波参数，能够有效地消除雷达影像中的噪声，并且保持图像的边缘和细节信息。
该方法适用于要求较高图像质量和清晰度的场景，能够有效地提高雷达影像的识别和分析能力。
四、多通道滤波多通道滤波是一种基于多波段信息的滤波方法，它可以利用雷达影像中不同波段的特性进行组合和处理。
在ENVi软件中，可以通过组合不同波段的滤波结果来增强图像的特定信息，比如边缘信息、纹理信息等。
该方法适用于需要综合利用多波段信息的场景，能够有效地提高雷达影像的数据融合和分析能力。
ENVi提供了多种有效的雷达影像滤波方法，包括中值滤波、高斯滤波、维纳滤波和多通道滤波等。
matlab中计算功率谱的4种方法
在MATLAB中，计算功率谱是信号处理和频谱分析中的重要任务。
功率谱可以帮助我们了解信号中不同频率成分的能量分布情况，对于理解信号特性和进行频谱分析都是至关重要的。
在MATLAB中，有多种方法可以用来计算功率谱，在本文中，我将介绍并比较其中的四种常用方法。
第一种方法是使用MATLAB中的`periodogram`函数。
`periodogram`函数可以直接计算信号的功率谱密度（PSD），它采用傅里叶变换的方法，将信号从时域转换到频域，并计算功率谱密度。
这种方法简单直接，适用于对功率谱快速估计的情况。
在使用`periodogram`函数时，我们可以指定窗函数和重叠比例等参数，来对功率谱的估计进行优化。
第二种方法是使用`pwelch`函数。
`pwelch`函数也可以用来计算信号的功率谱密度，它采用Welch方法，通过对信号进行分段，然后对每个段进行傅里叶变换，并对结果进行平均来估计功率谱密度。
Welch 方法可以减小估计的方差，得到更平滑和可靠的功率谱估计结果。
在使用`pwelch`函数时，同样可以指定窗函数和重叠比例等参数来优化估计结果。
第三种方法是使用`fft`函数和自行计算功率谱。
通过对信号进行傅里叶变换得到频谱，然后对频谱的幅度进行平方运算，即可得到功率谱。
这种方法的好处是灵活性高，可以根据具体需求对傅里叶变换和求平方的结果进行后续处理，比如进行平滑或滤波操作。
但是需要注意的是，自行计算功率谱需要对信号处理和频谱分析有较深的理解。
第四种方法是使用`cpsd`函数。
`cpsd`函数可以用来计算信号之间的交叉功率谱密度，适用于多信号系统中不同信号之间的频谱分析。
交叉功率谱密度可以帮助我们理解不同信号之间频率成分的相关性和影响程度，对于系统建模和故障诊断都是非常有帮助的。
MATLAB提供了多种方法来计算功率谱，每种方法都有其适用的场景和优势。
在具体应用中，我们可以根据信号特性和分析需求来选择合适的方法。
功率谱分析和频域滤波技术
功率谱分析和频域滤波技术在功率谱分析中，首先需要将时域信号进行离散化处理，通常使用傅里叶变换（Fourier Transform，FT）或者快速傅里叶变换（FastFourier Transform，FFT）来得到信号的频谱。
频谱表示信号在不同频率下的能量分布情况，通常以功率谱密度（Power Spectral Density，PSD）表示。
功率谱分析可以用于分析各种信号，例如音乐、语音、图像等。
对于音乐和语音信号，可以通过功率谱分析得到它们的音调、音量和谐波分布等特征。
对于图像信号，功率谱分析可以用于检测图像中的噪声、边缘和纹理等特征。
频域滤波技术利用功率谱信息对信号进行滤波处理，主要包括低通滤波、高通滤波、带通滤波和带阻滤波等。
在频域滤波中，通常将信号的频谱与一个滤波器的频谱相乘，以实现相应频率成分的增强或抑制。
通过调整滤波器的频率响应和幅度响应，可以实现不同滤波效果。
以音频去噪为例，可以通过功率谱分析和频域滤波技术去除音频信号中的噪声。
首先将音频信号进行离散化，然后通过FFT得到其功率谱密度。
根据噪声在频率上的特性，可以设计一个适当的频域滤波器来抑制噪声的频率成分。
最后，将滤波后的功率谱进行逆傅里叶变换，得到滤波后的音频信号。
在图像处理中，频域滤波技术可以应用于图像增强和去噪。
通过功率谱分析和频域滤波，可以提取图像中的高频特征，如边缘和纹理，以实现图像增强。
同时，通过频域滤波可以对图像中的噪声进行抑制，例如通过低通滤波器去除图像中的高频噪声。
总之，功率谱分析和频域滤波技术是信号处理中常用的方法。
通过功率谱分析可以将信号从时域转换到频域，以便更好地理解信号的特性。
频域滤波技术则可以利用功率谱信息对信号进行滤波处理，以实现信号的去噪、增强或频率选择。
这些方法在音频、语音和图像处理等领域都有广泛应用，并且不断有新的研究和改进方法出现，为信号处理提供了更多的工具和技术。
《功率谱估计》课件
目录
• 引言 • 功率谱估计的基本原理 • 常见功率谱估计方法 • 现代功率谱估计方法 • 功率谱估计的性能评估 • 实际应用案例分析
01
引言
功率谱估计的定义
功率谱估计是对信号的频率内容进行描述的方法，通过分析信号在不同频率的功 率分布情况，可以了解信号的特性。
功率谱估计可以分为非参数方法和参数方法两类，其中非参数方法包括傅里叶变 换、Welch方法等，而参数方法则包括AR模型、MA模型、和ARMA模型等。
非参数模型
不假设信号的功率谱具有特定参数形式，而是直接从数据中估计功率谱。
03
常见功率谱估计方法
直接法
定义
直接法是通过测量信号的样本值，利用离散 傅里叶变换（DFT）直接计算信号的频谱。
特点
计算简单，但容易受到频率偏移和相位失真的影响 。
应用场景
适用于信号频率稳定且对相位精度要求不高 的场合。
间接法
THANKS
感谢观看
分辨率与假峰率
分辨率（Resolution）
衡量功率谱估计中能够区分两个相近频率成分的能力。分辨率越高，说明估计的功率谱能够更好地分 辨出相近的频率成分。
假峰率（False Peak Rate）
衡量估计的功率谱中出现的虚假频率峰的概率。假峰率越低，说明估计的功率谱中虚假频率峰的出现 概率越小。
06
特点
能够减小频谱泄漏效应，提高频 谱分辨率。
应用场景
适用于信号持续时间较短或需要 高分辨率频谱分析的场合。
最大熵法
定义
最大熵法是一种基于信息论的方法，通过最 大化熵函数来估计信号的功率谱。
特点
能够提供平滑且连续的功率谱估计，但计算 复杂度较高。
谱估计(现代)
ak xx (m k ) Ex(n) (n m)
k 1
p
而
m0 0, E x(n) (n m) 2 , m 0
•Yule-Walker方程的推导
故
p a k xx (m k ) , m 0 k 1 xx (m) p a (k ) 2 , m 0 k xx k 1 或
p
2
需要推导AR参数与 xx (m)之间的关系。
3.1
• 估计方法
自回归模型法
2 与xx (m)乊间的关系 参数a1, a2, a3, …, ap及 ——Yule-Walker方程
已知：自相关函数 已知： 自相关函数
Yule-Walker方程
要求： AR模型的阶数p，以及p个AR 要求： AR模型的阶数p，以及p个 AR 参数a(i),激励源方差 2 参数a(k),激励源方差
3.2
最大熵谱估计法
• 基本思想——熵
代表一种不定度； 最大熵为最大不定度，即它的时间序列最随机， 它的PSD应是最平伏(最白色)。 Shannon对熵的定义： 当x的取值为离散的时，熵H定义为
H pi ln pi
i
pi：出现状态i 的概率。
当x的取值为连续的时，熵H定义为
p(x)：概率密度 函数
(n)
...
z-1 a1
z-1
z-1
a2
...
ap
3.1
自回归模型法
q
• MA(Moving Average)模型 ——全零点模型
x(n) bl (n l )
l 0
H ( z ) B( z ) 1 bl z k
高级计量分析(时间序列分解——趋势分解)
时间序列分解——趋势分解一、研究目的在季节调整案例中，介绍了如何对经济时间序列进行分解，但在季节调整方法中，趋势和循环要素视为一体不能分开。
本案例专门讨论如何将趋势和循环要素进行分解的方法。
测定长期趋势有多种方法，比较常用的方法有回归分析法、移动平均法、阶段平均法、HP （Hodrick-Prescott ）滤波算法和频谱滤波算法（frequency(bandpass)filter,BP 滤波）。
本案例主要介绍HP 滤波算法和BP 频谱滤波算法。
二、滤波算法的原理 1、HP 滤波算法设t Y 是包含趋势成分和波动成分的经济时间序列，T t Y 是趋势成分，c t Y 是波动成分。
则：,1,2,,T c t t t Y Y Y t T =+=HP 滤波算法就是从t Y 中将T t Y 分解出来。
一般的，时间序列t Y 中的可观测部分趋势Tt Y 常被定义为下面最小化问题的解：{}221min ()()TT T t t t t Y Y c L Y λ=⎡⎤-+⎣⎦∑ （1）其中，()c L 是延迟算子多项式：1()(1)(1)c L L L -=--- （2）将（2）式代入（1）式，则HP 滤波算法就是使得下面的损失函数最小，即：()()221111min ()T T T T T T Tt t t t t t t t Y Y Y Y Y Y λ+-==⎧⎫⎡⎤-+---⎨⎬⎣⎦⎩⎭∑∑ 最小化问题用2()Tt c L Y ⎡⎤⎣⎦来调整趋势的变化，并随着λ的增大而增大。
HP 滤波依赖于参数λ，该参数需要给定。
这里存在一个权衡问题，要在趋势要素对实际序列的跟踪程度和趋势光滑程度之间做一个选择。
当0λ=时，满足最小化问题的趋势序列与原序列相同；随着λ值的增加，估计的趋势越光滑；当λ趋于无穷大时，估计的趋势接近一条直线。
一般经验，λ的取值如下：100,1600,14400,λ⎧⎪=⎨⎪⎩年度数据季度数据月度数据2、BP 频谱滤波算法由于该方法的数学方法（傅立叶变换）较为复杂，这里我们只介绍其基本思想：该方法把时间序列看成是不同谐波的叠加，研究时间序列在频率域里的结构特征，由于这种分析主要是用功率谱的概念进行讨论，所以通常也称为谱分析。
welch谱估计的随机误差与置信度
welch谱估计的随机误差与置信度Welch谱估计法（Welch'smethod），也称为自由度不相等的t检验法，是一种常用的频谱估计方法，通过对信号进行分段并对每个段进行功率谱估计，最终得到整个信号的功率谱估计结果。
本文将探讨Welch谱估计法的随机误差和置信度。
一、Welch谱估计法的基本原理Welch谱估计法是一种常用的频谱估计方法，其基本原理是将信号分成多个重叠的段，对每个段进行功率谱估计，最终将估计结果加权平均得到整个信号的功率谱估计结果。
其基本步骤如下：1. 将信号分成M个长度为L的段，每个段之间有重叠部分。
2. 对每个段进行窗函数处理，以减少频谱泄漏。
3. 对每个段进行功率谱估计，得到M个功率谱估计结果。
4. 对M个功率谱估计结果进行加权平均，得到整个信号的功率谱估计结果。
二、Welch谱估计法的随机误差Welch谱估计法的随机误差可以分为两部分：第一部分是由于信号本身的随机性导致的误差，第二部分是由于功率谱估计方法本身的随机性导致的误差。
1. 信号本身的随机性导致的误差信号本身的随机性是指信号在时间和频率上的随机性。
由于Welch谱估计法将信号分成多个段进行功率谱估计，每个段的长度和位置都是随机的，因此信号本身的随机性会导致每个段的功率谱估计结果不同，最终整个信号的功率谱估计结果也会受到影响。
2. 功率谱估计方法本身的随机性导致的误差功率谱估计方法本身的随机性是指在进行功率谱估计时，采样数目有限，因此估计结果会受到采样噪声的影响。
同时，采用不同的窗函数和参数也会影响功率谱估计结果的随机性。
三、Welch谱估计法的置信度Welch谱估计法的置信度可以通过置信区间来衡量。
置信区间是指在给定置信水平下，估计结果的真实值有一定概率落在该区间内。
通常情况下，置信水平取为95%。
置信区间的计算方法如下：1. 根据置信水平和自由度计算t分布的上下界。
2. 根据t分布的上下界和估计结果的标准误差计算置信区间。
功率谱估计及其MATLAB仿真
功率谱估计及其MATLAB仿真一、本文概述功率谱估计是一种重要的信号处理技术，它能够从非平稳信号中提取有用的信息，揭示信号在不同频率上的能量分布特征。
在通信、雷达、生物医学工程、地震分析等领域，功率谱估计都发挥着至关重要的作用。
随着计算机技术的快速发展，功率谱估计的仿真研究也越来越受到重视。
本文将对功率谱估计的基本理论进行简要介绍，包括功率谱的概念、性质以及常见的功率谱估计方法。
随后，我们将重点探讨MATLAB 在功率谱估计仿真中的应用。
MATLAB作为一种功能强大的数值计算和仿真软件，为功率谱估计的研究提供了便捷的工具。
通过MATLAB，我们可以轻松地模拟出各种信号，进行功率谱估计，并可视化结果，从而更直观地理解功率谱估计的原理和方法。
本文旨在为读者提供一个关于功率谱估计及其MATLAB仿真的全面而深入的学习机会，帮助读者更好地掌握功率谱估计的基本原理和仿真技术，为后续的实际应用打下坚实的基础。
我们将通过理论分析和实例仿真相结合的方式，逐步引导读者深入了解功率谱估计的奥秘，探索MATLAB在信号处理领域的广泛应用。
二、功率谱估计的基本原理功率谱估计是一种在信号处理领域中广泛使用的技术，它旨在从时间序列中提取信号的频率特性。
其基本原理基于傅里叶变换，通过将时域信号转换为频域信号，可以揭示信号中不同频率分量的存在和强度。
功率谱估计主要依赖于两个基本概念：自相关函数和功率谱密度。
自相关函数描述了信号在不同时间点的相似程度，而功率谱密度则提供了信号在不同频率下的功率分布信息。
在实际应用中，由于信号往往受到噪声的干扰，直接计算功率谱可能会得到不准确的结果。
因此，功率谱估计通常使用窗函数或滤波器来减小噪声的影响。
窗函数法通过在时域内对信号进行分段，并对每段进行傅里叶变换，从而减小了噪声对功率谱估计的干扰。
而滤波器法则通过在频域内对信号进行滤波，去除噪声分量，得到更准确的功率谱。
MATLAB作为一种强大的数值计算和仿真软件，为功率谱估计提供了丰富的函数和工具。
第五章功率谱估计1-2节
经FFT变换,得:
ˆ ˆ ˆ Pxx (k ) FFT xx (m) xx (m)e
m0 L -1 -j 2 km L
k 0,1, 2, L -1
29/113
三、相关图法功率谱估计质量
用x(n)的N 个有限值得到 ˆ 自相关函数的估计 ( m)，
13/113
（a）间接法（BT法）
BT法又称为相关图法 对信号序列估计求其自相关函数值 对自相关函数的估计进行加权 对加权的自相关函数做傅里叶变换 获得功率谱估计。
直到1965年快速傅里叶变换算法（FFT） 问世以前，是最流行的谱估计方法。
14/113
（b）直接法（又称周期图 （periodogram)法）
对观测到的数据样本直接进行傅里叶变换 取模的平方，再除以N 得到功率谱估计。 不用估计自相关函数，且可以用FFT进行计算， 在FFT出现以后，周期图法才得到了广泛的应 用。
15/113
（2）现代谱估计
其基本思想是根据已有的观测数据，建 立信号所服从的模型，从而在观测不到 的区间上，信号的取值服从模型的分布 情况，不再认为是零。 主要讨论参数模型（AR、MA、ARMA） 法。
N
2 xx (l ) xx (l m)xx (l - m) (N - m - l )
N - m -1 2 l -( N - m -1) N - m -1 2 l -( N - m -1)
N - m
N
所以在实际中必须兼顾分辨率与方差的要求来适当选择信号仍然是均值为方差为的白噪声观察数据长度为了利用平均周期法估计其功率谱将它分成段分别按照平均周期图法估计其功率谱得到功率谱曲线如图从图中可以看出随着分段数的增加功率谱估计值在附近的幅度愈来愈小显示出分段平均对周期图方差减少有明显效果
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
清华大学《数字信号处理》期末作业2013 年 1 月第一题掌握去噪的方法1.1 题目描述MATLAB 中的数据文件noisdopp 含有噪声，该数据的抽样频率未知。
调出该数据，用你学过的滤波方法和奇异值分解的方法对其去噪。
要求：1．尽可能多地去除噪声，而又不损害原信号；2．给出你去噪的原理与方法；给出说明去噪效果的方法或指标；3．形成报告时应包含上述内容及必要的图形，并附上原程序。
1.2 信号特性分析MATLAB所给noisdopp信号极其频域特征如图1.1、图1.2。
图1.1含有噪声的noisdopp信号图1.2 noisdopp 信号频域特性其中横坐标f 采用归一化频率，即未知抽样频率Fs 对应2（与滤波器设计时参数一致）。
信号基本特性是一个幅值和频率逐渐增加的正弦信号叠加噪声，噪声为均匀的近似白噪声，没有周期等特点。
因为噪声信号能量在全频带均匀分布，滤波器截止频率过低则信号损失大，过高则噪声抑制小，认为频谱中含有毛刺较多的部分即为信噪比较小的部分，滤除这部分可以达到较好的滤波效果。
先给定去噪效果的评定指标。
信号开始阶段频率较高（如图1.3，红圈为信号值），一周期内采样点4~5个，即信号归一化频率达到0.4~0.5（Fs=2），难以从频域将这部分信号同噪声分离，滤波后信号损失较大，故对前128点用信噪比考察其滤波效果，定义：22()10lg (()())k k x k SNR y k x k =-∑∑其中，()x k 为原nosidopp 信号，()y k 为滤波后信号。
SNR 越大表示滤除部分能力越小，可以反映滤波后信号对原信号的跟踪能力，对前128点主要考察SNR ，越大滤波器性能越好。
对128点以后的点，信号能量集中在低频部分，滤波后效果不仅仅需考察SNR，同时考察滤波后信号的平滑程度。
定义平滑指标[1]22((1)())((1)())kky k y krx k x k+-=+-∑∑r越小表示滤波后信号相比原信号更加平滑，滤波效果更好。
图1.3信号前100点本文采用整体滤波、分时段滤波和奇异值分解的方法分别对信号进行去噪，并采用信噪比和平滑指数两个指标评定去噪效果。
1.3 滤波器设计根据信号的频谱，设计IIR滤波器和FIR滤波器对信号去噪，用重叠相加法对信号进行卷积运算（每128点截为16段）。
根据信号的频谱，均采用低通滤波方法去噪，将滤波器截止频率设为0.2。
1.3.1 IIR 滤波器滤波设计巴特沃斯低通滤波器，给定参数0.17,0.23,0.1,60p s p s R dB R dB ωω====，由buttord 得到阶数7N =。
用MATLAB 函数butter 设计滤波器，并得到70点频率响应如图1.4，滤波器冲击响应如图1.5。
图1.4 巴特沃斯低通滤波器幅频特性图1.5 巴特沃斯低通滤波器冲击响应采用重叠相加法计算nosidopp信号于巴特沃斯低通滤波器卷积结果。
nosidopp信号每128点分为一段，共16段，得到滤波后信号有1093点。
由图1.5可以得到，设计得到的巴特沃斯滤波器存在30点的延时，取卷积后信号31点开始的1024点为滤波后信号y。
得到滤波后信号如图1.6。
与原信号比较如图1.7。
图1.6滤波后信号图1.7滤波后信号与原信号直观看滤波后信号去除了大部分噪声，但开始高频部分信号损失较多（如图1.7），后尾部信号频率较低，滤波后噪声含量仍比较多。
图1.8信号起始部分前128点信噪比SNR为3.429，后896点信号SNR为12.715，全信号平滑指标为0.0642。
本题中设计的IIR滤波器的问题在于①不是线性相位，只能近似得到滤波后信号的延时，滤波后信号难以同原信号“对齐”②冲击响应无限长，有限长序列通过IIR系统得到的输出也是无限长，存在截断误差③分时段设计滤波器时计算较慢。
1.3.2 窗函数法设计FIR滤波器滤波选择性能较优的汉宁窗设计低通滤波器，截止频率0.2，N取34，得到35点长度的滤波器。
得到滤波器幅频响应和冲击响应如图1.9，图1.10图1.9低通滤波器幅频特性图1.10低通滤波器冲击响应该低通滤波器具有很好的线性相位，输出信号存在N/2=17点的延时，得到的滤波后信号如图1.11。
与原信号比较如图1.12。
图1.11滤波后信号图1.12滤波后信号与原信号直观感觉比巴特沃斯滤波器的滤波结果在信号尾部纹波减小，信号开始部分损失也比较大。
得到前128点信噪比SNR为6.448，后896点信号SNR为14.418，全信号平滑指标为0.0643。
由计算得到的SNR和平滑指标看，窗函数法得到的滤波器滤波效果比巴特沃斯滤波器要好。
采用矩形窗设计同样阶数的FIR滤波器，得到的滤波效果为：前128点信噪比SNR为10.303，后896点信号SNR为15.738，全信号平滑指标为0.184，如图1.13。
图1.13滤波后信号与原信号显然，全信号采用同样的滤波器参数滤波，信号高频分量保持较好，则会保留较多的噪声。
若在频域对信号直接加窗，相当于N=∞的窗函数滤波器，恢复信号有明显的吉布斯现象，如图1.14。
图1.14滤波后信号的吉布斯现象1.3.3 一致逼近法设计FIR滤波器滤波选择性切比雪夫一致逼近法设计低通滤波器，截止频率0.2，N取34，通带、阻带加权均取1，得到35点长度的滤波器。
滤波器幅频响应如图1.15。
图1.15低通滤波器幅频特性具有很好的线性相位，同样，输出信号存在N/2=17点的延时，得到的滤波后信号如图1.16。
与原信号比较如图1.17。
图1.16滤波后信号图1.17滤波后信号与原信号前128点信噪比SNR 为6.167，后896点信号SNR 为13.987，全信号平滑指标为0.0655，与窗函数法得到的滤波效果类似。
对全信号采用同样的滤波器参数滤波，因为没有考虑信号频率随时间的变化，不能得到更好的滤波效果。
1.4 分时段滤波用短时傅里叶变化（128,(128),127,2nfft window hanning overlap Fs ====）估计信号的时频特性大致如图1.18。
图1.18 Nosidopp 信号的时频特性同样将信号分为16段，每段64点，对每段信号做离散傅里叶变化求其频谱，得到幅值最大点频率为max f ，每两点间2/640.0313f ∆==，因为信号的包络为低频分量，需要保留，故设计max 3f f +∆为截止频率的低通滤波器。
采用切比雪夫一致逼近法，N 取34，通带、阻带加权均取1，输出信号存在N/2=17点的延时。
得到的滤波后信号如图1.19。
与原信号比较如图1.20。
图1.19滤波后信号图1.20滤波后信号与原信号图1.21滤波后信号与原信号前100点前128点信噪比SNR为8.144，后896点信号SNR为13.734，全信号平滑指标为0.0773，得到了兼顾信号高频分量和去噪的效果，滤波效果比不分时段要好。
1.5 奇异值分解方法滤波奇异值分解滤波的基本原理是通过矩阵分解方法，得到信号的特征值，特征值中较大的k个分量反应了信号低频成分，而较小的分量反应高频噪声，将较小的特征值置零再恢复信号，得到去噪的信号。
奇异值分解方法滤波需要考虑两个问题①一维信号构成用于奇异值分解的矩阵的方法②如何去除特征值。
本文采用教材式（9.5.18）给定的方式构造矩阵，对长度N的信号x(n)，定义矩阵X1为：1(0)(1)(1)(1)(2)()(1)()(1)x x x Mx x x Mx L x L x N-⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥--⎣⎦XLLM M M ML得到的矩阵为方阵可最大程度利用其特征值，对nosidopp信号，N=1024，故选择L=513，M=512构造矩阵。
利用svd分解得到特征值，如图1.21。
同时得到矩阵U，V。
图1.22 矩阵特征值可以看出特征值存在明显的转折点，将转折点（第30点）以后的特征值值为0。
恢复信号矩阵：1112111111ˆˆˆ(0)(1)(1)ˆˆˆ(1)(2)()ˆˆˆˆ(1)()(1)T x x x M x x x M x L x L x N -⎡⎤⎢⎥⎢⎥==⎢⎥⎢⎥--⎣⎦X USV L L M M M M L 此时得到的矩阵1ˆX 将不再对称，按上式中对应位置得到滤波后信号ˆ()y x n =为：1ˆˆ()(),0,1,...,1k kxi x i i N k ==-∑ 得到的滤波后信号如图1.22。
与原信号比较如图1.23。
图1.23 滤波后信号图1.24滤波后信号与原信号前128点信噪比SNR为6.145，后896点信号SNR为13.827，全信号平滑指标为0.0708。
直观观察滤波后信号质量好与滤波器滤波后的信号。
特别对信号前100点，如图1.24，直观观察要好于同样信噪比下滤波器去噪的效果。
图1.25滤波后信号与原信号前100点若采用将特征值中小于均值的数置零的方式滤波，得到滤波后信号如图1.26。
图1.26滤波后信号前128点信噪比SNR为10.83，后896点信号SNR为15.357，全信号平滑指标为0.315。
对噪声去除减少，对原信号特征保持较好。
构造矩阵时，减小矩阵的行数、增加矩阵列数，或采用文献[2]中提出的构造矩阵的方法：1(0)(1)(1)(1)(2)()(1)(0)(2)x x x Mx x x Mx N x x N-⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥--⎣⎦XLLM M M ML不能使得到的特征值拐点给清晰。
1.6 小结在以信噪比SNR和输出信号平滑指标两者作为滤波效果评定标准下，采用切比雪夫一致逼近法分时段设计低通滤波器得到了最好的滤波效果，最优指标为：前128点信噪比SNR为8.144，后896点信号SNR为13.734，全信号平滑指标为0.0773。
直观观察奇异值分解去噪后滤波效果最佳。
本文去噪方式存在的不足有①给定的评定滤波效果的标准不够合理，缺乏标准信号的情况下，信噪比作为标准难以区分估计误差和噪声，给定的标准不能完全反应直观观察出的滤波效果②缺乏对于信号所分段数L 的研究③奇异值分解方法中，对矩阵构造方式对滤波效果的研究不足。
第一题滤波器去噪部分MATLAB 程序见附录I ，奇异值分解去噪部分MATLAB 程序见附录II 。
第二题 掌握功率谱估计的方法2.1 试验数据的产生2.1.1 产生x(n)第一步，产生均值为0，功率为0.12的白噪声序列12(),()u n u n 。
利用教材所给出的式(1.10.2)，可以得到给定功率的白噪声序列，即rand 函数产生的500点均匀分布的白噪声序列，再乘以常数12*12*0.12=1.2a P ==得到12(),()u n u n 。
第二步，通过MATLAB 函数sos2tf 得到5个FIR 子系统级联组成的FIR 系统传递函数的系数firb 和fira ，firb 即为FIR 系统的冲击响应h(n)。