《可视化计算》课程作业2015年（冬）
姓名： 年级：
学 号： 学习中心：
————————————————————————————— 注意：所有图像的标题必须呈现足够你本人信息
1.（10分）求解下列线性方程组的解： u1.5v2x9y7z33.6v0.5x4y4z47u10v3x22y33z20 3u7v8.5x21y6z53u8v90y20z16
源程序： clear;
A=[1,1.5,2,9,7;0,3.6,0.5,-4,-4;7,10,-3,22,33;3,7,8.5,21,6;3,8,0,90,-20];
B=[3;-4;20;5;16]; C=A\\B C =
3.5056 -0.8979 -0.2745 0.1438
0.0137 图形：
2.（10分）信号y = 5*sin(pi*20*t)+3*cos(2*pi*50*t)幅度为1的白噪声的干扰，请画出此信号，然后进行傅立叶变换，画出变换后的频域信号。
源程序： clear;
t=0:0.001:0.6;
y=5*sin(pi*20*t)+3*cos(2*pi*50*t); y=y+randn(1,length(t)); subplot(2,2,1); plot(y(1:100));
xlabel('变换前的信号'); x=fft(y,256); subplot(2,2,2); x=real(x);
plot(x(1:256));
xlabel('变换后频域信号');
3.（10分）在空间有一个带正电的点电荷，请画出此点电荷的空间电位分布与电场强度的空间分布图。
源程序： 程序1： clear;
x=-2:0.1:2; y=-2:0.1:2;
[x,y]=meshgrid(x,y);
z=1./sqrt(x.^2+y.^2+0.01); mesh(x,y,z);
xlabel('电荷空间分布');
程序2： clear;clf; q=2e-6; k=9e9;
x=-6:0.6:6; y=x;
[X,Y]=meshgrid(x,y); r = sqrt(X.^2+Y.^2); V=q*k./r;
[Ex,Ey]=gradient(-V); AE=sqrt(Ex.^2+Ey.^2); Ex=Ex./AE; Ey=Ey./AE;
contour(X,Y,V,'r-'); hold on;
quiver(X,Y,Ex,Ey,.5,'b')
图形：
4.（10分）设计一个低通滤波器，从混合信号：
x(t)=sin(2*pi*10*t) + cos(2*pi*100*t) + 0.2*randn(size(t)) 中获取10Hz的信号（10分）。
源程序： clear; ws=1000;
t=0:1/ws:0.4;
x=sin(2*pi*10*t)+cos(2*pi*100*t)+0.2*randn(size(t)); wn=ws/2;
[B,A]=butter(10,30/wn); y=filter(B,A,x); plot(t,x,'b-') hold on
plot(t,y,'r.','MarkerSize',10) legend('Input','Output') xlabel('第四题');
图形：
可视化作业第四题-
5.（20分）设计一个程序，应用函数subplot(1,2,1)、subplot(1,2,2)分别显示您本人的二张照片，然后对二张照片分别进行傅立叶变换，并分别画出变换后的频域信号。再把2个频域信号相加，经傅立叶逆变换后，显示时域信号的图像。
源程序： clear;
ben=imread('ben.jpg'); ben2=imread('ben2.jpg'); figure(1);
subplot(1,2,1) imshow(ben);
title('照片1-'); subplot(1,2,2) imshow(ben2); title('照片2-'); F1=fft2(ben);
FS1=fftshift(F1); F2=fft2(ben2); FS2=fftshift(F2); figure(2);
subplot(1,2,1)
imshow(FS1,[]);title('傅立叶变换之后1-'); subplot(1,2,2)
imshow(FS2,[]);title('傅立叶变换之后2-');
k=min(min(FS1),min(FS2)):max(max(FS1),max(FS2));
s1=zeros(1,length(k));s2=s1s1(find((k>=min(FS1))&(k<=max(k1))==1))=FS1;
s2(find((k>=min(FS2))&(k<=max(k2))==1))=FS2; y=s1+s2; z=ifft(y);
imshow(z);title('逆变换后图像');
图形：
可视化作业第五题-
6.（30分）小论文
根据工作中的实际需要，请设计一个实际工程问题的可视化。可以选择以下之一：（1）工程动画的可视化；（2）大数据处理中的可视化；（3）算法与模型计算的可视化；（4）实际生产流程的可视化；（5）或其它有创新意义的可视化科学计算。要求： （1）题目有实际意义。 （2）有分析、算法描述 （3）程序源代码设计。
（4）问题结果有可视化显示。 （4）题目的问题有一定的新意。
小论文的字数不能少于2500字，格式由下列各部分组成：
中文题目
摘要：
中文关键词：
英文题目：
英文摘要：
英文关键词：
1.引言
2.算法基础 3.程序代码
4.结果分析（必须有可视化图） 5.结论
参考文献
基于MATLAB关于交混回响
的仿真研究
摘要： MATLAB是进行科学研究和数值分析的首选软件，本次对于交混回响的研究正是借助MATLAB的强大功能，通过编写程序实现了模拟语音信号的数字化处理，然后绘图直观的对混响现象进行仿真分析。
中文关键词：MATLAB；数字化；仿真
Simulation study of MATLAB
based on cross reverberation
Abstract：MATLAB is the first choice for scientific
research and numerical analysis. The research on the mixed reverberation is the powerful function of MATLAB. The digital processing of analog voice signal is realized by writing program.
Key words：MATLAB; digital; simulation
1.引言
MATLAB科学与工程计算的高级语言, 同时也是适用于科学和工程计算的数学软件系统。作为数值软件的MATLAB , 对计算数学领域的特定问题类型, 如数值线性代数、微分方程数值解、数值逼近、最优化方法、数值积分等进行求解, 可以得到问题的离散解或近似的解析表达式。针对每一类问题的求解, MATLAB 都能给出该类问题的各种高效算法。MATLAB 具有一个阵容强大、范围广泛的基本运算体系。同时，MATLAB基于数值运算的信号处理技术也使DSP的相关技术变得简单易行。此外,MATLAB 极其简洁的语言结构支持结构化程序设计, 且具有超过任一高级语言系统的调试技术。MATLAB 高度智能化的调试技术是用户开发MATLAB 程序的有力工具。
本次仿真研究利用MATLAB强大的绘图功能和声音输入输出的功能对.WAV的音频文件进行了数字化的采样，然后对数值进行了相应处理，通过绘图直观的显示了声音传播过程中的交混回响现象。
2.算法基础
由于本次仿真是对语音信号的数字化处理。下面介绍
MATLAB处理音频信号的一些理论依据和相应函数功能。
2.1 信号分类
按信号特点的不同，信号可表示成一个或几个独立变量的函数。例如，图像信号就是空间位置（二元变量）的亮度函数。一维变量可以是时间，也可以是其他参量，习惯上将其看成时间。信号有以下几种：
1．连续时间信号：在连续时间范围内定义的信号，但信号的幅值可以是连续数值，也可以是离散数值。当幅值为连续这一特点情况下又常称为模拟信号。实际上连续时间信号与模拟信号常常通用，用以说明同一信号。
2．离散时间信号：时间为离散变量的信号，即独立变量时间被量化了。而幅度仍是连续变化的。
3．数字信号：时间离散而幅度量化的信号。 本次仿真正是通过把模拟音频信号进行采样得到数字信号，再在此基础上借助MATLAB的一些工具完成的仿真。在数字化过程中，若是立体声，则以某采样频率采样产生n×2的数字信号序列；若是单声道，则采样生成n×1的数字信号序列。
2.2 函数介绍
MATLAB中对声音的处理以及矩阵运算和绘图提供了丰富的函数，下面介绍几个本次仿真过程中要用到的函数：
1. [y,Fs]=wavread(‘filename’,[N1 N2])； 该函数可以实现.wav格式音频文件的直接读取。其中y为返回的样值序列(如果是单声道，返回1×n序列；是立体声返回2×n序列)，另外y的取值范围要在[-1,1]内；Fs为采样率，需满
足取样定理(Fs大于信号最高频率的两倍)；N1和N2用来限定取样的起始和结束点。
2. zeros(m,n)；用来产生m行n列的全零矩阵。仿真过程中处理信号的延时传递效果，便是在序列前面加上适当长度的全零矩阵。
3. plot(x,y,’string’)，subplot(m,n,p)；为绘图函数，其中x表示绘图的横坐标，y表示绘图的纵坐标，string可以用来选择绘图的颜色，线条宽度等参数；subplot表示将界面分成m行n列，选中第p个选区为当前绘图区域。
4. *、+、^、’；为算术运算符，依次为：乘号、加号、乘方号、对矩阵转置运算符。
2.3 关于交混回响
建筑物里，由于距离较小，反射声波（回声）跟原来的声音混合在一起分不开，感觉效果是原来的声音加强并且延长了，这种现象称为交混回响。声源发出声音后，声强减弱到原来声强的百万分之一所需要的时间叫做建筑物的交混回响时间，是描述建筑物的声学特性的重要参数。适当的交混回响时间可以使声音变得浑厚、好听，增加感染力，教堂和音乐厅建筑都注意这一点。当交混回响时间太长会使音乐的节奏和语言模糊不清。交混回响时间一般以1秒-2秒为宜．交混回响时间与建筑物大小，墙壁、地面和天花板材料都由关系，还受到建筑物内人的多少以至服饰的影响。例如，首都剧场空座时的交响回响时间是3.3秒,满座时的交混回响时间大约是1.36秒。人民大会堂满座时的交混回响时间大约是1.8秒。
2.4 交混回响仿真实验
为了简化问题，假定墙面上有吸音材料，即忽略墙壁的反射。编程仿真的思想是：电信号在电线里传输的时间可以忽略，但是当相距若干米的扬声器发出的声音到达听者的位置时，不同扬声器的声音经过了不同的路程，路程的差别产生相位差，所以带有不同相位差的声音叠加就产生了交混回响的效果。
图3-1所示是交混回响仿真实验中的扬声器布置图，A、B、C是三个接到同一声源的扬声器，它距离听者的水平距离和垂直距离标注在图上。
通过在MATLAB上编写程序，可以通过图形象的展示出混响效果。程序中a.wav文件为一段下载的音频文件；听者听到各扬声器的声音处理方式为：计算扬声器到听者距离然后除以声速，得到传播时间，再将该时间乘以每秒取样点速率fs，可以得到离散点数目，再由zeros函数得到相应的1×n零矩阵，加到原始声音序列中合成即可。
图 3-1 交混回响仿真实验中的扬声器布置
3.程序代码
程序1-1：
x1=wavread('D:\\a.wav'); % 对a.wav进行取样，取样率为默认值，返回值存x1
y=x1(34100:44100)'; % 取34100到44100的采样点，转置后存y
fs=8000; % 取样速率 N1=10100;
v=330; % 声音速度
b1=(120^2+180^2)^.5; % 计算听者到扬声器B的距离 yb1=[zeros(1,(fs*b1/v)),y']'; % 听者听到扬声器B的声音 yb=yb1(1:N1); % 取1到N1离散值，存yb a1=(120^2+30^2)^.5; %计算听者到扬声器A的距离 ya1=[zeros(1,(fs*a1/v)),y']'; % 听者听到扬声器A的声音 ya=ya1(1:N1); % 取1到N1离散值，存ya c1=(200^2+200^2)^.5; % 计算听者到扬声器C的距离 yc1=[zeros(1,(fs*c1/v)),y']'; % 听者听到扬声器C的声音 yc=yc1(1:N1); % 取1到N1离散值，存yc ye=1.6*ya+2*yb+1.2*yc; % 听者听到三个扬声器声音叠加 subplot(211);plot(y,'r');grid % 作原声音时域波形图 subplot(212);plot(ye);axis([0,12000,-1,1]);grid % 作混叠声音时域波形图
4.结果分析
根据程序1-1在MATLAB中运行后，可以得到如图3-2的时域波形图。由图形可以看到，由于交混回响作用，听者要在声音发出一段时间后才能接受到信号；同时由于几个扬声器的混叠作用，使声音的各部分有放大或是减弱的效果。若是合理利用交混回响，可以在现实中使声音效果更丰富，更具空间感；同时由于交混回响作用，地震波在水层中的多次反射，也包括鸣震，造成有时在地震记录上以不同振幅
在多处出现。在海洋地震勘探中，它是一种极其严重的干扰。在陆地地震勘探中偶尔也可以遇到交混回响。
图 3-2 声音信号的时域图(上：原信号；下，听到的交混回响的信号)
5.结论
通过本次基于MATLAB关于交混回响的仿真研究，使我对MATLAB的基本操作有了很好的熟悉和掌握；在翻阅资料和学习课本知识的过程中，基本熟悉了通过MATALB软件仿真现实问题的方式方法。再次，由仿真研究的各个流程，可以充分体现MATLAB在数值运算，数字信号处理，绘图等方面的强大优势；其编程语言的通俗与简洁，也使其当之无愧成为进行科学研究和数值计算的首选软件。
参考文献
[1] 周品 何正风·MATLAB数值分析[M]·机械工业出版社，
2009.1
[2] 徐明远 邵玉斌·MATLAB仿真在通信与电子工程中的应用[M]·西安电子科技大学出版社，2009.9
[5] ·MATLAB在计算机辅助数值分析与计算中的应用 [EB/0L]·
http://wenku.baidu.com/view/feb66e4f767f5acfa1c7cda7.html
因篇幅问题不能全部显示，请点此查看更多更全内容