方向图计算

方向图计算（精选4篇）

方向图计算 篇1

0 引 言

天线阵方向图的综合问题是天线阵设计中的核心问题, 是指按规定的方向图要求, 用一种或多种方法来进行天线阵列的设计, 使其产生的方向图与所要求的方向图良好逼近。Dolph[1]利用切比雪夫多项式的性质导出了等距线阵的权值, 可以在给定最大副瓣电平的情况下使波束宽度最小, 或在给定波束宽度的情况下使副瓣电平最低;Villeneuve[2]把Taylor的方法应用于离散直线阵列的综合, 它可以在靠近主瓣的区域产生等幅度副瓣, 而在远处的副瓣电平逐渐减小。其他的用于解决直线阵列方向图综合问题的方法在文献[3]中有详细的描述。

二维矩形阵列天线方向图综合最简单的方法是推广直线阵列的方法。如果口径分布是二维可分离的, 则平面阵的方向图就等于两个正交直线阵的方向图的乘积。由此得到的方向图在相应主平面内达到所要求的副瓣电平, 而在其他平面内副瓣电平低于设计目标, 导致主瓣变宽[4]。

Tseng和Cheng提出了一种设计切比雪夫矩形平面阵的方法[5], 能使平面阵所产生的方向图在每一φ剖面内都是切比雪夫型最佳方向图, 但该方法要求行单元数等于列单元数, 并且不能把阵元的方向图特性考虑在内。

Stutzman和Coffey[6]提出了一种用于平面阵列综合的迭代抽样算法。首先根据其他算法或者实验测量数据得到一初始方向图, 利用正交函数sin (x) /x的性质对其进行采样修正, 使方向图在采样点与目标方向图匹配, 重复该采样迭代过程直到得到满意的方向图。该方法比较简单, 但综合的方向图在采样点之间会有波动, 也不能把单元的方向图特性考虑在内。

文献[7]运用一种改进的粒子群优化算法, 实现了不等幅激励的二维矩形平面阵列天线的方向图综合。通过采用对全局最优粒子微扰和跳变的惯性权重策略, 并使用粒子群算法本身对参数组合进行了优化选择, 很好地改善了优化速度和收敛精度, 使算法具有普遍的适用性, 取得了较好的结果。

该文把自适应算法应用于天线阵方向图的综合, 给出了一种二维矩形平面阵列的方向图综合方法。通过对副瓣区峰值电平的控制, 综合得到了满足要求的目标方向图, 可以把阵列单元的特性考虑在内, 具有较强的实际意义。

1 算法描述

1.1 自适应算法

自适应阵列综合算法是近年来在自适应信号处理理论的基础上发展起来的, 文献[8]比较系统地阐述了这种方法的原理:在干扰入射方向, 自适应阵列会产生零陷, 干扰越强零陷越深;干扰信号个数超过N-2时 (N是天线单元个数) , 自适应阵列将无法在各个干扰方向上形成零陷, 而是自适应地调节方向图使干扰信号功率在输出中最小。

在阵列天线方向图的副瓣区域, 峰值电平是最高的电平, 当峰值电平低于目标副瓣电平时, 副瓣区的电平自然就达到了设计要求。通过只对副瓣峰值电平进行控制, 把自适应算法应用于二维矩形平面阵列天线方向图的综合。

假定有大量的初始干扰信号施加到阵列的副瓣区, 由于干扰的分布特性是已知的, 根据最大信噪比准则求解初始自适应权值, 得到自适应算法产生的初始方向图, 干扰信号方向的电平会降低。在副瓣峰值方向将得到的结果与设计目标相比较, 根据自适应算法的原理调整干扰信号的功率, 如果副瓣太高, 就增加相应角度上的干扰功率, 反之, 就减少相应角度上的干扰功率, 然后再重新计算新的权值。重复进行这个迭代过程, 直到得到满足要求的目标方向图。流程图如图1所示。

1.2 二维矩形阵列的方向图综合

二维矩形天线阵列如图2所示, 图中黑点表示天线单元。此阵沿x轴方向有M个单元, 沿y方向有N个单元, dx和dy分别为行距和列距, (θ, φ) 表示主瓣的指向。天线阵列输出的信号方向矢量和加权值分别为:

