Matlab与Simulink仿真学习心得

班级：07610 学号：072016 姓名：吕天雄

一 Matlab学习心得体会与编程实践

<1>学习Matlab的心得体会

真正开始接触Matlab是大二上就开始了，到现在已经一年多了，在此之间，Matlab的确为我提供了很多便利。Matlab的确不愧成为是草稿纸上的语言。我们不必去为很简单的显示效果图形去找一些什么其他软件或者研究比较复杂的计算机图形学，一个plot或者别的函数往往就可以得到很满意的效果。

其实最初开始学习matlab的时候感觉这个东西和C没什么两样，但是后来具体到一些 东西，比如信号处理和数学建模上以后才感觉到使用matlab编写程序去验证结果比C要节 省很多时间，而且matlab写东西基本都是按照自己的思路平铺直叙很少去考虑什么函数的 嵌套调用或者指针等等很头疼的东西。

关于matlab的学习，我感觉其实百度和matlab自带的help基本能够解决绝大数问题， 而且一些比较好的论坛比如www.ilovematlab.cn都会为你产生很大的帮助，关键是在于多动手实践，多思考。但是matlab毕竟只是一个工具，原理和一些基本的编程素质还是必须有的，否则matlab最多也只能是验证一些别人的东西而已，根本帮不上什么忙的。

<2>遇到的一些问题的思考方式与解决办法

最开始用matlab的时候是在大物实验，实验要求去根据测量得到的数据作出图。但是 手动用铅笔去画确实很麻烦，所以用matlab确实可以省去很大的麻烦。但是第一次遇到问 的时候是有关极化坐标下的曲线拟合。

首先是一个物理实验的问题；在做一个关于光的偏振的实验的时候，最后的结果要在一个极化坐标下显示出来；因为数据是离散的，所以显示出来的图像是一个折来折去的一个东东； 然后很自然的想法是对这个曲线进行插值处理。

但是极化坐标下MATLAB并未提供插值处理的函数，interp1这个函数只能在笛卡尔坐标系，也就是直角坐标系下使用。

然后就想到把极坐标的数据转换的直角坐标系下，pol2cart可以实现这个想法，但是随后而来，也就是最后导致整个问题失败的关键也在这里。

pol2cart以后产生的一串数据中出现了重复的数据，那么interp1这个东西也就无能为力了，因为interp1不能处理一串数据中有重复出现的情况。 最后的处理办法是把这些数据c os，sin这些东西变换一下后，使其大致规则，然后再用polar画出极坐标下的图形。

接下来这个问题就有点超过我的范围了，可能会有点叙述不清楚。问题可以概要为：人脸网格插值。

这是一个用三角形网格表示的一个人脸模型。需要通过插值使其变得光滑，当然问题的复杂之处在于，插值会改变原本的网格结构。

对于这个问题许多人都给出了解决的办法，当然是一些关于人脸识别技术和运动图像处理的范畴之内的。

主要有两种，第一种是face—ls算法。这种算法是基于RBF（径向插值函数）和loop细分原则的一种算法，当然其精确程度比较差点。但是速度挺好。原理是：基于网格上的点，产生出顶点的迭代函数，从而产生出新的定点，进行细化。然后进行插值。也就是引进新的顶点。分为两个过程内部插值和边缘插值。具体可以去看兰州大学信息工程学院的学报。 第二种也就是基本的老办法RBF。首先解释一下什么叫做径向插值 径向基函数

主要是考虑多维空间的数据插值问题，径向基函数在三维图形的变形中常常被采用，用其来变分三维人脸的时候还应结合人脸的固有特征。确定该方法来实现三维人脸变形是可行的，但要构造好的基函数，以及解这些高维方程的解仍是难点。 插值方法：

假设大多数人脸的形状都可以由一个拓扑原型变化得来，那么，通过调整一个一般模型的构造参数可以建立不同的面部模型。但是，这种参数模型仅仅局限于那些构造参数已知的情况，并且对特定人脸参数的调整非常困难。在离散数据的多变量插值问题方面，径向基函数（radial basis function，RBF）插值方法是一个行之有效的工具，所以也适用于类似人脸这样高维曲面的近似或平滑插值。现有的许多方法使用了基于RBF的插值技术，将一般人脸网格变化到特定人脸的形状。这种方法的优点在于：（1）通过插值可以得到丢失的数据点，所以源网格和目标网格不需要相同数目的结点；（2）如果选择了合适的匹配点，数学上可以保证能够将源网格变形到目标网格。

