风螺旋公切线算法详解
好久不见,近来一切可好?2017年最后这段时间里,狂补了一把C#,希望未来能够从软件代码层面实现风螺旋算法的验证与推广。今天跟大家分享的这个话题的底图就是最近一段时间的学习成果:一个基于WPF架构的非常简单的绘图框架,以及对风螺旋的自动化绘制进行的实现。闲话少叙,开始今天的主题。
在掌握了风螺旋切线计算的基础上,就可以开始公切线算法的研究了。公切线的计算是飞行程序模板中非常关键的一项内容,因此,在开始模板算法分享之前,详细回顾一下公切线的算法是非常有必要的。
风螺旋公切线计算首先是基于以下4个基础设定:
/// 1. 同一程序模块内,转弯基础参数相同。无论是程序转弯、还是等待模板,在同一个ICAO标准模板中,转弯参数是相同的,即转弯半径是相同的。
/// 2. 风螺旋的切线与Esita方向相垂直。之前从原理及数学软件上对此项内容进行了证明,这是公切线的一个基础条件。
/// 3. Esita方向与对应的基准圆半径之间的夹角恒等于偏流角DA。ICAO DOC8168中的配图明确标注了此项内容DA= arcsin(w/v)。
/// 4. Esita的大小与sita角度可以互相换算。距离可以换算成旋转角度,反之旋转的角度可以换算为对应的距离。这是本公众号中多次提及的一个基本概念。风螺旋是假想的航空器的飞行轨迹,对这个轨迹整体进行旋转时,它的外扩规律符合按照Esita方向等距离外扩的特性,外扩的大小与旋转过的角度有直接的换算关系。
除以上4个基础设定以外,公切线计算中按照个人偏好使用了以下的标识来表示风螺旋的参数(如图1所示):
图1 风螺旋参数示意图
rotation: 从系统的零度旋转至风螺旋起点的角度。航空地图中采用的是磁北为零度,计算机系统中多数以水平向右为零度的起点,顺时针角度增大为默认条件。
offset:表示外扩距离的一个量。对于一个风螺旋来说,外扩和旋转角度可以互相换算,因此,最佳方案是对每一个风螺旋都设定一个初始的系统角度(rotation)和一个初始的外扩距离(offset)。
sita:用来表示从风螺旋起点开始实际的转弯角度。通过该角度可以换算出从起点开始一共外扩了多少距离,或者用来表示风螺旋中的某个具体位置点。
Esita:与sita角相对应的外扩距离。这个Esita的方向与对应半径的延长线之间的夹角为偏流角DA,与该位置处风螺旋切线的夹角为90度。
为了更清晰的表示计算关系,对图1中的关键点进行编号得到图2的内容:
图2 风螺旋公切线的角度关系
图2中绘制了两条实心的风螺旋线ws1和ws2,以及与它们对应的无外扩的虚线风螺旋。虚线风螺旋与圆周的交点视为风螺旋的实际起始角度。旋转与外扩可以相互换算,因此可知,实线风螺旋与虚线风螺旋之间是“等距平行”关系。
由于Esita与风螺旋的切线相垂直(基础设定2),因此,若两个风螺旋之间存在有公切线,则必定分别有两个Esita与公切线相垂直,如图1中的Esita1和Esita2,它们二者为平行关系。
由于两个Esita与对应的半径Radius之间夹角固定为偏流角DA(基础设定3),因此可知对应半径之间为平行关系。因为,两条风螺旋的基础参数是相同的,因此,两个半径是相等的。于是图2中的abcd是一个平行四边形。从d点向fc边做垂线,得到垂足e点,可知defg是一个矩形。
ws2风螺旋在公切点位置处的Esita可以分成两部分,EsitaA和EsitaB,前者与ws1中的外扩距离是相等的,后者产生的原因是两条风螺旋初始条件的差异所造成。若每弧度可以外扩的距离为Edist,则EsitaB可表示为offset2 - offset1 + (rotation1-rotation2) * Edist,即初始外扩距离之差,加上初始角度之差所换算成的外扩距离,就构成了图2中的EsitaB。图2中的EsitaA可视为两条风螺旋共同旋转过的角度所形成的距离,因为外扩的效率是相同的,相同的外扩角度可以得到相同的外扩距离。rotation1-rotation2这个公式是浪费了两天多的时间才最终纠错纠出来的,之前在Flash系统中不存在此问题,C#中的绘图系统真是复杂许多。rotation1-rotation2代表了以rotation2为基准角度,初始的角度差值。
具体的计算步骤如下:
1、放置风螺旋,计算dc线的角度。
当我们向屏幕中放置一条风螺旋时,它的中心点位置就可视为一个已知条件。设两个中心点坐标分别为c1(x,y),c2(x,y),那么,线段ab在软件系统中的角度为: arctan((c2.y-c1.y)/(c2.x-c1.x)),这里需要用到反正切函数,获得一个角度值,称之为基准线角度,用baseLineAngle来表示,圆心点间的距离用distance来表示。dc线与ab线平行且相等,因此,baseLineAngle就是dc线的方向。
2、计算dg线的角度。
图2中⊿dce为一个直角三角形,∠cde的值用angle来表示,它等于arcsin(EsitaB/distance)。图2中逆时针方向角度减小,因此可将线段dg的角度表示为:baseLineAngle – angle - Pi/2 (所有角度值均用弧度来表示),公切线的方向就等于baseLineAngle – angle。
若ws1的外扩距离大于ws2时,angle角的值为负值,位置关系如图3所示:
图3 angle角为负值的情况
此时,∠cdg的值可表示为Pi/2 + angle(因为angle小于零),线段dg的角度可以表示为baseLineAngle -(Pi/2 + angle),化简后,得到的公式为baseLineAngle – angle- Pi/2。与前面讨论的angle大于零时的公式是完全一致的,也就是说不论ws1或ws2外扩距离谁大,dg线的角度始终可以用baseLineAngle – angle - Pi/2来表示,公切线的角度可以用baseLineAngle– angle来表示。
3、计算dg线的长度。
从dg线的角度反推对应的半径的方向,可采用公式:baseLineAngle - angle - Pi/2 + DA(因为Esita与半径的夹角总等于DA)。得到圆周半径方向之后,分别与各自的rotation参数相减,就可以得到风螺旋真实的旋转角度,用这个旋转角度sita换算成外扩距离,再加上初始的外扩距离,就可以得到dg线的长度,角度、长度都知道以后,可以得到准确的切点位置。连接两个切点的坐标,得到的线段即是两条风螺旋线的公切线。
公切线的计算过程就是这样,您get到了没?
另外再附两张公切线的示意图,因为手工条件下真是很难想像公切线会出现在什么地方(以后这个问题将不是问题 :-))
这两张图与前面的图的区别在于:初始的 rotation 大小关系对调了一下,结论完全相同。这就是自动化计算的魅力所在,哈哈哈。
可以预见的未来中,风螺旋线将成为一门新兴的学科,它将与阿基米德螺旋、渐开线螺旋共同纳入 等距螺旋 的框架,成为几何课程的一部分,风螺旋相关计算将成为飞行程序设计入门课程的一部分。本公众号将以此为目标持续努力,全面开源共享。期待在不远的将来风螺旋线精确计算方法将成为DOC8168的一部分,并且成为飞行程序领域送给 几何研究领域的一份礼物!
祝大家:新年快乐!
![](https://i-blog.csdnimg.cn/blog_migrate/33c876aa47b2ab9120200861eb974dab.png)
微信扫一扫
关注该公众号