$\begin{array}{l}\mathit{V}{}_{2D}=\left[\begin{array}{cccc}v{}_{11}& v{}_{12}& \cdots & v{}_{1{\rm N}}\\ v{}_{21}& v{}_{22}& \cdots & v{}_{2{\rm N}}\\ ⋮& ⋮& ⋮& ⋮\\ v{}_{{\rm M}1}& v{}_{{\rm M}2}& \cdots & v{}_{{\rm M}{\rm N}}\end{array}\right](1)\\ \mathit{W}{}_{2D}=\left[\begin{array}{cccc}w{}_{11}& w{}_{12}& \cdots & w{}_{1{\rm N}}\\ w{}_{21}& w{}_{22}& \cdots & w{}_{2{\rm N}}\\ ⋮& ⋮& ⋮& ⋮\\ w{}_{{\rm M}1}& w{}_{{\rm M}2}& \cdots & w{}_{{\rm M}{\rm N}}\end{array}\right](2)\end{array}$

其中:vmn=fmnejkφmn, fmn是阵列单元 (m, n) 的方向图, k=2π/λ是波数, φmn=mdxsin (θ) cos (φ) +ndysin (θ) sin (φ) , 是各单元的空间相位差。阵列输出信号由每个vmn与相应的复数加权值wmn (即天线单元馈点电流) 相乘并求和得:

$S{}_{2D}={\displaystyle \sum _{m=1}^{{\rm M}}{\displaystyle \sum _{n=1}^{{\rm N}}w}}{}_{mn}v{}_{mn}={\displaystyle \sum _{m=1}^{{\rm M}}{\displaystyle \sum _{n=1}^{{\rm N}}w}}{}_{mn}f{}_{mn}\text{e}{}^{\text{j}k\phi {}_{mn}}(3)$

将上面两个二维阵列V2D和W2D表示成一维形式, 转化成直线阵来综合, V2D和W2D中元素要一一对应, M和N分别是阵元的行数和列数。

$\begin{array}{l}\mathit{V}=[v{}_{11}\cdots v{}_{{\rm M}1},v{}_{12}\cdots v{}_{{\rm M}2}‚\cdots ‚v{}_{1{\rm N}}\cdots v{}_{{\rm M}{\rm N}}]{}^{\text{Τ}}(4)\\ \mathit{W}=[w{}_{11}\cdots w{}_{{\rm M}1},w{}_{12}\cdots w{}_{{\rm M}2}‚\cdots ‚w{}_{1{\rm N}}\cdots w{}_{{\rm M}{\rm N}}]{}^{\text{Τ}}(5)\end{array}$

阵列输出信号可表示为:

$S={\displaystyle \sum _{m=1}^{{\rm M}}{\displaystyle \sum _{n=1}^{{\rm N}}w}}{}_{mn}v{}_{mn}=\mathit{W}{}^{\text{Τ}}\mathit{V}(6)$

天线方向图与V和W相关, 取其模值为:

$p(\theta ,\phi )=\left|\mathit{W}{}^{\text{Τ}}\mathit{V}\right|(7)$

假设阵列的副瓣区有A×B个二维干扰信号Xiab (a=1, 2, …, A, b=1, 2, …, B) , A和B远大于阵列的行数M和列数N, 其入射角为 (θia, φib) , 幅度为Aiab;Xn是各阵元中存在的功率为δ2的高斯热噪声, 取干扰噪声比为 (k表示第k次迭代) :

$\epsilon {}_{iab}(k)=A{}_{iab}(k){}^2/\delta {}^2,\epsilon \ge 0(8)$

天线阵列接收的非期望信号为:

$X{}_u=X{}_n+{\displaystyle \sum _{a=1}^A{\displaystyle \sum _{b=1}^{B}X}}{}_{iab}(9)$

假定第k次迭代时主瓣的最大值为P0 (k) , 如果要求的旁瓣电平低于主瓣电平R (dB) , 则目标副瓣电平为:

$d(k)={\rm P}{}_0(k)/10{}^{R/20}(10)$