当然我也不是很懂，只是理解了部分。这种方法的关键是找到一个很好的核函数来计算出新的顶点。

最后问题解决的方式是查着了一些图书馆的学报。然后而且在网上找到了部分代码。然后做出来的。这是我校大四一个同学毕设中的一个东西。

其实这个事件最好的解决办法是在3D-MAX中，不过因为没能找到MATLAB与3D-MAX是怎么接口的，所以作罢了。

第三次遇到的问题是我在数学建模课上遇到的一个关于矩阵LU分解的问题：

Matlab作出的结果和手算的结果竟然不一样。

10111 U=1 一个矩阵 A= 经过三角分解为L=10111因为1/为一个极大数 则1—1/可以看成—1/；

设=0.000000000000000000000001%10的负24次方

从而L*U=

则最后一位1没有了 变成了0

其中L=

U=

Matlab 程序为

epusino=0.000000000000000000000001; U=[epusino 1;0 -1/epusino]; L=[1 0;1/epusino 1]; L*U ans =

0.000000000000000 1.000000000000000 1.000000000000000 0

而直接用lu函数对矩阵A分解的结果为 A=[epusino 1;1 1] A =

0.0000 1.0000 1.0000 1.0000

>> [L,U]=lu(A) L =

0.0000 1.0000 1.0000 0 U =

1 1 0 1

我们老师当时的解释是这是由于我们所用matlab是破解版的缘故，破解版的计算精度 没正版的高所以导致运算结果的错误，后来我看了一些LU这个函数的help文档 psychologically lower triangular matrix\" (i.e. a product of lower triangular and permutation matrices) in L。翻译过来是L是一个心理上的下三角矩阵，其实是下三角矩阵和置换矩阵也就是最后的结论。 P*A = L*U.。由此可见这并不是什么精度的问题导致的结果只是LU分解的矩阵意义和书本上的不同而已。

以下是我自己写的一些程序的代码

1：网络随机拓扑图

目的是要生成一个度数随机、权值随机的拓扑图 。 首先的理解是利用邻接矩阵，先随机生成一个矩阵，该矩阵为一个对称矩阵。然后画出这个矩阵就行了。

随机拓扑图10026516398191 5908070160550440430 253035404550556065

%产生数组A用来存放表示两点之间权值的矩阵A，也就是临接矩阵，那么两点之间权值不为零元素的个数即为该点的度数

DEF=5;%设定一个东东 方便改变随机点的个数 A=rand(DEF,DEF);%产生DEF*DEF的随机矩阵 for i=1:DEF

A(i,i)=0%将对角线上的数置为0 end A=10*A;

A=floor(A);%向下去整 for i=1:DEF for j=1:i

A(j,i)=A(i,j)%将A矩阵变为一个上三角或者下三角矩阵 end end

x=100*rand(1,DEF);y=100*rand(1,DEF);%产生10个随机的点 plot(x,y,'r+');

for i=1:DEF

a=find(A(i,:)>0)%将A矩阵每行大于0的数的在该行的地址找出来放在a中

for j=1:length(a)

c=num2str(A(i,j)); %将A中的权值转化为字符型

if c~='0'%不显示为0的值 因为A矩阵为零代表两点不相连

text((x(i)+x(j))/2,(y(i)+y(j))/2,c,'Fontsize',18);%将权值显示在两点连线中间 end hold on;

line([x(i) x(a(j))],[y(i) y(a(j))]);%连线 end end

title('随机拓扑图'); e=num2str(DEF);

legend(e);%左上角显示节点的个数 for m=1:DEF A(m,m)=m; f=num2str(A(m,m));

hold on;text((x(m)+x(m))/2,(y(m)+y(m))/2,f,'Fontsize',18);%将权值显示在两点连线中间 end

