实验三 动画技术:电偶极子辐射的动态仿真

实验三 动画技术：电偶极子辐射的动态仿真

一、实验目的

物理过程或物理现象通常都是动态过程，因此对于物理过程的仿真或模拟应该也是动态。通过对物理过程的动态仿真能够近似地还原物理过程，帮助我们更好的理解物理现象和物理过程，揭示蕴藏其中的规律性东西。本次实验将以电偶极子天线的电磁波辐射动态仿真为例，介绍MA TLAB 的动画技术，以期实现如下目的：

1. 掌握两种MA TLAB 的动画制作的技术：影片动画和实时动画；

2. 掌握矢量场力线图的制作方法，并了解电偶极子辐射的规律，以便更好的理解《电磁场与电磁波》课程中的相关知识点，也为进一步学习其他专业课程(如天线原理、天线技术) 建立基础。

二、实验预备知识

1. MATLAB 动画技术

MA TLAB 提供了两种制作动画的方法：影片动画和实时动画。

(1) 影片动画

这种动画技术类似于电影的制作，其原理是首先对仿真的过程按时间次序进行“拍照”，获得一帧一帧的画面(称为帧) ，并将之存档，然后再按时间顺序以高于视觉暂留的帧频率播放帧，即可获得类似于电影的动画效果。这种动画技术适用于难以实时快速绘制的复杂画面，计算量大，占用内存较多。MA TLAB 提供了下列几种函数用于实现影片动画：

① moviein 函数 该函数将产生一个结构体数组（structure ，以下称帧结构体）来存放动画的帧（即所拍摄的一幅幅画面），每帧画面作为结构体的一个元素保存。调用格式 fmat = moviein (N)

产生一个能存放N 个帧的（1×N ）结构体数组fmat 。该结构体包含两个域cdata 和colormap ，前者存放帧的图像数据，后者存放帧使用的颜色表。

② getframe 函数 该函数作用是对当前的图像进行快照(“抓拍”) ，通常有两种使用格式: getframe “抓拍”当前坐标轴（一种图形对象）里的内容；

getframe(h) “抓拍”某个图形窗口或坐标轴里的内容，该图形窗口或坐标轴以句柄h 标识（图形窗口和坐标轴都是一种图形对象，每一种图形对象都有自己特有的句柄handle ，即标识，类似于“身份证”）。例如

>> fmat(1) = getframe (gcf)

抓拍当前图形窗口下的内容，并将该帧画面存放于帧结构体fmat 的第一个元素中；gcf 为Get handle to current figure的缩写，意思是获取当前图形窗口的句柄。在命令窗口中输入gcf 可显示当前图形窗口的句柄值，是个整数。再如

>> fmat(1) = getframe (gca)

抓拍当前坐标轴内的内容，并将该帧画面存放于帧结构体fmat 的第一个元素中；gca 的意思是获取当前坐标轴(axis)的句柄(是双精度数) 。

注意，两种格式所“拍摄”的画面一般是不一样的。

③ movie 函数 作用是按顺序回放帧结构体fmat 中的各帧画面。通用格式: movie (h, famt, n, fps, loc)

h 是播放动画的图形窗口或坐标轴的句柄，缺省时表示在当前的坐标轴中播放动画； famt 是帧结构体，不可省去； n 是重复播放的次数，缺省时，只播放一次； fps 代表每秒播放的帧数(帧频) ，应快于视觉暂留，缺省时系统默认fps=12。

loc 指定动画播放在图形窗口或坐标轴内的具体位置，为4元素的行向量[x,y,0,0]，头两个元素指定动画播放的相对位置(相对于由h 标识的图形窗口或坐标轴的左下角来计算) 。

在动画播放前，movie 函数首先将每帧图像的数据载入内存(此时图像会一帧一帧的显示出来) ，然后再按照用户设定的参数（重复次数n 、帧频fps 等）播放动画。

除了movie 函数，MA TLAB 还提供了一个函数movie2avi ，该函数能够在当前的工作目录下创建一个avi 视频格式的动画文件。一般的调用格式：

movie2avi (fmat, filename)

famt 为前述的帧结构体，filename 是字符串，指定avi 格式的文件名。

例1：模拟地球自转

[XE,YE,ZE] = sphere(200); % 产生单位球数据

surf(XE,YE,ZE); % 绘球体

shading interp; % 插值方式着色，并去掉格线

axis image vis3d; % 带参数vis3d 能避免旋转球体时，图形忽大忽小

axis off;

hold on

h=plot3([0 0],[0 0],[1.3 -1.3],'Linewidth',2,'Color','k'); % 地轴 rotate(h,[0 0],10); % 围绕x 轴将地轴转动10度，h 是图形对象plot3的句柄 delta_fi=12; % 每次旋转球体的角度，单位: 度

N=360/delta_fi; % 转动次数，也即拍摄次数

m=moviein(N); % 创建N 个元素的帧结构体

for i=1:N

view(-37.5+delta_fi*(i-1),30); % 通过改变视角，获得转动效果；

