随机信号分析仿真

随机信号分析原理大作业报告专业： 水声工程 姓名： xxx 学号：xxxxxxxxxx20题目要求： 给定一个白噪声信号，它的均值和方差自定。1. 设计一个线性滤波器，使该滤波器的输出为一个窄带信号。并给出该窄带信号在不同的3个典型中心频率和带宽时的波形。2. 对该滤波器输出的上述窄带信号，用莱斯表示法对其进行建模，画出和的波形。3. 计算上述3种窄带信号对应的瞬时频率和瞬时相位，并进行包络检测。1 窄带信号的生成1.1 高斯白噪声的产生若为一个具有零均值的平稳随机过程，其功率谱密度均匀分布在（-，+）整个频率区间，即 (1)其中，为一个正实常数，则称为白噪声。白噪声的自相关函数为 (2)白噪声在任意两个相邻时刻（不管这两个时刻多么的近）的取值都是不相关的，这意味着白噪声过程随时间的起伏很快，过程的功率谱极宽。这种形式定义的白噪声只是一种理想化的模型，实际上这种白噪声是不存在的，因为按照定义，白噪声的均方值为无限大，而物理上存在的随机过程，其均方值总是有限的。在实际工作中，当所研究的随机过程通过某一系统时，只要过程的功率谱密度在一个比系统宽大的多的频率范围内近似均匀分布，就可以把它作为白噪声来处理。现产生一均值为0，方差为3高斯白噪声，如图1所示：图1 白噪声时域波形图1.2 设计线性滤波器为了产生一个窄带信号，让白噪声通过一理想带通线性系统。设理想带通系统的频幅特性为 (3)若输入白噪声的物理谱,则输出的物理谱为 (4)输出相关函数为 (5)式中若，即(3)式所示的系统的中心频率远大于系统的带宽，则称这样的系统为窄带系统。这时，(4)式所表示的随机信号的功率谱分布在高频周围一个很窄的频域内。设计线性滤波器为窄带滤波器，使高斯白噪声通过其变成窄带信号，此线性滤波器可选择带通滤波器，采样频率fs=25000Hz。分别设置三种不同中心频率和带通宽度，如下：滤波器1：中心频率f0=1000Hz，带通：9001100 Hz。滤波器1的幅度相位图以及滤波器1输出信号波形如下：图2 滤波器1的幅度相位图图3滤波器1输出信号的时域波形图4 滤波器1输出信号的频域波形滤波器2：中心频率f0=700Hz，带通：650750 Hz。滤波器2的幅度相位图以及滤波器2输出信号波形如下：图5 滤波器2的幅度相位图图6 滤波器2输出信号的时域波形图7 滤波器2输出信号的频域波形滤波器3：中心频率f0=400Hz，带通：375425 Hz。滤波器3的幅度相位图以及滤波器3输出信号波形如下：图8 滤波器3的幅度相位图图9滤波器3输出信号的时域波形图10 滤波器3输出信号的频域波形2 窄带随机过程的莱斯（Rice）表示任何一个实平稳随即过程都可以表示为 (5)式中，是固定值，对于窄带随机过程来说，一般取窄带滤波器的中心频率。和是另外两个随机过程。 (6) (7)为的希尔伯特变换。对于本文的窄带信号，在经过不同滤波器后和波形分别如下为：图11 滤波器1输出信号和波形图11 滤波器2输出信号和波形图13 滤波器3输出信号和波形3 窄带随机信号瞬时频率、瞬时相位和包络检测窄带随机过程可以表示为 (8)式中是窄带随机过程的中心频率或称载波频率。其中， (9) (10)瞬时相位为： (11)瞬时角频率为： (12)由希尔伯特变换得： (13)则包络幅值为： (14)瞬时相位为：瞬时角频率为：三个滤波器输出信号瞬时频率和瞬时相位，以及窄带信号的包络检测如图14图22所示。图14 滤波器1输出信号瞬时频率图15 滤波器1输出信号瞬时相位图16 滤波器1输出信号包络检测图17 滤波器2输出信号瞬时频率图18 滤波器2输出信号瞬时相位图19 滤波器2输出信号包络检测图20 滤波器2输出信号瞬时频率图21 滤波器2输出信号相位图22 滤波器3输出信号包络检测 本文相关仿真所使用的软件为matlab，以滤波器1为例，仿真程序见附件一。附件一 滤波器1输出信号仿真程序clear allclose allclc%产生高斯白噪声N=25000; %序列长度my_var = 2;noise = sqrt(my_var)*randn(1,N);%均值为0，方差为2figure(1)plot(noise)title(均值为0方差为2的高斯白噪声)grid on fs = 25000;%采样频率f0 = 1000;%中心频率%滤波器f_pass = 900 1100;omega_pass = 2*f_pass/fs;b = fir1(192,omega_pass);figure(2)freqz(b,1,1024)%滤波器幅度和相位图像grid on%噪声通过窄带滤波器filter_outpu = filter(b,1,noise);figure(3)plot(filter_outpu)title(窄带信号在时域的波形)grid on%做fft变换Nfft = fs;fft_x = fft(filter_outpu,Nfft);ff = 0:fs/Nfft:fs-fs/Nfft;figure(4)plot(ff,20*log10(abs(fft_x)%窄带信号的频谱title(窄带信号的频谱)xlabel(频率 Hz)ylabel(幅度 dB)grid on%窄带信号在时域的波形X_t = filter_outpu;t = 0:1/fs:1-1/fs;figure(5)plot(t,X_t)title(窄带信号在时域的波形)xlabel(t / s)grid on%莱斯表示法h_X = hilbert(X_t,Nfft) ;%希尔伯特变换omega0 = 2*pi*f0;A_t = X_t.*cos(omega0*t)+h_X.*sin(omega0*t);B_t = -1*X_t.*sin(omega0*t)+h_X.*cos(omega0*t);figure(6)subplot(2,1,1);plot(t,A_t)grid onhold onsubplot(2,1,2);plot(t,B_t)grid on%瞬时频率 瞬时相位 theta_t = atan(h_X./X_t);xh1=unwrap(angle(h_X);omega_t=fs*diff(xh1)/(2*pi);figure(7)plot(omega_t);title(瞬时频率)