接下来是一个在信息安全课上写的一个关于256色图分层的程序，一副图像有m*n个像素然后每个像素是一个8bit的二进制数据换算为十进制是0-255之间。要做的就是把每bit的信息提取出来。

I=imread('bupt副本.bmp');%读入源图像

I=double(I);%将图像转换为double类型便于MATLAB中的计算 %确定图像的长宽 M=size(I,1);%长 N=size(I,2);%宽 l=M*N;%图像长宽之积 for n=1:8;%剔除每层 for i=1:M for j=1:N

B=numdec2bin(I(i,j),8);%先将每个像素转换为8位2进制序列 B(n)=[];%删除相应位置的元素

y(i,j)=numbin2dec(B);%将剩下的7个元素转换为10进制数放入y矩阵中 end end

y=uint8(y);%将10进制的矩阵转换为无符号整形

subplot(4,2,9-n);%显示n个图像，每层是相应的剔除该层的结果 imshow(y);%显示这个图像

title(strcat('去除第',num2str(9-n),'层后的结果'));%显示标题 end

附带的两个函数：

function y=numdec2bin(x,n);

%从函数将十进制数x转化为n位二进制 y=zeros(1,n); a=x;

i=n;

while a>0

y(i)=mod(a,2); a=floor(a/2); i=i-1; end

function y=numbin2dec(x) %将二进制数转换为十进制数 a=0;

lx=length(x); for i=1:lx

a=a+x(i)*2^(lx-i); end y=a;

去除第1层后的结果去除第2层后的结果去除第3层后的结果去除第4层后的结果去除第5层后的结果去除第6层后的结果去除第7层后的结果去除第8层后的结果

Matlab的RGB转换为YCbCr之后转换回来和原图不符

这个是做JPEG图像压缩的时候遇到的一个问题，当时的代码是用C写的但是结果不对，后来想到用matlab进行验证结果发现通过公式转换得到的图是错的。

原理：BMP图像压缩为JPEG的第一步是将RGB色彩空间通过这个公式映射到YCbCr空间上。

Y=0.299 R + 0.587 G + 0.114 B

Cb = - 0.1687R - 0.3313G + 0.5 B + 128 Cr = 0.5 R - 0.4187G - 0.0813 B + 128

然后再进行量化，DCT，编码等等步骤。JPEG解压时需要将YCbCr空间的图像又转化回来到RGB上。

R = Y + 1.402 (Cr - 128)

G = Y - 0.34414 (Cb- 128) - 0.71414 (Cr - 128)

B = Y + 1.772 (Cb - 128)

但是通过书上给的公式和网上大部分公式却发现根本转化不会来。下面是写的程序代码以及测试得到的结果图像。

clc,clear;

Source=imread('hl.jpg');%读入原始RGB图像

figure(1); subplot(1,2,1);

imshow(Source):title('original image');%显示图像 [r c d]=size(Source);%计算图像大小 %------计算红色分量并显示分解图------% R(:,:,1)=Source(:,:,1); R(:,:,2)=zeros(r,c); R(:,:,3)=zeros(r,c); R=uint8(R); whos; figure(2); subplot(1,3,1); imshow(R)

title('Red Component');

%-------计算绿色分量并显示分解图-------% G(:,:,2)=Source(:,:,2); G(:,:,1)=zeros(r,c); G(:,:,3)=zeros(r,c);

G=uint8(G); figure(2); subplot(1,3,2); imshow(G)

title('Green Component');

%--------计算蓝色分量并显示分解图-------% B(:,:,3)=Source(:,:,3); B(:,:,1)=zeros(r,c); B(:,:,2)=zeros(r,c); B=uint8(B); figure(2); subplot(1,3,3) imshow(B)

title('Blue Component'); %------------合成-------------% Comp(:,:,1)=R(:,:,1); Comp(:,:,2)=G(:,:,2); Comp(:,:,3)=B(:,:,3); figure(1); subplot(1,2,2);

imshow(Comp):title('composition image'); Y=0.229*R+0.587*G+0.114*B; Cb=0.5*B-0.1687*R-0.3313*G+128; Cr=0.5*R-0.4187*G-0.0813*B+128; red=Y+1.402*(Cr-128);