% 其中仰角保持默认的30度不变，方位角每次增加delta_fi m(i)=getframe(gcf); % 在当前图形窗口下抓拍帧

end

movie(gcf,m,10,15,[0 0 0 0]);

image(m(2).cdata); % 展示某一帧（第2帧）图像

程序中为了获得地球自转的效果，使用了view 函数，其作用是指定3D 图的观察角度，调用格式：view (az, el)，参数az 代表观察的方位角(Azimuth)，默认值az =－37.5o ; el 代表观察的仰角(Elevation)，默认值el=30o ；MA TLAB

中对方位角和仰角的定义及正方向规定如图示。

在上例中，若保持仰角不变，不断改变观察者的

方位角az ，则可获得球体转动的效果。

(2) 实时动画

实时动画的原理是保持图形窗口中绝大部分

的像素颜色不变，而只更新部分像素的颜色构成

运动图像。这种动画适用于每次变化少、图形精确度要求不高的场合。

如何显示新对象、擦除旧对象而又不破坏背景图案，在一般的编程语言中不太容易实现，而MA TLAB 图形系统提供的擦除属性EraseMode 可以使上述过程较容易实现。MA TLAB 的绘图命令都可以带擦出属性EraseMode ，其格式为：

Command (…, ’EraseMode’,’属性值’)

command 代表相应的绘图命令，比如plot 、mesh 等，EraseMode 属性的可选值如下表：

当新对象的属性设置后，应使用命令drawnow 立即刷新屏幕，使新对象立即显示出来。drawnow 命令迫使MA TLAB 暂停目前的任务序列而去刷新屏幕；若没有drawnow 命令，MA TLAB 要等任务序列执行完后才会去刷新屏幕。

例2：在实验二的基础上，进一步动态模拟卫星的轨道运动。

G=6.672e-11;

Me=5.97e24; % G 有引力常数，Me 地球质量

Vx0=0;Vy0=4000;x0=-4.2e7;y0=0; %设置初始条件

t0=0; tf=60*60*24*9;

tspan=t0:1600:tf; % tspan为等间隔的时间向量，则ode45输出对应时间的解 y0=[x0 y0 Vx0 Vy0]'; % 初始向量

options=odeset('Events','on');

[t Y Te Ye Ie]=ode45('odefile',tspan,y0,options,G,Me,y0); % 下面在 3D 坐标上画地球

[XE,YE,ZE] = sphere(30); % 产生单位球数据

RE=0.64e7*2; % 将地球半径放大一倍

XE=RE*XE;YE=RE*YE;ZE=RE*ZE;

surf(XE,YE,ZE);

shading interp;

colormap('hot');

hold on;

% 下面画轨道

plot(Y(:,1),Y(:,2),'b','Linewidth',2);

axis image off;