主瓣区不施加干扰信号, 即干扰信号功率为0。由上文可知, 第k+1次迭代时干扰信号强度为:

$\begin{array}{l}\epsilon {}_{iab}(k+1)=\{\begin{array}{l}0,(\theta {}_{ia},\phi {}_{ib})\in \text{主}\text{瓣}\text{区}\\ h(k),(\theta {}_{ia},\phi {}_{ib})\in \text{副}\text{瓣}\text{峰}\text{值}\text{方}\text{向}\\ \epsilon {}_{iab}(k),\text{e}\text{l}\text{s}\text{e}\end{array}(11)\\ h(k)=\max \{0,\epsilon {}_{iab}(k)+f[p(\theta {}_{ia},\phi {}_{ib},k)-d(k)]\}\end{array}$

其中:p (θia, φib, k) 为第k次迭代后方向图的实际电平;d (k) 为由式 (10) 求得的目标副瓣电平;f为迭代步长。目标副瓣电平越低, f取值越大;f取值较小可以保证收敛, 但会使迭代过程收敛速度比较慢, 增大f的值可以加快收敛速度, 但f过大会影响收敛过程的稳定性, 在实际应用中可参考文献[8]给出的副瓣电平分别为-30 dB, -35 dB, -40 dB时f的最大取值, 通过试探法确定合适的f值。

主瓣最大值P0 (k) 在迭代过程中是不断变化的, 每次迭代前需重新确定。副瓣峰值所在方向在迭代过程中也是不断变化的, 即需要调整的干扰信号的方向是不断变化的, 所以在每次迭代前也要重新确定, 可以通过比较相邻区域电平值的相对大小得到。

根据最大输出信噪比准则, 由文献[9]可知最佳权值为:

$\mathit{W}=\mu \mathit{\phi }{}_u^{-1}V{}_0^{*}(12)$

其中:V0为期望主瓣方向, 用于确定方向图主瓣的指向;φu为非期望信号的协方差矩阵;μ为任意常数。

$\mathit{\phi }{}_u=E[\mathit{X}{}_u^{*}\mathit{X}{}_u^{\text{Τ}}]=\delta {}^2\left[\mathit{{\rm I}}+{\displaystyle \sum _{a=1}^A{\displaystyle \sum _{b=1}^B\epsilon }}{}_{iab}(k)\mathit{V}{}_{im}^{*}\mathit{V}{}_{im}^{\text{Τ}}\right](13)$

其中:I为单位矩阵;Vim为干扰信号所在方向的方向矢量, 将 (θia, φib) 代入式 (4) 即可求得;*表示共轭, T表示转置。

可设定εiab的初始值为0, 把根据式 (12) 和式 (13) 求得的权值W0代入式 (7) 得到初始方向图, 将其与设计目标相比较, 根据式 (11) 调整干扰信号强度, 由式 (12) 和式 (13) 求得新的权值W1, 代入式 (7) 得到新的方向图, 重复上面的迭代过程, 直到得到满意的方向图, 具体流程如图1所示。

一般阵列方向图为阵元方向图和阵因子的乘积, 上述方法是在没考虑互耦的情况下, 对天线方向图进行的综合。在阵元间互耦较强需要考虑其影响时, 由文献[10]可知, 可以先按照本文算法综合出不考虑互耦时的激励, 再对耦合矩阵求逆, 得到所需要的实际激励。

2 仿真示例与分析

示例1:计算一个6×8矩形平面阵, 阵元行距和列距均为λ/2 (λ为波长) , 主瓣指向侧射方向, 目标副瓣电平取为-30 dB, f=2, 仿真结果如图3所示, 取坐标值u1=sin (θ) cos (φ) , u2=sin (θ) sin (φ) 。加权值以最大值为基准进行了归一化, 如表1所示。

示例2:计算一个7×7矩形平面阵, 阵元行距和列距均为λ/2 (λ为波长) , 阵列单元是电偶极子, 其方向图[11]为cos (θ) , 主瓣指向侧射方向, 目标副瓣电平取为-30 dB, f=2, 仿真结果如图4所示, 取坐标值u1=sin (θ) cos (φ) , u2=sin (θ) sin (φ) 。加权值以最大值为基准进行了归一化, 如表2所示。

3 结 语

该文把自适应算法应用于二维矩形天线阵方向图的综合, 通过只对副瓣峰值电平进行控制, 根据最大信噪比准则调整加权值, 得到了满足要求的目标方向图。由于在每次迭代时只调整副瓣峰值方向的干扰信号, 计算量较小, 使算法得到了简化, 加快了收敛速度, 在具体应用中更灵活, 更易于编程实现。此算法考虑了阵元方向图特性的影响, 在已知阵列的耦合系数矩阵时可以将互耦的影响考虑在内, 克服了传统算法的不足。可以拓展应用于阵元分布不规则的平面阵列的方向图综合问题, 同时副瓣区电平的包络可根据实际需要设定。仿真示例证实了该算法是可行的。

摘要：在干扰入射方向, 自适应天线阵方向图会产生零陷, 从而实现抗干扰的功能。基于最大信噪比准则, 将自适应算法应用于二维矩形天线阵方向图的综合问题, 假定有大量干扰信号施加到天线阵的副瓣区, 通过只对副瓣峰值电平进行控制调整加权值, 得到满足设计要求的目标方向图。该算法计算量较小, 收敛速度比较快, 可以将阵元的方向图特性的影响考虑在内。仿真结果证实了该算法是可行的。

关键词：二维矩形阵列,自适应算法,方向图综合,最大信噪比准则

参考文献

[1]Balanis C A.Antenna Theory Analysis and Design[M].3rdEdition.Hoboken:John Wiley&Sons Inc., 2005.

[2]Milligan T A.Modern Antenna Desigen[M].2nd Edition.Hoboken:John Wiley&Sons Inc., 2005.

[3]Mailloux R J.Phased Array Antenna Handbook[M].Nor-wood:Artech House Inc., 2005.

[4]汪茂光, 吕善伟, 刘瑞祥.阵列天线分析与综合[M].成都:电子科技大学出版社, 1989.

[5]林昌禄.天线工程手册[M].北京:电子工业出版社, 2002.

[6]Stutzman W L, Coffey E L.Radiation Pattern Synthesis ofPlanar Antennas Using the Iterative Sampling Method[J].IEEE Trans.on Antennas Propagat., 1 9 7 5, 2 3 (6) :7 6 4-7 6 9.

[7]袁智皓, 耿军平, 金荣洪, 等.基于改进的粒子群算法的二维阵列天线方向图综合技术[J].电子与信息学报, 2007, 29 (5) :1 236-1 239.

[8]Olen C A, Compton R T.A Numerical Pattern Synthesis Al-gorithm for Arrays[J].IEEE Trans.on Antennas Pro-pagat, 1990, 28 (10) :1 666-1 676.

[9]Applebaum S P.Adaptive Arrays[J].IEEE Trans.on An-tennas Propagat., 1976, 24 (5) :585-598.

[10]束咸荣, 何炳发, 高铁.相控阵雷达天线[M].北京:国防工业出版社, 2007.

[11]Kong J A.电磁波理论[M].吴季, 译.北京:电子工业出版社, 2003.

方向图计算 篇2

MIMO天线方向图的辐射鲁棒性设计

MIMO天线各单元传输的`基带信号可以完全不相关,天线上传输的信号可以是任一幅度和相位的组合,可能产生零辐射现象.MIMO天线零辐射缺陷可以用其单元方向图来判定,对于无损天线,该缺陷还可以用天线本地端口S参数判断,不必测量天线方向图.使用方向图鲁棒性设计方法,可以避免零辐射缺陷.针对MIMO系统基站侧使用阵列天线的情况,本文还给出了阵列天线盲点的判定准则.

作 者：颜罡 杜正伟 王蔷 龚克 YAN Gang DU Zheng-wei WANG Qiang GONG Ke  作者单位：颜罡,YAN Gang(清华大学电子工程系微波与数字通信国家重点实验室,北京,100084;联想研究院,上海,03)

杜正伟,王蔷,龚克,DU Zheng-wei,WANG Qiang,GONG Ke(清华大学电子工程系微波与数字通信国家重点实验室,北京,100084)

刊 名：微波学报  ISTIC PKU英文刊名：JOURNAL OF MICROWAVES 年，卷(期)： 23(6) 分类号：O4 关键词：MIMO   天线   方向图   零辐射   盲点  

测控天线方向图测试系统研究 篇3

天线的辐射方向图简称方向图,是天线的辐射参量随空间方向变化的图形表示。辐射参量包括:辐射的功率通量密度、场强、相位和极化。通常情况下,辐射方向图在远区测定并表示为空间方向坐标的函数[1]。根据天线方向图的图形可以确定天线的波瓣宽度、波束宽度、副瓣电平和方向系数等参数。目前船载测控天线大多采用手动测试,虽然可以通过频谱仪实时打印图形,但仍存在需多人配合和数据处理不便的问题。本文在研究现有测试技术的基础上,研究一种基于GPIB虚拟仪器方向图测试系统,完成对第一副瓣电平、3 dB点波瓣宽度以及10 dB点波瓣宽度的指标的自动化测试和计算,对所采集的数据进行图形化显示,并以可视化报告的形式打印出相关的测量曲线和特征点信息。

1 测试原理

天线在方向(θ,Φ)辐射的电场强度为:

undefined。 (1)

式中,A0为与方向无关的常数;f(θ,Φ)为场强方向函数。由式(1)可得:

undefined。 (2)

实际上常用场强归一化表示方向图,称归一化方向图,表示为:

undefined。 (3)

式中,EM为场强的最大值;F(θ,Φ)为归一化场强方向函数。所以要测归一化场强方向图,只要测得在距离天线的球面上各点的场强即可。

根据天线方向图的图形可以确定天线的波瓣宽度、波束宽度、副瓣电平和方向系数等。基于国军标的雷达方向图测量方法有3大类6种方法,其中源雷达高架的斜式测试场测量方法是在试验场常用的一种方法。

本文在现有测试技术的基础上研究一种基于GPIB及串口的虚拟仪器方向图测试系统,完成对第一副瓣电平、3 dB点波瓣宽度以及10 dB点波瓣宽度的指标的自动化测试和计算,对所采集的数据进行图形化显示,并以可视化报告的形式打印出相关的测量曲线和特征点信息[2,3]。

2 船载雷达方向图测试技术现状

传统的手动方向图测量方法属于静态逐点测试,数据的录取方向图的绘制以及参数的计算都是手工方式,操作复杂,工作量大,完成一次测量需要多人协作,测试与数据处理时间长,且精度很难提高。目前船载测控天线大多采用的半自动化测试,虽然可以通过频谱仪实时打印图形,但仍需多人配合、数据处理不便的问题[4,5]。

下面以船载测控雷达方向图测试系统为例,介绍测试原理和测试方法。测试系统采用源雷达高架斜式测量法。被测天线为接收状态,可在平面内自由转动。源天线为发射状态,固定不动。收、发均为同向线极化或圆极化[6]。

2.1 测试要求

源雷达架设在高塔上,被测雷达及定位装置架设在地面上,源雷达垂直方向图的第一零点指向几何反射点,架设高度应满足远场测试要求以及信标塔有足够的高度,以消除地面反射波的影响[7,8]。

2.2 测试方案及流程

被测天线为接收状态,可在平面内自由转动。源天线为发射状态,固定不动。收、发均为同向线极化或圆极化。在被测天线转动过程中,不断记录可得到该天线的方向图[9]。

① 对准目标:

将源天线安装在标校塔中。按测试要求配置信号,源天线发出测试频率的信号。转动被测天线使其波束中心对准信标塔源天线,微调天线至接收信号最强。

② 设置频谱仪:

设置频谱仪中心频率和带宽。扫描时间和视频带宽视具体情况而定,扫描过程中天线匀速转动。

③ 方位测试:

完成上述步骤后可正式进行测试。被测天线左、右连续偏转3°,在频谱仪上得到±3°范围的包括主瓣和旁瓣在内的方位方向图。

④ 俯仰测试:

使被测天线回到中心点处,并再次精调被测天线使其对准信标塔源天线,以此为零点使天线向上下偏转2°,即在频谱仪上得到±2°范围的包括主瓣和旁瓣的俯仰方向图。

⑤ 数据处理:

根据天线方向图计算出天线增益和其他相关指标[10]。

3 船载雷达方向图自动化测试

3.1 系统概述

为缩短开发时间,提高开发效率,系统采用虚拟仪器中的仪器控制技术。将多个可以完成固定功能的独立仪器组合在一起进行信息交换、协调工作,完成独立仪器的系统集成。独立仪器的系统集成需要2个条件:① 独立仪器自身配置了行业通用的接口或总线的硬件,并且提供与外界通信的代码;② 有一个统一的平台能够与各个相关的独立仪器的接口或总线进行通信。总线技术在系统集成中占有非常重要的作用。传统测试方法与自动化测试方法比较如表1所示。

3.2 系统组成及功能

3.2.1 系统硬件构成

自动化测试系统建立在已有测试系统的基础上。整个系统由笔记本电脑、高速GPIB通信卡、电平转换器、便携式打印机、基于GPIB通用总线接口和RS232接口的虚拟仪器技术测试软件系统组成,如图1所示。

通过虚拟仪器软件将角度控制和测试频谱显示在一个界面同时显示,提高工作效率。系统可以实现多种抛物面雷达和方向图、差方向图、雷达轴比、雷达增益 、雷达差零深、差斜率、差线性度、第一副瓣电平、3 dB点和10 dB点波瓣宽度的指标自动化测试与计算,系统软件分为以下3个模块[11]:

① 被测天线控制模块。通过RS232接口,控制测角监控微机,进而完成对天线转动的实时远程控制。

② 源天线控制模块。通过光纤远程控制标校塔信号源,实现开关机控制、点频设置和电平控制。

③ 测试模块。通过GPIB接口,控制测试仪器。设置测试频率,数据存储路径,测试数据显示处理,完成方向图绘制。

3.2.2 软件实现

雷达天线方向图测试靠测试软件来完成。实现各种雷达天线方向图测试软件的程序流程图如图2所示。整个程序为用户提供了一个智能化、自动化的测试环境和测试界面[12]。

3.2.3 测试方法比较

天线方向图自动化测试系统分别与对应设备的数据通讯后,完成一体化的天线方向图测试。与传统的测试方法相比,系统优点如下:

① 接口编程方便,可实现各种自动标准、多次测量平均等要求,测量精度高;

② 对所采集的数据进行自动分析处理,以图形化显示,以可视化报告的形式打印相关的测量曲线和特征点信息。

③ 便于扩展传统仪器的功能,将多台带有串口及GPIB接口的仪器组合起来,形成较大的自动测试系统可以完成雷达指标的综合测试。高效灵活地完成各种不同的测试任务,而且组建和拆散灵活,使用方便。

4 测试结果分析

以船载测控雷达为测试对象,系统联通后,在方向图自动化测试系统中设置方位偏转角度为-3°~ +3°,俯仰的偏转角度为-2°~+2°,可以得到天线的方向图曲线并及时记录保存。天线俯仰角为90°固定、方位角在中心点(90°)附近正负偏转3°时的天线方向图测试曲线如图3所示[13]。

5 结束语

方向图计算 篇4

！本程序用于图3－10连杆机构运动分析，输出G点的位置、速度、加速度！本TB程序见“PRINT3-1.TRU”

!图3－10连杆机构的动画程序见“图3-10连杆机构动画程序.DOC”！若要绘制从动件的位移、速度、加速度曲线，见“在TB中绘制SVA三曲线.DOC”以及“在ACAD中绘制SVA三曲线.DOC”

！若要将图3－10连杆机构从动件的位移、速度、加速度具体数值打印出来，见“如何将TRUE BASIC的输出数值打印出来.DOC”

OPTION NOLET

FOR I =0 TO 360 STEP 2

CALL LINK(0, 0, 0, 0, 0, 0, I*PI/180, 10, 0, 14, XD, YD, VDX, VDY, ADX, ADY)CALL RRR(XD, YD, VDX, VDY, ADX, ADY, 41, 0, 0, 0, 0, 0, 39, 28, QDA, W3, E3，QBA, W2, E2)

QDC=QDA+35*PI/180 CALL LINK(XD, YD, VDX, VDY, ADX, ADY, QDC, W3, E3, 15, XC, YC, VCX, VCY，ACX, ACY)

CALL RPR(0, 0,-34, 0, 0, 0, 0, XC, YC, VCX, VCY, ACX, ACY, 0, QFG, W4, E4)

CALL LINK(0,-34, 0, 0, 0, 0, QFG, W4, E4, 55, XG, YG, VGX, VGY, AGX, AGY)

PRINT I, QFG*180/PI, W4, E4, XG, YG, VGX, VGY, AGX, AGY NEXT I END

SUB LINK(XA, YA, VAX, VAY, AAX, AAY, QAB, W, E, L, XB, YB, VBX, VBY, ABX, ABY)

XB=XA+L*COS(QAB)

YB=YA+L*SIN(QAB)

VBX=VAX-L*SIN(QAB)*W

VBY=VAY+L*COS(QAB)*W

ABX=AAX-L*COS(QAB)*W^2-L*SIN(QAB)*E

ABY=AAY-L*SIN(QAB)*W^2+L*COS(QAB)*E END SUB

SUB RRR(XA, YA, VAX, VAY, AAX, AAY, XC, YC, VCX, VCY, ACX, ACY, LAB, LCB, QAB，WAB, EAB, QCB, WCB, ECB)

LAC=SQR((XC-XA)^2+(YC-YA)^2)

COSQAC=(XC-XA)/LAC

SINQAC=(YC-YA)/LAC

QAC=ANGLE(COSQAC,SINQAC)

COSQCBA=(LAB^2+LAC^2-LCB^2)/(2*LAB*LAC)

SINQCBA=SQR(1-COSQCBA^2)

QCBA=ANGLE(COSQCBA,SINQCBA)

QAB=QAC-QCBA

XB=XA+LAB*COS(QAB)

YB=YA+LAB*SIN(QAB)

COSQCB=(XB-XC)/LCB

SINQCB=(YB-YC)/LCB

QCB=ANGLE(COSQCB,SINQCB)

WAB=((VAX-VCX)*COSQCB+(VAY-VCY)*SINQCB)/LAB/SIN(QAB-QCB)

WCB=((VAX-VCX)*COS(QAB)+(VAY-VCY)*SIN(QAB))/LCB/SIN(QAB-QCB)

G=AAX-ACX-LAB*COS(QAB)*WAB^2+LCB*COSQCB*WCB^2

F=AAY-ACY-LAB*SIN(QAB)*WAB^2+LCB*SINQCB*WCB^2

EAB=(G*COSQCB+F*SINQCB)/LAB/SIN(QAB-QCB)

ECB=(G*COS(QAB)+F*SIN(QAB))/LCB/SIN(QAB-QCB)END SUB

SUB RPR(M, XA, YA, VAX, VAY, AAX, AAY, XC, YC, VCX, VCY,ACX,ACY, LAB,QBD,W,E)

LAC=SQR((XC-XA)^2+(YC-YA)^2)

COSQAC=(XC-XA)/LAC

SINQAC=(YC-YA)/LAC

QAC=ANGLE(COSQAC,SINQAC)

LBC=SQR(LAC^2-LAB^2)

QACB=ATN(LAB/LBC)

QBD=QAC+M*QACB

DELTA=-(YC-YA)*SIN(QBD)-(XC-XA)*COS(QBD)

DELTAW=(VCX-VAX)*SIN(QBD)-(VCY-VAY)*COS(QBD)

DELTAV=-(YC-YA)*(VCY-VAY)-(XC-XA)*(VCX-VAX)

W=DELTAW/DELTA

VLBC=DELTAV/DELTA

T1=(ACX-AAX)+(VCY-VAY)*W+SIN(QBD)*W*VLBC

T2=(ACY-AAY)-(VCX-VAX)*W-COS(QBD)*W*VLBC

DELTAE=T1*SIN(QBD)-T2*COS(QBD)