green=Y-0.34414*(Cb-128)-0.71414*(Cr-128); blue=Y+1.772*(Cb-128); Comp2(:,:,1)=red(:,:,1); Comp2(:,:,2)=green(:,:,2); Comp2(:,:,3)=blue(:,:,3);

figure(3);imshow(Comp2);title('RGB转换为YCrCb后又转换为RGB的图像'); RD=R(:,:,1)-red(:,:,1); GD=G(:,:,2)-green(:,:,2); BD=B(:,:,3)-blue(:,:,3);

figure(4);subplot(1,3,1);imshow(RD);title('红色分量差异'); subplot(1,3,2);imshow(GD);title('绿色分量差异'); subplot(1,3,3);imshow(BD);title('蓝色分量差异');

然而可以通过figure(4)的图像清楚看到R 和G分量转换后恢复 与以前的差别很大。而且看了一下matlab自带的rgb2ycbcr和ycbcr2rgb并且用了这两个函数测试后图像是一样的。很是不解

实验得到的图为:

original imagecomposition image

将源图像分为RGB三个部分。

Red ComponentGreen ComponentBlue Component

将源图像的RGB分量转换到YCbCr上后又转化为RGB得到的三个分量差异的图。

红色分量差异绿色分量差异蓝色分量差异

公式转换的y分量公式转换的cb分量公式转换的Cr分量matlab函数转换的y分量matlab函数转换的Cb分量matlab函数转换的Cr分量

从最后的结果可以看出原图转换到YCbCr空间之后又转化回来得到的图像与原图不符….

RGB转换为YCrCb后又转换为RGB的图像

以上只是我自己用Matlab做过的部分问题，还有一些比如场声源定位Music算法仿真，DES加密算法。神经网络滤波器等等，鉴于代码量太大，就不再敖述。

二 Simulink学习心得与编程实践

SIMILINK模块库按功能进行分类，包括以下8类子库： Continuous（连续模块） Discrete（离散模块）

Function&Tables（函数和平台模块） Math（数学模块）

Nonlinear（非线性模块）

Signals&Systems（信号和系统模块） Sinks（接收器模块） Sources（输入源模块）

对于simulink来说其实没什么好说的，需要什么就把什么拖到Model里面，连线之后设定参数。不是很麻烦的。但是重要的是对于系统的设计与仿真思想才是最重要的。而且里面的S-Function也是比较好的一个设计，用S-Function可以自定义一些组件，使自己的仿真看上去更加清晰与精简。

Simulink的工作方式：

（1） 模块内的参数值首先会送到Matlab中进行计算，得到的参数值会用来当做以

后需要调用的参数。

（2） 模型系统中的各个层级将被平展开来，每一个子系统将被相应的模块所代替。 （3） 模块按被处理的顺序排列，此时代数回路结构也被检查出来，此种排列产生一

个列表，以确保具有代数回路的模块驱动输入的模块被更新后才更新。

（4） 检查块之间的链接，是否每一个块的输出端口与它所连接的模块输入端口有相

同的信号宽度。

现在可以准备执行仿真操作，仿真时使用数值迭代求的的结果，每种数值积分模型提供的连续状态的微分能力。

Simulink中的模型都是分级的，因此可以通过自上而下或者自下而上的方式建立模型。定义了一个模型以后，就可以通过Simulink的菜单或者在Matlab的Command中输入命令进行仿真。

关于学习的心得就写到这里了，接下来是一些自己做的仿真。

1：观察一个信号与积分之后的区别，目的是了解Scope的功能和用法。

仿真结果，左边为原始信号，右边为积分后的信号。

用XYGraph看到的对比

Rossler吸引子产生仿真：

Rossler 吸引子产生是服从下面这个方程的

23 1

212

331

仿真图

(xx)xxaxxbx(xc)xab0.2c5.7

仿真结果

LMS自适应滤波器： 仿真图：

遇到的问题

问题的解决办法：

将Spectrum Scope中的buffer input打勾就行了。

仿真结果图形：

1：Spectrum Scope的输出图形

2：Vector Scope的输出结果

3：Scope的输出结果图形