matlab加特定频率的噪声,如何在Matlab/Simulink中生成指定的白噪声和有色噪声
在进行系统仿真时,经常需要利用白噪声和有色噪声作为系统输入,在Matlab和Simulink环境下提供了多种方式生成白噪声或者有色噪声。在这里对自己常用的一些方法进行简单总结。欢迎一起交流并提供更多的思路。
白噪声和有色噪声基础知识
在Matlab/Simulink中生成单位功率谱密度的白噪声
在Matlab/Simulink中利用成形滤波器生成有色噪声
功率谱密度分析
1. 白噪声和有色噪声基础知识
在学术声,白噪声和有色噪声的定义如下:
白噪声(white noise)是指功率谱密度在整个频域内是常数的噪声。 所有频率具有相同能量密度的随机噪声称为白噪声。
有色噪声( coloured noise)是指功率谱密度函数不平坦的噪声。大多数的噪声的频谱主要都是非白色频谱,通过信道的白噪声受信道频率的影响而变为有色的。
在物理世界中,任何系统都是有限带宽的。而根据白噪声的定义,其带宽是无穷大的,意味着其能量是无限的,这显然是不现实的。所以,在数学分析或者仿真时,我们常认为在有限带宽内具有平坦功率谱密度的信号就是白噪声。
式中,
为
的脉冲响应函数,则输入信号
和输出信号
的功率谱密度具有如下对应关系:
式中,
为传递函数
的频率响应。因此,可以通过将白噪声信号
通过线性传递函数
,则可以得到有色噪声信号
,如下所示:
为
的复共轭。
就是成形滤波器的传递函数形式。
下图给出了白噪声和有色噪声的功率谱密度对比图:
下面结合实际代码,介绍如何在Matlab/Simulink中分别生成指定功率谱密度的白噪声和有色噪声。
2. 在Matlab/Simulink中生成单位功率谱密度的白噪声
功率谱密度与频谱图一样,关于
对称,即
。因此对于频率
的功率谱密度,称之为双边功率谱密度(double-side power spectrum density);对于频率
则称为单边功率谱密度密度(single-side power spectrum density)。由双边功率谱密度转换为单边功率谱密度时需要乘以系数2:
下面分别介绍两种功率谱密度在Matlab/Simulink中的生成方式。
2.1 双边功率谱密度为1的白噪声
在Matlab中可以用随机数函数生成白噪声,例如randn()函数。根据白噪声的定义,均值为0,需要确定的只剩随机序列的方差。在前述理论基础中已经介绍了,在工程上,我们假设在有限带宽上具有常值功率谱密度的信号就是白噪声,因此有:
式中,
表示双边功率谱密度;
表示单边功率谱密度;
表示系统带宽。
若系统的采样频率$f_s=10\text{ Hz}$,根据奈奎斯特采样定理,系统最高频率或带宽则为$F = f_s/2=5 \text{ Hz}$。所以对于双边功率谱密度为$S_{xx}(f)=1$的信号,对应的噪声序列的方差为$\sigma^2 = 2\times5\times1=10$。对应的Matlab代码如下:
noise_power=1;
fs=10;
ts=1/fs;
t1_m = 0:1/fs:100;
x1_m = sqrt(2*fs/2*noise_power)*randn(length(t1_m),1); %fs/2 is Nyquist frequency
在Simulink中则提供了更加方便的方法。采用有限带宽白噪声模块(Band-limited White Noise)可以很方便的生成指定功率谱密度的白噪声。该模块需要设置的参数和定义如下:
Noise power:白噪声功率谱密度的幅值
Sample time:采样时间
Seed:随机数种子
需要注意的是该模块的中的功率谱密度的定义是双边功率谱密度。因此可以在该模块中赋如下值:
两种方法生成的生成的白噪声的仿真结果对比如下,可以看出所得到的结果是一致的。
利用mean()函数计算功率谱密度的平均值,得到的结果分别为:
和
。可见仿真结果与预期是一致的。
2.2 单边功率谱密度为1的白噪声
有了前述分析,单边功率谱密度为1的白噪声就可以类似得到了。
首先,在Matlab中,唯一要确定的就是噪声序列的方差。同样,若系统的采样频率
,根据奈奎斯特采样定理,系统最高频率或带宽则为
。所以对于单边功率谱密度为
的信号,对应的噪声序列的方差为
。对应的Matlab代码如下:
noise_power=1;
fs=10;
ts=1/fs;
t1_m = 0:1/fs:100;
x1_m = sqrt(fs/2*(noise_power))*randn(length(t1_m),1); %fs/2 is Nyquist frequency
在Simulink中,需要注意有限带宽白噪声模块的noise power值。由于在Simulink中,该值的定义为双边功率谱密度,根据双边功率谱密度和单边功率谱密度的关系,noise power需要设定为0.5,如下所示:
两种方法生成的生成的白噪声的仿真结果对比如下,可以看出所得到的结果是一致的。
利用mean()函数计算功率谱密度的平均值,得到的结果分别为:
和
。可见仿真结果与预期是一致的。
3. 在Matlab/Simulink中利用成形滤波器生成有色噪声
根据基础知识,将单位功率谱密度的白噪声通过成形滤波器可以得到有色噪声。成型滤波器可以通过传递函数建模得到,白噪声的建立已经分析过了,因此剩余的工作就十分简单了。
利用2.1中生成的双边功率谱密度为1的白噪声,在Matlab中定义传递函数如下所示的成形滤波器。
可以看出该成形滤波器也可以看作是一个低通滤波器。再利用lsim()函数则可以仿真得到有色噪声了,代码如下:
noise_power=1;
fs=10;
ts=1/fs;
t1_m = 0:1/fs:10000;
x1_m = sqrt(2*fs/2*noise_power)*randn(length(t1_m),1); %fs/2 is Nyquist frequency
s = tf('s');
shapeFilter = (1e0*s+1)^2/(1e2*s+1)^2;
y = lsim(shapeFilter,x1_m,t1_m);
在simulink中可以更加直观的进行仿真。搭建如下模块,则可以进行仿真得到有色噪声:
有色噪声仿真结构-水印.png
对比仿真结果如下所示。可以看出,两种方法得到的结果一致。仿真得到的有色噪声也与成形滤波器符合,证明了仿真的准确性。同时,中间的图也表明输入白噪声的高频变化被成形滤波器的低通滤波特性所滤去了。
4. 功率谱密度分析
功率谱密度是随机过程中无法绕开的概念。关于功率谱密度的相关介绍可以参考这篇文章--如何理解随机振动的功率谱密度? 。这篇文章非常形象生动,我就是从这篇文章开始理解了什么是功率谱密度,在此也表示对作者深深地感谢。
在 Matlab中(simulink中也有相应模块,可以自行搜索power spectrum density)提供了三种计算功率谱密度的方法:
利用快速傅立叶变换,在帮助文件中搜索Power Spectral Density Estimates Using FFT,可以学习相关算法
pwelch()函数
periodogram()函数
本次所有的功率谱密度计算均是采用periodogram()函数。下一篇文章计划结合实际例子和代码,总结功率谱密度的相关概念和知识,欢迎关注。
matlab加特定频率的噪声,如何在Matlab/Simulink中生成指定的白噪声和有色噪声相关推荐
- 如何在Angular 10中生成QR码
In this tutorial, we'll learn how to generate QR codes in Angular 10 by building a simple example ap ...
- 白噪声,有色噪声的定义、特性及其MATLAB仿真
一.白噪声 白噪声(white noise)是指功率谱密度在整个频域内是常数的噪声. 所有频率具有相同能量密度的随机噪声称为白噪声.白噪声是指在较宽的频率范围内,各等带宽的频带所含的噪声功率谱密度相等 ...
- 画出序列的图形matlab,江恩时间序列怎么画,如何在matlab上绘制基于时间序列的图形...
Q1:如何在matlab上绘制基于时间序列的图形 ..flag.. Q2:怎样用spss软件画出时间序列图 第一步:定义时间.步骤:数据-定义日期.有许多种日期模式,依实际情况定. 第二步:创建模型. ...
- MATLAB计算杨氏模量,四阶弹性模量Cijkl如何在matlab里表示啊? - 计算模拟 - 小木虫 - 学术 科研 互动社区...
matlab 四元数运算计算包就可以了吧 Matlab 四元数操作函数 2012-06-03 21:02:55| 分类: MATLAB&Mathemati | 标签:四元数 quater ...
- 如何在ASP.NET中生成HTML5离线Web应用
传统的Web应用程序有一个很大的症结是当用户的网络连接不好时,应用会加载失败,为了 解决这一问题,HTML5中引入了Web的离线工作的功能.离线功能使得Web应用程序类似于本机应用程序,当断开网络连接 ...
- java aes密钥生成_如何在Java(Android)中生成与.Net中相同的AES密钥?
我需要从.Net WebService提供的salt和密码生成 Java( Android)中的AES密钥.我需要使用与.net生成的密钥相同的密钥和相同的密码和盐(使用Rfc2898DeriveBy ...
- 使用matlab播放特定频率的声音
A=2;%振幅 f_0=397%声音频率 fs=10000; %采样频率 N=3000; % 信号样点数,播放时长 y=A*sin(2*pi*f_0*(0:N-1)/fs); %单频信号 sound( ...
- 在access窗体中加图片_如何在Access窗体中显示指定路径的图片
在Access中,如果把图形对象以OLE格式的字段保存,那么在窗体中可以直接显示出图片来.但是这样做有以下不足:一.需要将图片逐一插入到表中,工作量太大.二.使数据库文件变得庞大.三.相同的图片文件, ...
- c语言随机函数怎么循环,如何在C ++循环中生成不同的随机数?
千万不能使用rand(); 使用新的C ++ 11的设施(如std::mt19937,std::uniform_int_distribution等)来代替. 你可以使用这样的代码(在Ideone上住这 ...
- 如何在bat脚本中列出指定目录下的所有文件信息
今天在研究windows系统下的bat脚本,想要去遍历一个目录并输出该目录下的所有文件的详细信息.经过一番查找和摸索,使用以下代码即可实现: @echo offset source_dir=" ...
最新文章
- Dagger2的使用
- IIS7.5配置对PHP的支持
- How to Fix “Username is not in the sudoers file. This incident will be reported” in Ubuntu
- pb的webserver增加的方法发布后没有显示_震惊!!!Diboot 2.0.5 发布,让开发工作又快又爽...
- 软件测试缺陷定义和管理
- 基于麻雀算法优化的Tsallis相对熵图像多阈值分割 -附代码
- delphi mysql ado_delphi2010利用ADO连接MySQL数据库
- 凸优化读书笔记01(仿射集合、仿射组合,仿射包)
- 聊聊Hadoop DistCp的数据切分处理方式
- webp格式以及工具介绍
- 渐进式Express源码学习2-道士下山
- AndroidStudio的强大搜索功能介绍!全局搜索、搜索文件、搜索类、搜索文本、搜索一切
- 如何优化我的世界服务器,我的世界服务器太卡怎么办 MC服务器优化攻略
- 计算机图形学——游戏方向 第一章 计算机图形学概述
- Access to script at ‘file:///D:/html/vue-%E6%A8%A1%E5%9D%97%E5%8C%96%E5%BC%80%E5%8F%91/js/aaa.js‘ fr
- Windows10系统登陆界面“出现问题,PIN不可用”-解决方法
- 运行 python manage.py runserver
- 阿里云的这群疯子- 文/史中
- Cannot access ‘state‘ before initialization
- 深入剖析Spring Web源码(十九) - 整理的文档和日志的索引(第一版)
热门文章
- Matlab数据线性化
- Delphi学习第一课
- 实现LAYERED窗口
- 虚拟拨号服务器怎么用,windows实现虚拟拨号服务器
- Android自定义View之仿金山词霸加载效果
- 盒子综合案例——德云社十八愁与宠物知识栏
- php使用 163邮箱接口,G. PHP发送邮件功能实现(使用163邮箱)
- 什么是CDN,网站被攻击时该怎么防
- php error unexpected,PHP错误syntax error unexpected T-FUNCTION的解决方案-深圳做网站-创络...
- WIFI WDS不同应用模式简介