title('地球卫星轨道运行动态模拟’,’color’,’r,’ ,'FontSize',15,... 'FontName',' 隶书');

% 标记远地点和近地点

plot(Ye(2,1),Ye(2,2),'rv','Linewidth',2,'MarkerSize',8)

text(1.05*Ye(2,1),Ye(2,2),'远地点','FontWeight','Bold');

plot(Ye(3,1),Ye(3,2),'r^','Linewidth',2,'MarkerSize',8)

text(Ye(3,1)-3.0e7,Ye(3,2)-2e7,'近地点','FontWeight','Bold'); % 下面画卫星

[XWE,YWE,ZWE] = sphere(10); % 产生单位球数据

RE=0.64e7;

XWE=RE*XWE;YWE=RE*YWE;ZWE=RE*ZWE;

h=surf(XWE+Y(1,1),YWE+Y(1,2),ZWE,'EraseMode','normal' );

% 擦除模式为normal, 绘图命令surf 返回句柄h

shading interp;

% 下面实现卫星动态运动

m=size(Y,1); % 计算解矩阵Y 的行数

for k=1:4 % 重复动画k 次

for i=2:m

set(h,'xdata',XWE+Y(i,1),'ydata',YWE+Y(i,2));

% 通过其句柄h 更新绘图函数surf 绘制卫星时的x 和y 坐标。卫星的运动由其坐标x 和 % y的变化来反映。xdata 、ydata 和zdata 属性控制图形对象的x 、y 和z 坐标。 drawnow; % 立即刷新屏幕

end

end

上述程序绝大部分取自实验二附录中的satellite.m 文件（出于篇幅的原因，ode45指令调用的odefile 函数没有给出，可参见实验二附录中的odefile.m ），其中红色粗体部分是为实现动态模拟所专门设计的。下面说明程序中的第5条命令：tspan= t0:1600:tf。在ode 求解器中，通常输入参数tspan=[t0,tf]表示积分区间，此时ode 求解器在积分迭代时，采用的是变步长，也即输出的不是等间隔的时间点上的解。若tspan=[t0,t1,...,tf]为行向量（其中各个元素须按升序或降序排列），则输出由tspan 指定的时间点处的解。因此输入第5条命令后，ode45指令将输出等间隔(间隔1600秒) 的时间点上的解。为什么要这样做呢？因为我们不仅需要让卫星动起来，而且要能如实反映其运动速度的不同。而速度的差别可由等时间间隔内的位置的差别来反映。上述程序运行后可观察到卫星在远地点的速度明显慢于近地点。

2. 电偶极子辐射的理论分析

当发射电磁波的天线的长度远小于辐射波长时，天线可视为电偶极子。电偶极子的电磁波辐射分析是线形天线工程计算的基础。关于电偶极子天线辐射的详细理论分析过程可参见谢处方编著的第4版《电磁场与电磁波》8.2节所述。下面仅给出有关的结果，电偶极子辐射场的复数形式的电场分布：

E r =

E θ=⎡1j ⎤-jkr cos(θ) ⎢-e 23⎥4πωε0(kr ) (kr ) ⎣⎦2Ilk ⎡j 1j ⎤-jkr sin(θ) ⎢+-e (1)23⎥4πωε0(kr ) ⎦⎣kr (kr ) Ilk 33

E φ=0

其中I 为电偶极子天线的电流。令电偶极子的电量Q =Q 0e

I =dQ

dt =j ωQ 0e j ωt j ωt ，则

代入(1)式，并取实部得到真实的电场分布：

E r =2A cos θ[cos(ωt -kr ) -kr sin(ωt -kr )]/(kr )

2233E θ=A sin θ[(1-k r ) cos(ωt -kr ) -kr sin(ωt -kr )]/(kr ) (2)

E φ=0

其中A =P 0k 3/4πε0, P 0=Q 0l ，P 0为偶极子的电偶极矩的大小。按《电磁场与电磁波》的

1.4.1节所述，可算出球坐标中电偶极子的电力线(矢量线) 满足的微分方程为：

E r

E θd θ=1r dr (3)

代入(2)式并积分，得到电力线方程：

sin θ⋅[cos(ωt -kr ) -kr sin(ωt -kr )]/(kr ) =K (4) 2

式中K 为积分常数，K 取不同的值则得到不同的电力线。因此由(4)式可绘制出电偶极子的电力线族。在绘图时，需要将球坐标还原成直角坐标：

r = θ=cos (z /-1 φ=tan (y /x ) -1

由于电场分布与φ角无关，故电场分布关于z 轴对称，因此可以只考虑某个过z 轴的平面(如xoz 平面) 上电力线图，对于xoz 平面，y=0，因此(4)式中球坐标

r = θ=cos (z /-1 且x 、z 的取值范围均为[-r , r ]。

显然，(4)式可以写成 u (x , z ) =K 的形式，这其实是标量函数u(x,z)的等值线方程，因此电偶极子的电力线方程就是函数u(x,z)的等值线方程。MA TLAB 提供了一个专门的函数用于绘制标量函数u 的等值线(或称等高线) 图：

[c, h] = contour (X, Y, Z, V)

其中，X ，Y ，Z 为同维的矩阵，X ，Y 指定平面上点的x 、y 坐标，可由meshgrid 命令取得，在本例中：

x=–r:0.1:r; y=-r:0.1:r; [X,Y]=meshgrid(x, y);

Z 是函数u(x,y)在坐标X,Y 上的值，V 是向量，指定各条等高线的高度值（例如(4)式中的K 值）。h 是返回的句柄值。

三、实验内容

1. 以影片动画的方式仿真电偶极子辐射过程

要模拟电偶极子辐射场的动态过程，首先要绘制各个时刻的电力线图，即使用contour 函数在ωt 取不同值的情况下绘制电力线方程(4)式。绘制电力线图时应注意下面几个环节： ① 适当选取每个画面上电力线的根数，太多连成一片，太少没有真实感。有2个参数控制

电力线的根数，一为K 值，K 每取一个值代表一条电力线（环形线，见附图），K 的值越多则电力线越多，一组K 值对应一套电力线(族) ；另一个是波数k ，k 越大，电力系将越密，每幅画面将包含更多的电力线数。

② 每个周期内，画面的个数，即适当选取ωt 以及∆ωt 的值，应以感觉画面连续为准。 ③ 最大辐射半径r max 的选取，即x 、z 的范围。r max 越大，x 、z 的范围越大，所画电力线也

越多。其值的选取应以感觉向无限远处传播出去为宜。

根据经验，上述参数可参照下列值：

K=[-2.0,-1.5,-0.8,-0.4,-0.2,0.2,0.4,0.8,1.5,2.0];

k=1；

rmax=3.5*pi;

ωt=n*pi/N，N=24，n=0,1,2,...,N-1，即∆ωt=π/24。N 实际就是“拍照”次数，也是帧结构体的长度，N 越小，动画速度越快。

2. 以实时动画的方式仿真电偶极子辐射过程

由于实时动画以实时擦除的方式实现动画效果，故画面不宜太复杂，否则效果不好。因此电力线数目不宜过多，上述参数需要重新调整。

四、思考题

请根据仿真结果，尝试找出电偶极子的辐射有何特点？

附图：各个时刻的电力线图(红色表示电偶极子天线)

本次实验的程序radiation_movie.m 和radiation_erase.m见附件。