CN118050814A - 低空光学定姿三分量磁测***及测定磁场矢量数据的方法 - Google Patents
低空光学定姿三分量磁测***及测定磁场矢量数据的方法 Download PDFInfo
- Publication number
- CN118050814A CN118050814A CN202410450963.5A CN202410450963A CN118050814A CN 118050814 A CN118050814 A CN 118050814A CN 202410450963 A CN202410450963 A CN 202410450963A CN 118050814 A CN118050814 A CN 118050814A
- Authority
- CN
- China
- Prior art keywords
- axis
- theta
- data
- camera
- magnetic
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000005259 measurement Methods 0.000 title claims abstract description 88
- 238000000034 method Methods 0.000 title claims abstract description 52
- 230000003287 optical effect Effects 0.000 title claims abstract description 51
- 230000002829 reductive effect Effects 0.000 claims abstract description 27
- 238000013016 damping Methods 0.000 claims description 38
- 238000004364 calculation method Methods 0.000 claims description 37
- 238000012937 correction Methods 0.000 claims description 32
- 230000005540 biological transmission Effects 0.000 claims description 18
- 230000033001 locomotion Effects 0.000 claims description 8
- 238000006243 chemical reaction Methods 0.000 claims description 5
- 238000006073 displacement reaction Methods 0.000 claims description 5
- 238000012935 Averaging Methods 0.000 claims 1
- 230000001131 transforming effect Effects 0.000 claims 1
- 239000000969 carrier Substances 0.000 abstract description 4
- 238000005265 energy consumption Methods 0.000 abstract description 2
- 230000000694 effects Effects 0.000 description 13
- 230000006870 function Effects 0.000 description 12
- 230000005358 geomagnetic field Effects 0.000 description 10
- 230000000007 visual effect Effects 0.000 description 10
- 238000005516 engineering process Methods 0.000 description 9
- 238000012360 testing method Methods 0.000 description 7
- 238000010586 diagram Methods 0.000 description 6
- 238000011161 development Methods 0.000 description 4
- 238000005070 sampling Methods 0.000 description 4
- 239000004065 semiconductor Substances 0.000 description 4
- 238000004891 communication Methods 0.000 description 3
- 239000000463 material Substances 0.000 description 3
- 238000004377 microelectronic Methods 0.000 description 3
- 230000001681 protective effect Effects 0.000 description 3
- 230000003068 static effect Effects 0.000 description 3
- XEEYBQQBJWHFJM-UHFFFAOYSA-N Iron Chemical compound [Fe] XEEYBQQBJWHFJM-UHFFFAOYSA-N 0.000 description 2
- 230000002159 abnormal effect Effects 0.000 description 2
- 238000012512 characterization method Methods 0.000 description 2
- 230000007423 decrease Effects 0.000 description 2
- 230000009977 dual effect Effects 0.000 description 2
- 230000008030 elimination Effects 0.000 description 2
- 238000003379 elimination reaction Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 239000012530 fluid Substances 0.000 description 2
- 238000011835 investigation Methods 0.000 description 2
- 230000002035 prolonged effect Effects 0.000 description 2
- 230000000717 retained effect Effects 0.000 description 2
- 239000000523 sample Substances 0.000 description 2
- 206010034719 Personality change Diseases 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 230000004313 glare Effects 0.000 description 1
- 238000002329 infrared spectrum Methods 0.000 description 1
- 229910052500 inorganic mineral Inorganic materials 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 229910052742 iron Inorganic materials 0.000 description 1
- 230000005389 magnetism Effects 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- 239000011707 mineral Substances 0.000 description 1
- 238000002360 preparation method Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000000758 substrate Substances 0.000 description 1
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Length Measuring Devices By Optical Means (AREA)
Abstract
本发明公开了一种低空光学定姿三分量磁测***及测定磁场矢量数据的方法,属于地球物理勘探领域,在2台小型低空航空载体上分别安装定位激光光源,在另一航空载体上安装微型图像传感器的相机组件与磁场矢量测量***,使用相机瞬时观测到的图像数据,解算相机本身的空间姿态,与相机同一载体的磁场矢量测量元件也就获得了准确可控的姿态信息,最终解算真地理坐标环境下的地磁矢量数据。本申请规避了昂贵惯导器件的定姿方式,***成本较低,体积较小,空中负荷轻,耗能较小,适用于较长时间的空中飞行;可以靠近地面测量,得到地下一定深度地质体的明确地磁信息;大幅降低了对磁矢量元件的电磁干扰,测量误差优于7.5nT。
Description
技术领域
本发明涉及地球物理勘探领域,具体地讲是涉及一种利用光学方法进行航空载体姿态确定进行大地磁场矢量测量的低空光学定姿三分量磁测***及测定磁场矢量数据的方法。
背景技术
磁三分量测量是一种测量磁场的方法,它能够提供磁场三个分量的信息。这种测量方法在多个领域有着广泛的应用,例如地质调查、矿产资源勘探、海洋地磁研究等。
在航空地球物理勘探中,航磁三分量测量技术能够获得包括磁偏角、磁倾角、垂直及水平分量在内的地磁场全要素信息,这对于对磁异常进行定性分析和定量解释非常有利。在地球物理勘探领域中,航空磁测技术目前以地磁场总场强度或总场梯度测量为主。与总场测量相比,矢量测量可同时获取地磁场模量大小和方向信息,有效减少反演中的多解性,有助于对磁性体的定量解释,提高地下矿体探测分辨率和定位精度,使得矢量磁测成为航空磁测的主要发展方向之一。
磁三分量为矢量数据,需要确定测试设备自身的姿态方位,根据设备自身的姿态方位和测量得到的磁三分量数据,解算得到真地理坐标环境下的磁分量数值。磁测仪器在运动中的姿态控制,成为该项技术的关键点。由于运动中的高精度姿态控制的技术困难,制约了地球物理勘探中磁场矢量测量技术的发展与应用。
在地面磁三分量测量中,传统的方法是使用经纬仪进行测量基线定位,并将磁矢量测试元件安装于无磁经纬仪器之上,在静止状态下逐个测量基线上每个测点的磁矢量数据。该方法定位精度较高,但属于静态观测,速度较慢,效率较低。
已有的卫星导航长基线定姿方式,是通过解算安置在运动载体上的多条天线间基线矢量信息,再通过基线矢量求解旋转矩阵来获得载体的三维姿态。精度尚可,性能稳定。由于卫星导航中的卫星天线需要一定距离才可提高定姿精度,该技术一般适用于大型船只,不适用于较小尺寸的航空载体。
微型惯性器件,即IMU,已在民用技术领域大量应用,如可用于无人机、AR/VR设备以及光学图像防抖等,IMU为主的组合传感器可做到微型芯片化,体积重量较小,但该类器件不能达到磁场矢量测量的定姿精度要求。
目前确实需要一种能够快速高精度实现大地磁场矢量测试的设备和方法。
发明内容
本发明的目的在于克服上述传统技术的不足之处,提供一种融合光学、卫星导航通讯,机器视觉、微电子、空间姿态控制,磁矢量元件、软件等技术,达到快速高精度实现大地磁场矢量测试功能的低空光学定姿三分量磁测***。
本发明的目的是通过以下技术措施来达到的:
低空光学定姿三分量磁测***,其特征在于:包括地面主机和若干台航空机,所述地面主机以无线传输方式给航空机发送轨迹指引、测量参数、卫星导航基站差分数据,所述航空机包括飞行器、主机、卫星导航模块、无线传输模块、电源、VCSEL光带环、阻尼杆、若干台相机、IMU模块、倾角传感器和磁分量传感器,所述主机包括处理器和存储器,所述卫星导航模块包括卫星导航信号天线,所述IMU模块包括微型惯性器件,所述VCSEL光带环由VCSEL光点组成,所述卫星导航信号天线安装于飞行器上端,所述主机和电源安置于飞行器上,所述VCSEL光带环安置于飞行器下部,所述阻尼杆位于VCSEL光带环下部,所述相机和微型惯性器件通过安装平台安置于阻尼杆下部,所述磁分量传感器通过刚性杆安装于阻尼杆下部,所述磁分量传感器上安装有温度传感器,所述主机与卫星导航模块、无线传输模块、VCSEL光带环、相机、IMU模块、温度传感器和磁分量传感器电连接。磁分量传感器为磁三分量磁通门元件,包括磁场矢量X元件、磁场矢量Y元件和磁场矢量Z元件,分别获得X、Y、Z分量数据。
地面主机实现地面台站与观测磁场日变功能。地面主机可以采用航空机相同结构,地面主机的区别在于可以使用附加电源与较稳定的固定支架。相同结构的地面主机和航空机可以简化设备类型,现场工作时可以随机配用,且使用附加电源可以延长测量作业时间。
本申请中选择所述飞行器为旋翼式飞行器,但是选择其他类型飞行器也同样适用。所述飞行器机身涂装有与背景亮度接近的粗糙近红外光线反射材料。
阻尼杆为具阻尼作用的连接杆,其长度取0.8-1.0m,阻尼杆具有三个功能:1.适当加长长度,减少主机与旋翼式飞行器对磁三分量元件的电磁干扰;2.具有阻尼作用,使旋翼式飞行器的震动传递到下部的震动减弱减小,特别是大幅减少高频分量;3.起飞前支撑。阻尼是指任何振动***在振动中,由于外界作用,如流体阻力、摩擦力等,和/或***本身固有的原因引起的振动幅度逐渐下降的特性,以及此一特性的量化表征。阻尼的力学模型在物理学和工程学上应用广泛,根据现有技术即可构造具有阻尼作用的连接杆,不赘述。
相机与磁分量传感器之间的刚性杆长度约0.4m。通过刚性的刚性杆连接相机、微型惯性器件和磁分量传感器,既可以使相机与磁分量传感器间保持适当距离,尽可能减少相机对X、Y、Z元件的电磁干扰,还可以使三者间避免弹性震动误差,提高精度。但是在航机倾斜时,VCSEL光点处位置与卫星导航信号处差别较小,但相机位置通过阻尼杆连接因此影响较大,可以显著缩短阻尼杆,甚至阻尼杆缩短为零,通过适当加长刚性杆以减少航机对磁矢量X、Y、Z元件的干扰,同时减少相机与卫星导航信号误差。
飞行器具有在磁分量传感器之下设置的三条支撑杆,用于起飞前在地面起支撑作用。在飞行器旋翼外设简易防护条,防止降落时倾倒损坏。如进一步地减少干扰,将阻尼杆加长至1.0m或以上时,应在阻尼杆最上部设置铰链,在起飞、降落时将相机、磁分量探头、刚性杆等元件收起,在主机处另外设置起降支撑杆。
若干台所述相机环绕围绕安装平台均匀排布。所述相机数量为6台或者8台。相机采用6台时,对应相机为水平视角64°近红外感光相机,相机采用8台时为水平视角49°近红外感光相机;相机像素可用1280×720像素,尺寸适中,既能满足测量精度需要,又可减少数据计算,提高效率。相机采用固定光圈镜头,焦距固定在无穷远处。相机镜头与感光元件之间安装半带宽15nm的窄带滤光片,滤光片中心波长取943-945nm。在光源为频谱中心为940nm时,由于镜头与感光元件之间最边缘处的光线与主轴夹角可能接近10°,使得窄带滤光片的窄带滤光特征改变,峰值向左偏移低于940nm,且VCSEL光点光源在使用较长时间后,温度上升,会使光线光谱中心右移,滤光片中心波长取943-945nm最优。
本申请中镜头的后焦加长,使感光元件之前接收到的光线尽可能接近垂直方向地通过窄带滤光片,兼采用较小半带宽的滤光片,进一步地减少杂散光干扰。
所述相机帧率为20fps以上,快门为全局快门,以适应动态下的拍摄,快门速度低于1/1000s,快门速度与每一帧图像的拍摄可由外部信号控制。
所述VCSEL光带环为VCSEL球笼型光带环,VCSEL球笼型光带环由三排密集均匀排布的VCSEL光点组成,三排所述VCSEL光点中的一排位于球笼型的中间位置,另两排所述VCSEL光点对称分布在球笼型的上下位置且与中间位置的VCSEL光点位置交错设置。所述VCSEL光点即垂直腔面发射激光器,也可以称为VCSEL半导体激光器,VCSEL光线中心波长940nm。本***使用波长940nm的近红外激光作为标志点光源,在地表或低空环境中,940nm的自然背景光强相对较弱,从而减少杂散光干扰。以地球仪为例,赤道位置排布36只以上VCSEL光点,在北纬和南纬6-8°线上排布36只以上VCSEL光点,并与赤道位置光点错开排布,每只VCSEL光点功率3-5w,由主机通过脉冲电源控制光点,一般接通时间2ms以下,整个VCSEL光带环的平均功率4-6w。VCSEL是一种在与半导体外延片垂直方向上形成光学谐振腔、发出的激光束与衬底表面垂直的半导体激光器结构。光点即为VCSEL半导体激光器发射激光形成的发光点。
磁分量传感器之上安装无磁的温度传感器,监控元件温度变化。
主机包括处理器和存储器,所述处理器接收卫星导航模块、IMU模块、倾角传感器和磁分量传感器的电信号,所述处理器通过无线传输模块与地面主机通讯,所述处理器通过脉冲电源控制VCSEL光带环中的VCSEL光点点亮,所述处理器通过相机拍摄照片并解算两帧照片的帧间差分情况,使用当时IMU模块的数据复算即时抖动量,优化帧间差分计算结果,帧间差分计算后利用灰度质心法求解光点像素坐标,处理器与磁分量传感器得到的X、Y、Z分量数据进行真地理坐标转换,得到设定地理坐标环境下的磁场矢量数据,将该磁场矢量数据与位置、高度信息一起经无线传输模块传输给地面主机并转入存储器。
可以采用双处理器和存储器,以加快解算速度,具体可以采用:主机包括处理器一、存储器一、处理器二和存储器二,所述处理器一接收卫星导航模块、IMU模块、倾角传感器和磁分量传感器的电信号,所述处理器一通过无线传输模块与地面主机通讯,所述处理器一通过脉冲电源控制VCSEL光带环中的VCSEL光点点亮并将相机拍摄照片送入处理器二,所述处理器二解算两帧照片的差分情况,并使用当时IMU模块的IMU数据复算即时抖动量,优化差分结果,差分后利用灰度质心法求解像素坐标,处理器二解算得到的像素坐标回传给处理器一,处理器一与磁分量传感器得到的X、Y、Z分量数据进行真地理坐标转换,得到设定地理坐标环境下的磁场矢量数据,将该磁场矢量数据与位置、高度信息一起经无线传输模块传输给地面主机并转入存储器一。
所述处理器二将每一帧照片暂存于存储器二,用于停飞后取出复算核对。
一种地磁矢量数据测定方法,其特征在于:在至少2台定位航空载体上分别安装定位激光光源和卫星定位装置,在测量航空载体上安装相机组件、卫星定位装置和磁场矢量测量***,使用相机瞬时观测到的定位航空载体上定位激光光源图像数据和该定位航空载体上卫星定位装置的空间坐标、结合测量航空载体自身卫星定位装置的空间坐标解算得到相机本身的空间姿态,该测量航空载体上的磁场矢量测量元件获得磁场矢量三分量数据,结合空间姿态信息和磁场矢量三分量数据解算得到真地理坐标环境下的地磁矢量数据。
GPS定位中的RTK方式可以实时提供高精度的位置信息。
本申请的低空光学定姿三分量磁测***综合利用了以上技术发展,融合了光学、GPS通讯,机器视觉、微电子、空间姿态控制,磁矢量元件、软件等技术,达到快速高精度实现大地磁场矢量测试的功能。
应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,其特征在于:采用上述的低空光学定姿三分量磁测***通过以下步骤测定磁场矢量数据:
S1:航空机主机将时间信号传递给相机执行拍摄,同时给予远端航空机信号,使远端航空机的VCSEL光点在相机拍摄第一帧照片时同步点亮,并同步采集一次磁分量传感器的磁三分量数据;
S2:按照设定帧率,远端航空机VCSEL光点关停,同步拍摄第二帧图像,将第二帧图像参照姿态数据进行姿态改正后,与第一帧图像进行帧间差分计算,得到光点像素坐标;
S3:航空机主机***程序执行中多台相机同步工作,主机对各相机拍摄到的光点像素坐标数据进行判别,对于相机像素边缘部位的光点位置与相邻相机拍摄得到的光点位置进行比对,舍去靠近边缘的像素点,并将该位置信息反馈给处理器,控制航空机调整飞行轨迹,使相机得到的像素光点位置远离像素框边缘,提高解算精度;
S4:初步判定光点像素坐标位置后,航空机保持相对稳定飞行状态时,锁定一个像素范围,减少下一步解算帧间差分图像的计算范围;
S5:得到光点像素坐标后,解算球坐标环境下的θ、ψ数据,结合提供VCSEL光点的远端航空机得到的θ′、ψ′共同解算本机的姿态数据,引入前一时刻测量得到的磁三分量数据,解算真地理坐标环境下的磁分量数值,即得磁场矢量数据,将磁场矢量数据与本机位置、高程信息进行组合后加以存储,并/或发送到地面主机或地面台站。
其中,帧间差分计算方法如下:
使VCSEL光点点亮时相机拍摄第一帧照片,其像素灰度表示为P1(i,j),然后VCSEL光点关停,相机拍摄第二帧照片,其像素灰度表示为P2(i,j);
其中,P1(i,j)是VCSEL光点点亮时相机拍摄第一帧图像;P2(i,j)是VCSEL光点关停时相机拍摄第二帧图像;P表示像素灰度数值,i表示左右方向的像素坐标,由左至右数值增大,i取值范围为1-n;j表示上下方向的像素坐标,由上至下数值增大,j取值范围为1-m;
在P2(i,j)中,根据IMU模块中的微型惯性器件给出的第一帧与第二帧之间的姿态变动数据,分别求相机平台在水平、垂直方向位移与转动后像素的新位置:
上述算式中,noij是受相机水平抖动与自转影响、水平变动后图像像素坐标错动量,计算数值精度为小数点后3位以上;moij是受相机竖直抖动与自转影响、垂直变动后图像像素坐标错动量,计算数值精度为小数点后3位以上;Plev是水平抖动变化量,计算数值精度为小数点后3位以上;Pver是垂直抖动变化量,计算数值精度为小数点后3位以上;θrot是相机沿镜头轴线自转抖动的角度,数值范围一般为0.2°以下;
得出应有的校正量:
上述算式中,nij是为消除相机水平抖动与自转影响、水平变动后图像像素坐标修正量,计算数值精度为小数点后3位以上;mij是为消除相机竖直抖动与自转影响、垂直变动后图像像素坐标修正量,计算数值精度为小数后3位以上;
由于两帧之间时间间隔较小,其校正量是一个微小的数值,一般不超过相当的3-5个像素角度,所以由于像素自身转动引起每一个像素方格倾斜量的校正可忽略,只进行像素点垂直与水平位移的影响计算。
将P2(i,j)的图像像素灰度投影到新的坐标图像中,得到P3(i,j),然后与第一帧图像进行直接差分,得到差分图像P4(i,j):
上述算式中,P3(i,j)是对第二帧P2(i,j)图像进行抖动校正后的图像;P4(i,j)是差分后的光点灰度图像。
在差分图像中,相似的背景数据减差后基本为0值,光点数据就得到了明显的突出表示,但噪音数据不能完全消除。
由于噪音影响或较大灰度值的异常像素点可能存在,不能使用单点的像素灰度值搜索最大亮度点Pmax的位置,可采用以下方法消除噪音数据以获得准确像素位置:
使用9个网格数据组合形成较大的搜索框P4max(i,j),并降低最大灰度值标准为Psea:
上述算式中,P5(i,j)是缩小搜索范围的光点灰度图像;Pmax是在P4(i,j)中单点的像素灰度中最大的灰度数值;Psea是在P4(i,j)中进行组合框搜索时使用的降低标准的灰度数值;KL是降低最大灰度标准的系数,一般取值0.4-0.7;P4max(i,j)是组合框搜索范围;
对整个图像P4(i,j)进行搜索,当检索到P4max(i,j)大于Psea时,就立刻记录所在像素点的位置,记为I1、J1,不再进行搜索。当然,在全部像素位置中,可能存在多个像素点满足搜索要求,本步骤只是按照约定搜索方向找到第一个符合条件的像素点。
设定亮点像素周边范围,设置应用半径Pr,在P4(i,j)中,进行如下判别:
上述算式中,I1是第一次组合框搜索得到的亮点中心的水平像素坐标;J1是第一次组合框搜索得到的亮点中心的竖直像素坐标;Pr是设定亮点像素周边范围的应用半径;
如上式成立,则将所在像素点灰度数据恢复为0值,形成一个围绕光点中心的圆形范围。再对圆形范围内的背景噪音进行判别:
上述算式中,PK是背景噪音消除限定值,一般取3-5;
如上式成立,则将圆形范围所在像素点灰度数据恢复为0值。如上式不成立,则保留原亮点数据,得到缩小范围的亮点图像P5(i,j),像素坐标值仍保持不变。
在图像P5(i,j)中,使用灰度质量中心法求解亮度中心点的像素坐标:
上述算式中,I2是图像亮度中心点左右方向的坐标;J2是图像亮度中心点竖直方向的坐标;
将第一次灰度质心法得到的I2、J2数据再次代入到原I1、J1位置,重新界定一个新的数据缩小范围,使用相同算法再次进行灰度质心求解,得到像素坐标I3、J3,就可提高亮点中心的定位精度,获得定位精度可达0.1倍图像像素尺度的高精度的光点中心位置,也即像素位置。除灰度质心法外,也可使用其他方法求取光点中心位置。I3是第二次计算得到的图像亮度中心点左右方向的精确坐标;J3是第二次计算得到的图像亮度中心点竖直方向的精确坐标。
将各个相机上得到的高精度的像素坐标数值按照校正公式转化为球坐标体系下的视角数据θ、ψ,就可与卫星导航定位数据一起进行定姿计算。
坐标体系运动定位计算方法如下:
建立三维坐标系,并同时建立球坐标系,设球坐标模值为1构成圆球面,倾角θ以三维坐标的X、Y平面为零点,方位角ψ以三维坐标的X、Z平面为零点。1光点、2光点表示相机与卫星导航控制的两台远端航空机光点方位。在相机控制下的1光点、2光点的方位与磁矢量X、Y、Z元件方向同时固定在一个姿态之中,在1光点、2光点受到卫星导航控制产生方位变化时,磁矢量X、Y、Z元件方向随之同步变化。磁矢量X、Y、Z元件即磁分量传感器。
在从球坐标系圆心到X轴方向上,通过1光点与X轴垂直的平面与圆球面相交构成X1相交圆,X1相交圆的半径为:
上述算式中,θ1是相机像素计算的1光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ1是相机像素计算的1光点的方位角弧度数据,数值范围介于0°至360°。
视线方向沿X轴从外向内观测,以X轴左侧为起点1光点与该X1相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到Y轴方向上,通过1光点与Y轴垂直的平面与圆球面相交构成Y1相交圆,Y1相交圆的半径为:
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,1光点与Y1相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到X轴方向上,通过2光点与X轴垂直的平面与圆球面相交构成X2相交圆,X2相交圆的半径为:
上述算式中,θ2是相机像素计算的2光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ2是相机像素计算的2光点的方位角弧度数据,数值范围介于0°至360°;
视线方向沿X轴从外向内观测,以X轴左侧为起点2光点与该X2相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到Y轴方向上,通过2光点与Y轴垂直的平面与圆球面相交构成Y2相交圆,Y2相交圆的半径为:
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,2光点与Y2相交圆圆心的连线与基准水平面的夹角为:
原始数据θ1、ψ1和θ2、ψ2与运动后的θ1'、ψ1'和θ2'、ψ2'共8个观测数据,设θx、θy为整个平台分别围绕X轴、Y轴转动的角度,是未知量,由θ1、ψ1和θ2、ψ2可得出θrx1、rx1、θry1、ry1、θrx2、rx2、θry2、ry2的解析表达式,建立二元方程组:
上述算式中,θ1'是卫星导航方式计算的1光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ1'是卫星导航方式计算的1光点的方位角弧度数据,数值范围介于0°至360°;θ2'是卫星导航方式计算的2光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ2'是卫星导航方式计算的2光点的方位角弧度数据,数值范围介于0°至360°;rx1是过1光点与X轴垂直的平面与球坐标圆球面的交线的圆的半径,rx1取值范围为0-1;θrx1是圆球面上1光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度,数值范围介于﹣180°至﹢180°;ry1是过1光点与Y轴垂直的平面与球坐标圆球面的交线的圆的半径,取值范围为0-1;θry1是圆球面上1光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度,数值范围介于﹣180°至﹢180°;rx2是过2光点与X轴垂直的平面与球坐标圆球面的交线的圆的半径,rx2取值范围为0-1;θrx2是圆球面上2光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度,数值范围介于﹣180°至﹢180°;ry2是过2光点与Y轴垂直的平面与球坐标圆球面的交线的圆的半径,取值范围为0-1;θry2圆球面上2光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度,数值范围介于﹣180°至﹢180°;θx是X轴的假设旋转角度,数值范围介于﹣30°至﹢30°;θy是Y轴的假设旋转角度,数值范围介于﹣30°至﹢30°;上式中所述球坐标圆球面的半径均为1。
该式求解可得到θx、θy。由于三角函数内都有θx、θy,左侧分别为θrx1、θrx2、θry1、θry2等,不能视为同一变量,所以该方程组不能用线性法直接求解。应适当变换,将等式左侧强制分解,使θrx1、θrx2、θry1、θry2固定,等式右侧只取变化量。将不同参数对θx、θy的影响变化取系数,就可将θx、θy作为独立项列出:
由于上式左侧:
可以认为如下形式成立:
其中a'、b'为系数,a'受到θrx1、θx影响,b'受到θry1、θy影响。
引入如下系数:
Kx、Ky为θx、θy变化时,弧长与弦长的比例关系。Krx1、Krx2为θx变化的中点处,受θx变化影响,sinθx的变化幅度系数,Kry1、Kry2为θy变化的中点处,受θy变化影响,sinθy的变化幅度系数:
上述算式中,Kx是当X轴旋转θx时其弦长与弧长的对比系数,在θx取值范围﹣30°至﹢30°时,计算精度可达10﹣5度;Ky是当Y轴旋转θy时其弦长与弧长的对比系数,在θy取值范围为﹣30°至﹢30°时,计算精度可达10﹣5度;Krx1是在X轴旋转时1光点在初始与第二次观测之间的角度位置系数;Kry1是在Y轴旋转时1光点在初始与第二次观测之间的角度位置系数;Krx2是在X轴旋转时2光点在初始与第二次观测之间的角度位置系数;Kry2是在Y轴旋转时2光点在初始与第二次观测之间的角度位置系数;
令:
上式中a、b、c、d、e、f是二元方程组参数;
建立方程组:
行列式为:
解析解为:
求解θx、θy:第一步,首先将上式a--d中含有的未知θx、θy置零,求解第一次结果,将得到θx、θy值再置入方程组,进行第二次迭代,所得θx、θy值再置入方程组,进行第三次迭代,以此方法重复求解到4-20步时,直至得到的θx、θy精度达到或者超过千分之一度,满足观测精度要求。
在θ1、θ2相差较小时,ψ1、ψ2处于非0°或非180°时,可在4步求出较高精度结果;在θ1、θ2方向同时变化时,在ψ1、ψ2相差180°左右,如160°或200°时,须20步才可求得百分之一度的精度的结果,在150°-210°时,20步可求得千分之一精度的结果;在θ1、θ2反向变化时,在ψ1、ψ2相差0°左右,如20°、-20°时,须20步达到百分之一;30°或-30°时,20步可达千分之一度的结果。
在ψ1、ψ2夹角为0°或180°时无法定姿,方程式无意义,实际应用时,ψ1、ψ2夹角一般大于30°,小于150°。
求得θx、θy后,设***平台初始的绕Z轴转动量为ψ0,ψ0与ψxy1、ψxy2为独立关系,共同叠加后形成观测量ψ1'和ψ2',可利用正算式求出1光点和2光点由θx、θy转动引起的绕Z轴相应的微小转动量ψxy1和ψxy2:
上述算式中,X1'是1光点受θx、θy影响后在X方向的新位置;Y1'是1光点受θx、θy影响后在Y方向的新位置;X2'是2光点受θx、θy影响后在X方向的新位置;Y2'是2光点受θx、θy影响后在Y方向的新位置;ψxy1是1光点由θx、θy转动引起的绕Z轴相应的微小转动弧度量,数值范围为0°至5°;ψxy2是2光点由θx、θy转动引起的绕Z轴相应的微小转动弧度量,数值范围为0°至5°;ψxy1、ψxy2一般在5°以内。
由于:
故:
上述算式中,ψ01是1光点数据计算得到的水平转动角度,数值范围介于0°至360°;ψ02是2光点数据计算得到的水平转动角度,数值范围介于0°至360°。
ψ01与ψ02理论上相同,但实际中观测误差存在,ψ01、ψ02接近,应取其平均值,得到ψ0:
上述算式中,ψ0是整个平台的水平转动角度,数值范围介于0°至360°。
磁矢量测量元件X、Y、Z的方位,在球坐标中可用3点、4点、5点表示,与1光点、2光点处于同一体系下,也同样受到θx、θy转动、ψ0转动的影响,为计算磁分量数值,将3、4、5点的位置由球坐标转化为三维坐标:
上述算式中,X3'是3点受θx、θy影响后在X方向的新位置;Y3'是3点受θx、θy影响后在Y方向的新位置;Z3'是3点受θx、θy影响后在Z方向的新位置;
同样可得X4'、Y4'、Z4':
上述算式中,X4'是4点受θx、θy影响后在X方向的新位置;Y4'是4点受θx、θy影响后在Y方向的新位置;Z4'是4点受θx、θy影响后在Z方向的新位置。
同样可得X5'、Y5'、Z5':
上述算式中,X5'是5点受θx、θy影响后在X方向的新位置;Y5'是5点受θx、θy影响后在Y方向的新位置;Z5'是5点受θx、θy影响后在Z方向的新位置。
令:
上述算式中,a11、a12、a13、a21、a22、a23、a31、a32、a33为三元行列式参数。
XT0'、YT0'、ZT0'测量元件得到的数据,应提前对三个分量元件进行误差项校正,以保证磁矢量元件反映的数值接近理想状态。具体如下:
上述算式中,XT0'、YT0'、ZT0'是磁矢量元件校正后测量得到的数据,计量单位为nT;XT、YT、ZT是磁矢量元件校正前测量得到的数据,计量单位为nT;Kxt、Kyt、Kzt是各个传感器的斜率系数;bxt、byt、bzt是各个传感器零点误差校正值,计量单位为nT;Kxc、Kyc、Kzc是温度飘移系数,计量单位为℃-1;T是磁矢量元件温度,来自磁分量传感器上安装的温度传感器,计量单位为℃。
将整个磁矢量测量元件直接测得的数据XT0'、YT0'、ZT0'输入后,建立方程组行列式:
上述算式中,XT'、YT'、ZT'是真地理坐标环境下的磁矢量数据,计量单位为nT。
利用线性法进行求解,即可得到真地理坐标环境下的磁场的三个正交分量的数值XT'、YT'、ZT'。
由于采用了上述技术方案,与现有技术相比,本发明的优点是:本发明公开了一种低空光学定姿三分量磁测***,规避了昂贵惯导器件的定姿方式,体积质量较小,空中负荷轻,空气动力条件较好,***耗能较小,适用于较长时间的空中飞行,常规的天气条件下均能适用;本发明***成本较低,对相机元件的像素要求不高,使用普通的像素数量,可以完成所需精度的定姿功能;程序简化,占用计算资源少,运行速度快,满足高精度密集采样要求;采用合理的***布局与较小的飞行载体,大幅降低了对磁矢量元件的电磁干扰;本发明光学定姿原理也可推广至高空或地面磁场矢量测量模式,本***在航空磁测、铁矿勘探、多金属矿勘探、地质调查、军事等领域有较大的应用前景。
下面结合附图和具体实施方式对本发明作进一步说明。
附图说明
附图1是本发明低空光学定姿三分量磁测***中的航空机结构示意图。
附图2是本发明低空光学定姿三分量磁测***中的航空机中相机组件安装结构示意图。
附图3是本发明低空光学定姿三分量磁测***中各组件功能联系示意图。
附图4是本发明的***程序执行流程示意图。
附图5是本发明中P1—P5光点处的灰度图。
附图6是本发明中三维坐标系、球坐标系示意图。
附图7是本发明中地磁场矢量模值误差散点图。
图中:1、飞行器;2、主机;3、卫星导航模块;4、VCSEL光带环;5、阻尼杆;6、相机;7、IMU模块;8、磁分量传感器;9、刚性杆;10、支撑杆;11、简易防护条。
具体实施方式
实施例1:低空光学定姿三分量磁测***,包括地面主机和若干台航空机,所述地面主机以无线传输方式给航空机发送轨迹指引、测量参数、卫星导航基站差分数据,如附图1中航空机结构所示,所述航空机包括飞行器1、主机2、卫星导航模块3、无线传输模块、电源、VCSEL光带环4、阻尼杆5、若干台相机、IMU模块7、倾角传感器和磁分量传感器8,所述主机2包括处理器和存储器,所述卫星导航模块3包括卫星导航信号天线,所述IMU模块7包括微型惯性器件,所述VCSEL光带环4由VCSEL光点组成,所述卫星导航信号天线安装于飞行器1上端,所述主机和电源安置于飞行器1上,所述VCSEL光带环4安置于飞行器1下部,所述阻尼杆5位于VCSEL光带环4下部,所述相机6和微型惯性器件通过安装平台安置于阻尼杆5下部,所述磁分量传感器8通过刚性杆9安装于阻尼杆5下部,所述磁分量传感器8上安装有温度传感器,所述主机2与卫星导航模块3、无线传输模块、VCSEL光带环4、相机6、IMU模块7、温度传感器和磁分量传感器8电连接。磁分量传感器8为磁三分量磁通门元件,包括磁场矢量X元件、磁场矢量Y元件和磁场矢量Z元件,分别获得X、Y、Z分量数据。温度传感器用以监测磁分量传感器8的温度。
地面主机实现地面台站与观测磁场日变功能。地面主机可以采用航空机相同结构,地面主机的区别在于可以使用附加电源与较稳定的固定支架。相同结构的地面主机和航空机可以简化设备类型,现场工作时可以随机配用,且使用附加电源可以延长测量作业时间。
所述飞行器1为旋翼式飞行器。所述飞行器1机身涂装有与背景亮度接近的粗糙近红外光线全反射材料。
阻尼杆5为具阻尼作用的连接杆,其长度取0.8-1.0m,阻尼杆5具有三个功能:1.加长长度,减少主机与旋翼式飞行器对磁三分量元件的电磁干扰;2.具有阻尼作用,使旋翼式飞行器的震动传递到下部的震动减弱减小,特别是大幅减少高频分量;3.起飞前支撑。可以根据功能选取适宜材质制作,实现阻尼作用,隔离震动。阻尼是指任何振动***在振动中,由于外界作用,如流体阻力、摩擦力等,和/或***本身固有的原因引起的振动幅度逐渐下降的特性,以及此一特性的量化表征。阻尼的力学模型在物理学和工程学上应用广泛,根据现有技术即可构造具有阻尼作用的连接杆,不赘述。
相机6与磁分量传感器8之间的刚性杆9长度约0.4m。通过具有较好刚性的刚性杆9连接相机6、微型惯性器件7和磁分量传感器8,既可以使相机6与磁分量传感器8间保持适当距离,尽可能减少相机6对X、Y、Z元件的电磁干扰,还可以使三者间避免弹性震动误差,提高精度。在航机倾斜时,VCSEL光点处位置与卫星导航信号处差别较小,但相机6位置通过阻尼杆5连接因此影响较大,解决办法可以显著缩短阻尼杆5,甚至阻尼杆5缩短为零,通过适当加长刚性杆9以减少航机对磁矢量X、Y、Z元件的干扰,同时减少相机6与卫星导航信号误差。
飞行器1具有在磁分量传感器8之下设置的三条支撑杆10,用于起飞前在地面起支撑作用。在飞行器1旋翼外设简易防护条11,防止降落时倾倒损坏。如进一步地减少干扰,将阻尼杆5加长至1.0m或以上时,应在阻尼杆5最上部设置铰链,在起飞、降落时将相机、磁分量探头、刚性杆等元件收起,在主机处另外设置起降支撑杆。
如附图2所示,若干台所述相机6环绕围绕安装平台均匀排布。所述相机数量为6台或者8台。相机6采用6台时,相机6为水平视角64°近红外感光相机;相机6采用8台时,相机6为水平视角49°近红外感光相机;相机6像素可用1280×720像素,尺寸适中,既能满足测量精度需要,又可减少数据计算,提高效率。相机6采用固定光圈镜头,焦距固定在无穷远处。相机6镜头与感光元件之间安装半带宽15nm的窄带滤光片,滤光片中心波长取943-945nm。在光源频谱中心为940nm时,由于镜头与感光元件之间最边缘处的光线与主轴夹角可能接近10°,使得窄带滤光片的窄带滤光特征改变,峰值向左偏移低于940nm,且VCSEL光点光源在使用较长时间后,温度上升,会使光线光谱中心右移,滤光片中心波长取943-945nm最优。
本申请中镜头的后焦加长,使感光元件之前接收到的光线尽可能接近垂直方向地通过窄带滤光片,兼采用较小半带宽的滤光片,进一步地减少杂散光干扰。
所述相机6帧率为20fps以上,快门为全局快门,以适应动态下的拍摄,快门速度低于1/1000s,快门速度与每一帧图像的拍摄可由外部信号控制。
所述VCSEL光带环4为VCSEL球笼型光带环,VCSEL球笼型光带环由三排密集均匀排布的VCSEL光点组成,三排所述VCSEL光点中的一排位于球笼型的中间位置,另两排所述VCSEL光点对称分布在其上下位置且与中间VCSEL光点交错设置。所述VCSEL光点即垂直腔面发射激光器,VCSEL光线中心波长940nm。本***使用波长940nm的近红外激光作为标志点光源,在地表或低空环境中,940nm的自然背景光强相对较弱,从而减少杂散光干扰。以地球仪为例,赤道位置排布36只以上VCSEL光点,在北纬和南纬6-8°线上排布36只以上VCSEL光点,并与赤道位置光点错开排布,每只VCSEL光点功率3-5w,由主机通过脉冲电源控制VCSEL光点,一般接通时间2ms以下,整个VCSEL光带环的平均功率4-6w。
磁分量传感器之上安装无磁的温度传感器,监控元件温度变化。
主机2包括处理器和存储器,所述处理器接收卫星导航模块、IMU模块、倾角传感器和磁分量传感器的电信号,所述处理器通过无线传输模块与地面主机通讯,所述处理器通过脉冲电源控制VCSEL光带环中的VCSEL光点点亮,所述处理器通过相机拍摄照片并解算两帧照片的帧间差分情况,使用当时IMU模块的数据复算即时抖动量,优化帧间差分计算结果,帧间差分计算后利用灰度质心法求解光点像素坐标,处理器与磁分量传感器得到的X、Y、Z分量数据进行真地理坐标转换,得到设定地理坐标环境下的磁场矢量数据,将该磁场矢量数据与位置、高度信息一起经无线传输模块传输给地面主机并转入存储器。
参看附图3***中各组件功能联系示意图,可以采用双处理器和存储器,以加快解算速度,具体可以采用:主机包括处理器一、存储器一、处理器二和存储器二,所述处理器一接收卫星导航模块、IMU模块、倾角传感器和磁分量传感器的电信号,所述处理器一通过无线传输模块与地面主机通讯,所述处理器一通过脉冲电源控制VCSEL光带环中的VCSEL光点点亮并将相机拍摄照片送入处理器二,所述处理器二解算两帧照片的差分情况,并使用当时IMU模块的IMU数据复算即时抖动量,优化差分结果,差分后利用灰度质心法求解像素坐标,处理器二解算得到的像素坐标回传给处理器一,处理器一与磁分量传感器得到的X、Y、Z分量数据进行真地理坐标转换,得到设定地理坐标环境下的磁场矢量数据,将该磁场矢量数据与位置、高度信息一起经无线传输模块传输给地面主机并转入存储器一。
所述处理器二将每一帧照片暂存于存储器二,用于停飞后取出复算核对。
实施例2:一种地磁矢量数据测定方法,在至少2台定位航空载体上分别安装定位激光光源和卫星定位装置,在测量航空载体上安装相机组件、卫星定位装置和磁场矢量测量***,使用相机瞬时观测到的定位航空载体上定位激光光源图像数据和该定位航空载体上卫星定位装置的空间坐标、结合测量航空载体自身卫星定位装置的空间坐标解算得到相机本身的空间姿态,该测量航空载体上的磁场矢量测量元件获得磁场矢量三分量数据,结合空间姿态信息和磁场矢量三分量数据解算得到真地理坐标环境下的地磁矢量数据。
GPS定位中的RTK方式可以实时提供高精度的位置信息。
本申请的低空光学定姿三分量磁测***综合利用了以上技术发展,融合了光学、GPS通讯,机器视觉、微电子、空间姿态控制,磁矢量元件、软件等技术,达到快速高精度实现大地磁场矢量测试的功能。
实施例3:参看附图4本发明***执行流程图,应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,采用实施例1中的低空光学定姿三分量磁测***通过以下步骤测定磁场矢量数据:
S1:航空机主机将时间信号传递给相机执行拍摄,同时给予远端航空机信号,使远端航空机的VCSEL光点在相机拍摄第一帧照片时同步点亮,并同步采集一次磁分量传感器的磁三分量数据;
S2:按照设定帧率,远端航空机VCSEL光点关停,同步拍摄第二帧图像,将第二帧图像参照姿态数据进行姿态改正后,与第一帧图像进行帧间差分计算,得到光点像素坐标;
S3:***程序执行中,多台相机同步工作,主机对各相机拍摄到的光点像素坐标数据进行判别,对于相机像素边缘部位的光点位置与相邻相机拍摄得到的光点位置进行比对,舍去靠近边缘的像素点,并将该位置信息反馈给处理器,控制航空机调整飞行轨迹,使相机得到的像素光点位置远离像素框边缘,提高解算精度;
S4:初步判定像素位置后,航空机保持相对稳定飞行状态时,锁定一个较小的像素范围,减少下一步解算帧间差分图像的计算范围,减轻程序计算的任务量,提高运行速度;
S5:得到光点像素坐标后,解算球坐标环境下的θ、ψ数据,结合提供VCSEL光点的远端航空机的卫星导航信号得到的θ′、ψ′共同解算本机的姿态数据,引入前一时刻测量得到的磁三分量数据,解算真地理坐标环境下的磁分量数值,与本机位置、高程信息进行组合后加以存储,并/或直接发送到地面工作的台站。
其中,帧间差分计算方法如下:
使VCSEL光点点亮时相机拍摄第一帧照片,其像素灰度表示为P1(i,j),然后VCSEL光点关停,相机拍摄第二帧照片,其像素灰度表示为P2(i,j);
其中,P1(i,j)是VCSEL光点点亮时相机拍摄第一帧图像;P2(i,j)是VCSEL光点关停时相机拍摄第二帧图像;P表示像素灰度数值,i表示左右方向的像素坐标,由左至右数值增大,i取值范围为1-n;j表示上下方向的像素坐标,由上至下数值增大,j取值范围为1-m;
在P2(i,j)中,根据IMU模块中的微型惯性器件给出的第一帧与第二帧之间的姿态变动数据,分别求相机平台在水平、垂直方向位移与转动后像素的新位置:
上述算式中,noij是受相机水平抖动与自转影响、水平变动后图像像素坐标错动量,计算数值精度为小数点后3位以上;moij是受相机竖直抖动与自转影响、垂直变动后图像像素坐标错动量,计算数值精度为小数点后3位以上;Plev是水平抖动变化量,计算数值精度为小数点后3位以上;Pver是垂直抖动变化量,计算数值精度为小数点后3位以上;θrot是相机沿镜头轴线自转抖动的角度,数值范围一般为0.2°以下;
得出应有的校正量:
上述算式中,nij是为消除相机水平抖动与自转影响、水平变动后图像像素坐标修正量,计算数值精度为小数点后3位以上;mij是为消除相机竖直抖动与自转影响、垂直变动后图像像素坐标修正量,计算数值精度为小数后3位以上;
由于两帧之间时间间隔较小,其校正量是一个微小的数值,一般不超过相当的3-5个像素角度,所以由于像素自身转动引起每一个像素方格倾斜量的校正可忽略,只进行像素点垂直与水平位移的影响计算。
将P2(i,j)的图像像素灰度投影到新的坐标图像中,得到P3(i,j),然后与第一帧图像进行直接差分,得到差分图像P4(i,j):
上述算式中,P3(i,j)是对第二帧P2(i,j)图像进行抖动校正后的图像;P4(i,j)是差分后的光点灰度图像。
在差分图像中,相似的背景数据减差后基本为0值,光点数据就得到了明显的突出表示,但噪音数据不能完全消除。
由于噪音影响或较大灰度值的异常像素点可能存在,不能使用单点的像素灰度值搜索最大亮度点Pmax的位置,可采用以下方法消除噪音数据以获得准确像素位置:
使用9个网格数据组合形成较大的搜索框P4max(i,j),并降低最大灰度值标准为Psea,算式如下:
上述算式中,P5(i,j)是缩小搜索范围的光点灰度图像;Pmax是在P4(i,j)中单点的像素灰度中最大的灰度数值;Psea是在P4(i,j)中进行组合框搜索时使用的降低标准的灰度数值;KL是降低最大灰度标准的系数,一般取值0.4-0.7;P4max(i,j)是组合框搜索范围。
对整个图像P4(i,j)进行搜索,当检索到P4max(i,j)大于Psea时,就立刻记录所在像素点的位置,记为I1、J1,不再进行搜索。当然,在全部像素位置中,可能存在多个像素点满足搜索要求,本步骤只是按照约定搜索方向找到第一个符合条件的像素点。
设定亮点像素周边范围,设置应用半径Pr,在P4(i,j)中,进行如下判别:
上述算式中I1是第一次组合框搜索得到的亮点中心的水平像素坐标;J1是第一次组合框搜索得到的亮点中心的竖直像素坐标;Pr是设定亮点像素周边范围的应用半径;
如上式成立,则将所在像素点灰度数据恢复为0值,形成一个围绕光点中心的圆形范围。再对圆形范围内的背景噪音进行判别:
上述算式中PK是背景噪音消除限定值,一般取3-5。
如上式成立,则将圆形范围所在像素点灰度数据恢复为0值。如上式不成立,则保留原亮点数据,得到缩小范围的亮点图像P5(i,j),像素坐标值仍保持不变。在图像P5(i,j)中,使用灰度质量中心法求解亮度中心点的像素坐标:
上述算式中I2是图像亮度中心点左右方向的坐标;J2是图像亮度中心点竖直方向的坐标;
将第一次灰度质心法得到的I2、J2数据再次代入到原I1、J1位置,重新界定一个新的数据缩小范围,使用相同算法再次进行灰度质心求解,得到像素坐标I3、J3,就可提高亮点中心的定位精度,获得定位精度可达0.1倍图像像素尺度的高精度的光点中心位置,也即像素位置。除灰度质心法外,也可使用其他方法求取光点中心位置。不同阶段的光点位置灰度图见图5。
将各个相机上得到的高精度的像素坐标数值按照校正公式转化为球坐标体系下的视角数据θ、ψ,就可与卫星导航定位数据一起进行定姿计算。
坐标体系运动定位计算方法如下:
参看图6三维坐标系、球坐标系示意图所示,建立三维坐标系,并同时建立球坐标系,设球坐标模值为1构成圆球面,倾角θ以三维坐标的X、Y平面为零点,方位角ψ以三维坐标的X、Z平面为零点。1光点、2光点表示相机与卫星导航控制的两台远端航空机光点方位。在相机控制下的1光点、2光点的方位与磁矢量X、Y、Z元件方向同时固定在一个姿态之中,在1光点、2光点受到卫星导航控制产生方位变化时,磁矢量X、Y、Z元件方向随之同步变化。磁矢量X、Y、Z元件即磁分量传感器。
在从球坐标系圆心到X轴方向上,通过1光点与X轴垂直的平面与圆球面相交构成X1相交圆,X1相交圆的半径为:
上述算式中θ1是相机像素计算的1光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ1是相机像素计算的1光点的方位角弧度数据,数值范围介于0°至360°。
视线方向沿X轴从外向内观测,以X轴左侧为起点1光点与该X1相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到Y轴方向上,通过1光点与Y轴垂直的平面与圆球面相交构成Y1相交圆,Y1相交圆的半径为:
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,1光点与Y1相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到X轴方向上,通过2光点与X轴垂直的平面与圆球面相交构成X2相交圆,X2相交圆的半径为:
上述算式中θ2是相机像素计算的2光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ2是相机像素计算的2光点的方位角弧度数据,数值范围介于0°至360°。
视线方向沿X轴从外向内观测,以X轴左侧为起点2光点与该X2相交圆圆心的连线与基准水平面的夹角为:
在从球坐标系圆心到Y轴方向上,通过2光点与Y轴垂直的平面与圆球面相交构成Y2相交圆,Y2相交圆的半径为:
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,2光点与Y2相交圆圆心的连线与基准水平面的夹角为:
原始数据θ1、ψ1和θ2、ψ2与运动后的θ1'、ψ1'和θ2'、ψ2'共8个观测数据,设θx、θy为整个平台分别围绕X轴、Y轴转动的角度,是未知量,由θ1、ψ1和θ2、ψ2可得出θrx1、rx1、θry1、ry1、θrx2、rx2、θry2、ry2的解析表达式,建立二元方程组:
上述算式中θ1'是卫星导航方式计算的1光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ1'是卫星导航方式计算的1光点的方位角弧度数据,数值范围介于0°至360°;θ2'是卫星导航方式计算的2光点的倾角弧度数据,数值范围介于﹣30°至﹢30°;ψ2'是卫星导航方式计算的2光点的方位角弧度数据,数值范围介于0°至360°;rx1是过1光点与X轴垂直的平面与球坐标圆球面的交线的圆的半径,球坐标圆球面半径为1,rx1取值范围为0-1;θrx1是圆球面上1光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度,数值范围介于﹣180°至﹢180°;ry1是过1光点与Y轴垂直的平面与球坐标圆球面的交线的圆的半径,球坐标圆球面半径为1,取值范围为0-1;θry1是圆球面上1光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度,数值范围介于﹣180°至﹢180°;rx2是过2光点与X轴垂直的平面与球坐标圆球面的交线的圆的半径,球坐标圆半径为1,rx2取值范围为0-1;θrx2是圆球面上2光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度,数值范围介于﹣180°至﹢180°;ry2是过2光点与Y轴垂直的平面与球坐标圆球面的交线的圆的半径,球坐标圆半径为1,取值范围为0-1;θry2圆球面上2光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度,数值范围介于﹣180°至﹢180°;θx是X轴的假设旋转角度,数值范围介于﹣30°至﹢30°;θy是Y轴的假设旋转角度,数值范围介于﹣30°至﹢30°。
该式求解可得到θx、θy。由于三角函数内都有θx、θy,左侧分别为θrx1、θrx2、θry1、θry2等,不能视为同一变量,所以该方程组不能用线性法直接求解。应适当变换,将等式左侧强制分解,使θrx1、θrx2、θry1、θry2固定,等式右侧只取变化量。将不同参数对θx、θy的影响变化取系数,就可将θx、θy作为独立项列出:
由于上式左侧:
可以认为如下形式成立:
其中a'、b'为系数,a'受到θrx1、θx影响,b'受到θry1、θy影响。
引入如下系数:Kx、Ky为θx、θy变化时,弧长与弦长的比例关系。Krx1、Krx2为θx变化的中点处,受θx变化影响,sinθx的变化幅度系数,Kry1、Kry2为θy变化的中点处,受θy变化影响,sinθy的变化幅度系数:
上述算式中Kx是当X轴旋转θx时其弦长与弧长的对比系数,在θx取值范围﹣30°至﹢30°时,计算精度可达10﹣5度;Ky是当Y轴旋转θy时其弦长与弧长的对比系数,在θy取值范围为﹣30°至﹢30°时,计算精度可达10﹣5度;Krx1是在X轴旋转时1光点在初始与第二次观测之间的角度位置系数;Kry1是在Y轴旋转时1光点在初始与第二次观测之间的角度位置系数;Krx2是在X轴旋转时2光点在初始与第二次观测之间的角度位置系数;Kry2是在Y轴旋转时2光点在初始与第二次观测之间的角度位置系数;
令:
上式中a、b、c、d、e、f是二元方程组参数;
建立方程组:
行列式为:
解析解为:
求解θx、θy:第一步,首先将上式a-d中含有的未知θx、θy置零,求解第一次结果,将得到θx、θy值再置入方程组,进行第二次迭代,所得θx、θy值再置入方程组,进行第三次迭代,以此方法重复求解到4-20步时,直至得到的θx、θy精度小于千分之一度,满足观测精度要求。
在θ1、θ2相差较小时,ψ1、ψ2处于非0°或非180°时,可在4步求出较高精度结果;在θ1、θ2方向同时变化时,在ψ1、ψ2相差180°左右,如160°或200°时,须20步才可求得百分之一度的精度的结果,在150°-210°时,20步可求得千分之一精度的结果;在θ1、θ2反向变化时,在ψ1、ψ2相差0°左右,如20°、-20°时,须20步达到百分之一;30°或-30°时,20步可达千分之一度的结果。
在ψ1、ψ2夹角为0°或180°时无法定姿,方程式无意义,实际应用时,ψ1、ψ2夹角一般大于30°,小于150°。
求得θx、θy后,设***平台初始的绕Z轴转动量为ψ0,ψ0与ψxy1、ψxy2为独立关系,共同叠加后形成观测量ψ1'和ψ2',可利用正算式求出1光点和2光点由θx、θy转动引起的绕Z轴相应的微小转动量ψxy1和ψxy2:
上述算式中,X1'是1光点受θx、θy影响后在X方向的新位置;Y1'是1光点受θx、θy影响后在Y方向的新位置;X2'是2光点受θx、θy影响后在X方向的新位置;Y2'是2光点受θx、θy影响后在Y方向的新位置。ψxy1是1光点由θx、θy转动引起的绕Z轴相应的微小转动弧度量,数值范围为0°至5°;ψxy2是2光点由θx、θy转动引起的绕Z轴相应的微小转动弧度量,数值范围为0°至5°;ψxy1、ψxy2一般在5°以内。
由于:
得:
上述算式中,ψ01是1光点数据计算得到的水平转动角度,数值范围介于0°至360°;ψ02是2光点数据计算得到的水平转动角度,数值范围介于0°至360°。
ψ01与ψ02理论上相同,但实际中观测误差存在,ψ01、ψ02接近,应取其平均值得到ψ0:
上述算式中,ψ0是整个平台的水平转动角度,数值范围介于0°至360°。
磁矢量测量元件X、Y、Z的方位,在球坐标中可用3点、4点、5点表示,与1光点、2光点处于同一体系下,也同样受到θx、θy转动、ψ0转动的影响,为计算磁分量数值,将3点、4点、5点的位置由球坐标转化为三维坐标:
上述算式中,X3'是3点受θx、θy影响后在X方向的新位置;Y3'是3点受θx、θy影响后在Y方向的新位置;Z3'是3点受θx、θy影响后在Z方向的新位置;
同样可得X4'、Y4'、Z4':
上述算式中,X4'是4点受θx、θy影响后在X方向的新位置;Y4'是4点受θx、θy影响后在Y方向的新位置;Z4'是4点受θx、θy影响后在Z方向的新位置。
同样可得X5'、Y5'、Z5':
上述算式中X,X5'是5点受θx、θy影响后在X方向的新位置;Y5'是5点受θx、θy影响后在Y方向的新位置;Z5'是5点受θx、θy影响后在Z方向的新位置。
令:
上述算式中a11、a12、a13、a21、a22、a23、a31、a32、a33为三元行列式参数。
XT0'、YT0'、ZT0'测量元件得到的数据,应提前对三个分量元件进行误差项校正,以保证磁矢量元件反映的数值接近理想状态。具体如下:
上述算式中XT0'、YT0'、ZT0'是磁矢量元件校正后测量得到的数据,计量单位为nT;XT、YT、ZT是磁矢量元件校正前测量得到的数据,计量单位为nT;Kxt、Kyt、Kzt是各个传感器的斜率系数;bxt、byt、bzt是各个传感器零点误差校正值,计量单位为nT;Kxc、Kyc、Kzc是温度飘移系数,计量单位为℃-1;T是磁矢量元件温度,来自磁分量传感器上安装的温度传感器,计量单位为℃。
将整个磁矢量测量元件直接测得的数据XT0'、YT0'、ZT0'输入后,建立方程组行列式:
上述算式中XT'、YT'、ZT'是真地理坐标环境下的磁矢量数据,计量单位为nT。
利用线性法进行求解,即可得到真地理坐标环境下的磁场的三个正交分量的数值XT'、YT'、ZT'。
本***主机的现场应用说明:
将地面接收台站与日变测量机设置于野外平整之处,根据测量区域范围,布设测量基线,设置航空机的飞行轨迹、采集数据密度等参数。航空机高度相同,预设测线方向飞行并采集磁场矢量数据。3台航空机起飞后排布为正三角形位置,正三角形一边与测线方位垂直,一次飞行就可对3条测线进行测量,每条测线均有1台航空机进行测量。在航空机之间正三角形排布情况下,每一台测量机接收到的目标光源视角ψ1、ψ2相差60°左右,在姿态结算中可以保证得到准确的定位结果。
4台航空机时,将其位置排布设置为矩形,如果矩形一边与测线方位一致,则一次就对2条测线进行测量,每条测线均有2台航空机进行重复的测量。每一台测量机接收到的目标光源视角ψ1、ψ2相差90°左右,其姿态结算更为精确。但二者之间存在另一目标光源,其视角在ψ1、ψ2之间,与另外2台航机目标光源之间的夹角为45°,每台相机的视野中就会出现2个光源点,需要设定判别程序,识别不同位置的光源点信息。同样原理,5台与5台以上航空机时,易引起光源点的干扰现象,不推荐使用。
测量航机起飞状态时,测线测量前,航机缓慢旋转,利用不同位置的相机数据计算当前位置的磁场矢量数据,检测相机与磁矢量元件的校正情况。在进行测线测量前与测量结束后,均应利用多个航空机的与日变机之间的交互测试功能,对地面日变机进行静态核对,以使日变数据准确可靠。
卫星导航定位中的RTK方式采样频率较低,可根据微型IMU的姿态变动数据与飞行速度数据,在RTK方式采样的定位数据之间进行插值定位,提高磁场矢量测量采集密度。
航空机之间应保持200m以上的距离,以保证定位精度;在空间光线透视性较好的情况下,航空机之间距离可增大至500m-600m。对于测线间距小于200m时,可按照测线间距要求往返飞行测量。
根据天气情况,应记录太阳或较强的包含近红外光谱的人工光源情况,防止相机图像中出现眩光干扰,影响定姿精度。在后期检查数据时,结合存储器2导出的图像数据分析,局部测点的磁场矢量测量数据偏差较大时,应予以剔除。
本发明低空光学定姿三分量磁测***的精度情况如下:
相机6采用6台时,相机6为水平视角64°、水平像素1280的近红外感光相机,具有12800以上的视角分辨率,所以相机光点定位精度为0.005°,航空机之间距离200m时,GPS定位误差优于0.005°,磁三分量X、Y、Z元件校正后的方位误差0.01°,磁三分量X、Y、Z元件电气误差5.0nT,大地磁场值52000nT、磁倾角55°、磁倾角5°。模拟姿态解算、磁场矢量求解。真地理坐标下的X方向磁矢量模值与真实背景值的精度情况见图7地磁场矢量模值误差散点图,其平均误差约为7.5nT。由于***每个单次测量采样的数据都是独立的,且测量采样频率可达20hz,在航机的航迹稳定的情况下,对测量数据进行多点平均,可较大幅度地提高测量精度。
关于地磁场强度的单位,在国际单位制SD中B的单位为特,符号T。在高斯单位制CGSM中,B的单位为高斯,符号Gs。两单位制之间有1Gs=10-4T。由于特单位太大,常用更小的单位nT来表示,1nT=10-9T。
Claims (8)
1.低空光学定姿三分量磁测***,其特征在于:包括地面主机和若干台航空机,所述地面主机以无线传输方式给航空机发送轨迹指引、测量参数、卫星导航基站差分数据,所述航空机包括飞行器、主机、卫星导航模块、无线传输模块、电源、VCSEL光带环、阻尼杆、若干台相机、IMU模块、倾角传感器和磁分量传感器,所述主机包括处理器和存储器,所述卫星导航模块包括卫星导航信号天线,所述IMU模块包括微型惯性器件,所述VCSEL光带环由VCSEL光点组成,所述卫星导航信号天线安装于飞行器上端,所述主机和电源安置于飞行器上,所述VCSEL光带环安置于飞行器下部,所述阻尼杆位于VCSEL光带环下部,所述相机和微型惯性器件通过安装平台安置于阻尼杆下部,所述磁分量传感器通过刚性杆安装于阻尼杆下部,所述磁分量传感器上安装有温度传感器,所述主机与卫星导航模块、无线传输模块、VCSEL光带环、相机、IMU模块、温度传感器和磁分量传感器电连接。
2.根据权利要求1所述的低空光学定姿三分量磁测***,其特征在于:
若干台所述相机环绕围绕安装平台均匀排布,所述相机采用固定光圈镜头,所述相机镜头的焦距固定在无穷远处,所述相机镜头与感光元件之间安装半带宽15nm的窄带滤光片。
3.根据权利要求1所述的低空光学定姿三分量磁测***,其特征在于:
所述VCSEL光带环是由三排均匀排布的VCSEL光点组成球笼型的VCSEL光带环,三排所述VCSEL光点中的一排位于球笼型的中间位置,另两排所述VCSEL光点对称分布在球笼型的上下位置且与中间位置的VCSEL光点位置交错设置,所述主机通过脉冲电源控制VCSEL光点。
4.根据权利要求1所述的低空光学定姿三分量磁测***,其特征在于:
所述处理器接收卫星导航模块、IMU模块、倾角传感器和磁分量传感器的电信号,所述处理器通过无线传输模块与地面主机通讯,所述处理器通过脉冲电源控制VCSEL光带环中的光点点亮,所述处理器通过相机拍摄照片并解算两帧照片的帧间差分情况,使用当时IMU模块的数据复算即时抖动量,优化帧间差分计算结果,帧间差分计算后利用灰度质心法求解光点像素坐标,处理器与磁分量传感器得到的X、Y、Z分量数据进行真地理坐标转换,得到设定地理坐标环境下的磁场矢量数据,将该磁场矢量数据与位置、高度信息一起经无线传输模块传输给地面主机并转入存储器。
5.应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,其特征在于:采用权利要求1所述的低空光学定姿三分量磁测***通过以下步骤测定磁场矢量数据:
S1:航空机主机将时间信号传递给相机执行拍摄,同时给予远端航空机信号,使远端航空机的VCSEL光点在相机拍摄第一帧照片时同步点亮,并同步采集一次磁分量传感器的磁三分量数据;
S2:按照设定帧率,远端航空机VCSEL光点关停,航空机同步拍摄第二帧图像,将第二帧图像参照姿态数据进行姿态改正后,与第一帧图像进行帧间差分计算,得到光点像素坐标;
S3:航空机主机***程序执行中,多台相机同步工作,主机对各相机拍摄到的光点像素坐标数据进行判别,对于相机像素边缘部位的光点位置与相邻相机拍摄得到的光点位置进行比对,舍去靠近边缘的像素点,并将该位置信息反馈给处理器,控制航空机调整飞行轨迹,使相机得到光点位置的像素远离像素框边缘,提高解算精度;
S4:初步判定光点像素坐标位置后,航空机保持相对稳定飞行状态时,锁定一个像素范围,减少下一步解算帧间差分图像的计算范围;
S5:得到光点像素坐标后,解算球坐标环境下的θ、ψ数据,结合提供VCSEL光点的远端航空机得到的θ′、ψ′共同解算本机的姿态数据,引入前一时刻测量得到的磁三分量数据,解算真地理坐标环境下的磁分量数值,即得磁场矢量数据,将磁场矢量数据与本机位置、高程信息进行组合后加以存储,并/或发送到地面主机。
6.根据权利要求5所述的应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,其特征在于:所述帧间差分计算方法如下:
使VCSEL光点点亮时相机拍摄第一帧照片,其像素灰度表示为P1(i,j),然后VCSEL光点关停,相机拍摄第二帧照片,其像素灰度表示为P2(i,j);
其中,P1(i,j)是VCSEL光点点亮时相机拍摄第一帧图像;P2(i,j)是VCSEL光点关停时相机拍摄第二帧图像;P表示像素灰度数值,i表示左右方向的像素坐标,由左至右数值增大,i取值范围为1-n;j表示上下方向的像素坐标,由上至下数值增大,j取值范围为1-m;
在P2(i,j)中,根据IMU模块中的微型惯性器件给出的第一帧与第二帧之间的姿态变动数据,分别求相机平台在水平、垂直方向位移与转动后像素的新位置:
,
上述算式中,noij是受相机水平抖动与自转影响、水平变动后图像像素坐标错动量;moij是受相机竖直抖动与自转影响、垂直变动后图像像素坐标错动量;Plev是水平抖动变化量;Pver是垂直抖动变化量;θrot是相机沿镜头轴线自转抖动的角度;
得出应有的校正量:
,
上述算式中,nij是为消除相机水平抖动与自转影响、水平变动后图像像素坐标修正量;mij是为消除相机竖直抖动与自转影响、垂直变动后图像像素坐标修正量;
将P2(i,j)的图像像素灰度投影到新的坐标图像中,得到P3(i,j),然后与第一帧图像进行直接差分,得到差分图像P4(i,j):
,
上述算式中,P3(i,j)是对第二帧P2(i,j)图像进行抖动校正后的图像;P4(i,j)是差分后的光点灰度图像。
7.根据权利要求6所述的应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,其特征在于:帧间差分计算后降噪并获得光点像素坐标的方法如下:
使用9个网格数据组合形成的搜索框P4max(i,j),并降低最大灰度值标准为Psea:
,
上述算式中,P5(i,j)是缩小搜索范围的光点灰度图像;Pmax是在P4(i,j)中单点的像素灰度中最大的灰度数值;Psea是在P4(i,j)中进行组合框搜索时使用的降低标准的灰度数值;KL是降低最大灰度标准的系数;P4max(i,j)是组合框搜索范围;
对整个图像P4(i,j)进行搜索,当检索到P4max(i,j)大于Psea时,就立刻记录所在像素点的位置,记为I1、J1,设定亮点像素周边范围应用半径Pr,在P4(i,j)中,进行如下判别:
,
上述算式中,I1是第一次组合框搜索得到的亮点中心的水平像素坐标;J1是第一次组合框搜索得到的亮点中心的竖直像素坐标;Pr是设定亮点像素周边范围的应用半径;
如上式成立,则将所在像素点灰度数据恢复为0值,形成一个围绕光点中心的圆形范围;再对圆形范围内的背景噪音进行判别:
,
上述算式中,PK是背景噪音消除限定值;
如上式成立,则将圆形范围所在像素点灰度数据恢复为0值;如上式不成立,则保留原亮点数据,得到缩小范围的亮点图像P5(i,j);在图像P5(i,j)中,使用灰度质心法求解亮度中心点的像素坐标:
,
上述算式中,I2是图像亮度中心点左右方向的坐标;J2是图像亮度中心点竖直方向的坐标;
将第一次灰度质心法得到的I2、J2数据再次代入到原I1、J1位置,重新界定一个新的数据缩小范围后使用相同算法再次进行灰度质心求解,得到像素坐标I3、J3;I3是第二次计算得到的图像亮度中心点左右方向的精确坐标;J3是第二次计算得到的图像亮度中心点竖直方向的精确坐标。
8.根据权利要求5所述的应用低空光学定姿三分量磁测***测定磁场矢量数据的方法,其特征在于:所述解算真地理坐标环境下的磁分量数值的方法:
建立三维坐标系,并同时建立球坐标系,设球坐标系模值为1处构成圆球面,倾角θ以三维坐标的X、Y平面为零点,方位角ψ以三维坐标的X、Z平面为零点;1光点、2光点表示相机与卫星导航控制的两台远端航空机光点方位;在相机控制下的1光点、2光点的方位与磁分量传感器方向同时固定在一个姿态之中,在1光点、2光点受到卫星导航控制产生方位变化时,磁分量传感器方向随之同步变化;
在从球坐标系圆心到X轴方向上,通过1光点与X轴垂直的平面与圆球面相交构成X1相交圆,X1相交圆的半径为:
,
上述算式中,θ1是相机像素计算的1光点的倾角弧度数据;ψ1是相机像素计算的1光点的方位角弧度数据;
视线方向沿X轴从外向内观测,以X轴左侧为起点1光点与该X1相交圆圆心的连线与基准水平面的夹角为:
,
在从球坐标系圆心到Y轴方向上,通过1光点与Y轴垂直的平面与圆球面相交构成Y1相交圆,Y1相交圆的半径为:
,
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,1光点与Y1相交圆圆心的连线与基准水平面的夹角为:
,
在从球坐标系圆心到X轴方向上,通过2光点与X轴垂直的平面与圆球面相交构成X2相交圆,X2相交圆的半径为:
,
上述算式中,θ2是相机像素计算的2光点的倾角弧度数据;ψ2是相机像素计算的2光点的方位角弧度数据;
视线方向沿X轴从外向内观测,以X轴左侧为起点2光点与该X2相交圆圆心的连线与基准水平面的夹角为:
,
在从球坐标系圆心到Y轴方向上,通过2光点与Y轴垂直的平面与圆球面相交构成Y2相交圆,Y2相交圆的半径为:
,
视线方向沿Y轴从外向内观测,以Y轴左侧为起点,2光点与Y2相交圆圆心的连线与基准水平面的夹角为:
,
原始数据θ1、ψ1和θ2、ψ2与运动后的θ1'、ψ1'和θ2'、ψ2'共得8个观测数据,设θx、θy为整个平台分别围绕X轴、Y轴转动的角度,由θ1、ψ1和θ2、ψ2可得出θrx1、rx1、θry1、ry1、θrx2、rx2、θry2、ry2的解析表达式,建立二元方程组:
,
上述算式中,θ1'是卫星导航方式计算的1光点的倾角弧度数据;ψ1'是卫星导航方式计算的1光点的方位角弧度数据;θ2'是卫星导航方式计算的2光点的倾角弧度数据;ψ2'是卫星导航方式计算的2光点的方位角弧度数据;rx1是过1光点与X轴垂直的平面与球坐标系圆球面的交线的圆的半径,rx1取值范围为0-1;θrx1是圆球面上1光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度;ry1是过1光点与Y轴垂直的平面与球坐标系圆球面的交线的圆的半径;θry1是圆球面上1光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度;rx2是过2光点与X轴垂直的平面与球坐标系圆球面的交线的圆的半径;θrx2是圆球面上2光点的相对角度,以X轴为旋转轴心,以Y轴位置的方向为零点,从X轴向里看,右旋为正增加弧度,左旋为负增加弧度;ry2是过2光点与Y轴垂直的平面与球坐标系圆球面的交线的圆的半径;θry2圆球面上2光点的相对角度,以Y轴为轴心,以X轴位置的方向为零点,从Y轴向里看,左旋为正增加弧度,右旋为负增加弧度;θx是X轴的假设旋转角度;θy是Y轴的假设旋转角度;上式中球坐标系圆球面的半径均为1;
上式求解可得到θx、θy;变换等式,将等式左侧强制分解,使θrx1、θrx2、θry1、θry2固定,等式右侧只取变化量,将不同参数对θx、θy的影响变化取系数,将θx、θy作为独立项列出,得:
,
由上式左侧:
,
得到下式:
,
其中a'、b'为系数,a'受到θrx1、θx影响,b'受到θry1、θy影响;
引入如下系数:Kx、Ky为θx、θy变化时弧长与弦长的比例关系,Krx1、Krx2为θx变化的中点处受θx变化影响sinθx的变化幅度系数,Kry1、Kry2为θy变化的中点处受θy变化影响sinθy的变化幅度系数,得到:
,
上述算式中,Kx是当X轴旋转θx时其弦长与弧长的对比系数;Ky是当Y轴旋转θy时其弦长与弧长的对比系数;Krx1是在X轴旋转时1光点在初始与第二次观测之间的角度位置系数;Kry1是在Y轴旋转时1光点在初始与第二次观测之间的角度位置系数;Krx2是在X轴旋转时2光点在初始与第二次观测之间的角度位置系数;Kry2是在Y轴旋转时2光点在初始与第二次观测之间的角度位置系数;
令:
,
上式中a、b、c、d、e、f是二元方程组参数;
建立方程组:
,
行列式为:
,
解析解为:
,
求解θx、θy:第一步,首先将上式a--d中含有的未知θx、θy置零,求解第一次结果,将得到θx、θy值再置入方程组,进行第二次迭代,所得θx、θy值再置入方程组,进行第三次迭代,以此方法重复求解4-20步直至得到的θx、θy精度达到或者超过千分之一度;
求得θx、θy后,设***平台初始的绕Z轴转动量为ψ0,1光点由θx、θy转动引起的绕Z轴相应的微小转动量为ψxy1,ψ0与ψxy1为独立关系,共同叠加后形成观测量ψ1',2光点由θx、θy转动引起的绕Z轴相应的微小转动量为ψxy2,ψ0与ψxy2为独立关系,共同叠加后形成观测量ψ2',以下式求解微小转动量ψxy1和ψxy2:
,
上述算式中,X1'是1光点受θx、θy影响后在X方向的新位置;Y1'是1光点受θx、θy影响后在Y方向的新位置;X2'是2光点受θx、θy影响后在X方向的新位置;Y2'是2光点受θx、θy影响后在Y方向的新位置;
由于:
,
得:
,
上述算式中,ψ01是1光点数据计算得到的水平转动角度;ψ02是2光点数据计算得到的水平转动角度;
ψ01与ψ02取平均值得到ψ0:
,
磁矢量测量元件所测X、Y、Z的方位在球坐标中用3点、4点、5点表示,与1光点、2光点处于同一体系下,将3点、4点、5点的位置由球坐标转化为三维坐标:
,
上述算式中,X3'是3点受θx、θy影响后在X方向的新位置;Y3'是3点受θx、θy影响后在Y方向的新位置;Z3'是3点受θx、θy影响后在Z方向的新位置;
,
上述算式中,X4'是4点受θx、θy影响后在X方向的新位置;Y4'是4点受θx、θy影响后在Y方向的新位置;Z4'是4点受θx、θy影响后在Z方向的新位置;
,
上述算式中,X5'是5点受θx、θy影响后在X方向的新位置;Y5'是5点受θx、θy影响后在Y方向的新位置;Z5'是5点受θx、θy影响后在Z方向的新位置;
令:
,
上述算式中,a11、a12、a13、a21、a22、a23、a31、a32、a33为三元行列式参数;
XT0'、YT0'、ZT0'为三个磁矢量测量元件直接测得的数据,对三个磁矢量测量元件进行误差项校正,算式如下:
,
上述算式中,XT0'、YT0'、ZT0'是磁矢量元件校正后测量得到的数据,计量单位为nT;XT、YT、ZT是磁矢量元件校正前测量得到的数据,计量单位为nT;Kxt、K yt、Kzt是各个传感器的斜率系数;bxt、byt、bzt是各个传感器零点误差校正值,计量单位为nT;Kxc、Kyc、Kzc是温度飘移系数,计量单位为℃-1;T是磁矢量元件温度,计量单位为℃;
将数据XT0'、YT0'、ZT0'输入后,建立方程组行列式:
,
上述算式中,XT'、YT'、ZT'是真地理坐标环境下的磁矢量数据,计量单位为nT;
利用线性法进行求解,即可得到真地理坐标环境下的磁场的三个正交分量的数值XT'、YT'、ZT'。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410450963.5A CN118050814B (zh) | 2024-04-16 | 2024-04-16 | 低空光学定姿三分量磁测***及测定磁场矢量数据的方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202410450963.5A CN118050814B (zh) | 2024-04-16 | 2024-04-16 | 低空光学定姿三分量磁测***及测定磁场矢量数据的方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN118050814A true CN118050814A (zh) | 2024-05-17 |
CN118050814B CN118050814B (zh) | 2024-07-12 |
Family
ID=91052252
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202410450963.5A Active CN118050814B (zh) | 2024-04-16 | 2024-04-16 | 低空光学定姿三分量磁测***及测定磁场矢量数据的方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN118050814B (zh) |
Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20090326816A1 (en) * | 2006-05-30 | 2009-12-31 | Choon Bae Park | Attitude correction apparatus and method for inertial navigation system using camera-type solar sensor |
CN102520455A (zh) * | 2011-12-14 | 2012-06-27 | 吉林大学 | 航空地磁矢量检测装置 |
CN102621584A (zh) * | 2011-12-14 | 2012-08-01 | 吉林大学 | 航磁矢量数据处理方法 |
CN105094138A (zh) * | 2015-07-15 | 2015-11-25 | 东北农业大学 | 一种用于旋翼无人机的低空自主导航*** |
CN106989744A (zh) * | 2017-02-24 | 2017-07-28 | 中山大学 | 一种融合机载多传感器的旋翼无人机自主定位方法 |
CN107014380A (zh) * | 2017-05-26 | 2017-08-04 | 西安科技大学 | 基于飞行器的视觉导航与惯性导航的组合导航方法 |
CN109916394A (zh) * | 2019-04-04 | 2019-06-21 | 山东智翼航空科技有限公司 | 一种融合光流位置和速度信息的组合导航算法 |
CN110736457A (zh) * | 2019-11-12 | 2020-01-31 | 苏州工业职业技术学院 | 一种基于北斗、gps和sins的组合导航方法 |
US20200327696A1 (en) * | 2019-02-17 | 2020-10-15 | Purdue Research Foundation | Calibration of cameras and scanners on uav and mobile platforms |
CN112444245A (zh) * | 2020-11-17 | 2021-03-05 | 大连理工大学 | 一种基于偏振光、光流矢量、双目视觉传感器的仿昆虫视觉组合导航方法 |
CN114216454A (zh) * | 2021-10-27 | 2022-03-22 | 湖北航天飞行器研究所 | 一种gps拒止环境下基于异源图像匹配的无人机自主导航定位方法 |
CN117739972A (zh) * | 2024-02-18 | 2024-03-22 | 中国民用航空飞行学院 | 一种无全球卫星定位***的无人机进近阶段定位方法 |
-
2024
- 2024-04-16 CN CN202410450963.5A patent/CN118050814B/zh active Active
Patent Citations (12)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20090326816A1 (en) * | 2006-05-30 | 2009-12-31 | Choon Bae Park | Attitude correction apparatus and method for inertial navigation system using camera-type solar sensor |
CN102520455A (zh) * | 2011-12-14 | 2012-06-27 | 吉林大学 | 航空地磁矢量检测装置 |
CN102621584A (zh) * | 2011-12-14 | 2012-08-01 | 吉林大学 | 航磁矢量数据处理方法 |
CN105094138A (zh) * | 2015-07-15 | 2015-11-25 | 东北农业大学 | 一种用于旋翼无人机的低空自主导航*** |
CN106989744A (zh) * | 2017-02-24 | 2017-07-28 | 中山大学 | 一种融合机载多传感器的旋翼无人机自主定位方法 |
CN107014380A (zh) * | 2017-05-26 | 2017-08-04 | 西安科技大学 | 基于飞行器的视觉导航与惯性导航的组合导航方法 |
US20200327696A1 (en) * | 2019-02-17 | 2020-10-15 | Purdue Research Foundation | Calibration of cameras and scanners on uav and mobile platforms |
CN109916394A (zh) * | 2019-04-04 | 2019-06-21 | 山东智翼航空科技有限公司 | 一种融合光流位置和速度信息的组合导航算法 |
CN110736457A (zh) * | 2019-11-12 | 2020-01-31 | 苏州工业职业技术学院 | 一种基于北斗、gps和sins的组合导航方法 |
CN112444245A (zh) * | 2020-11-17 | 2021-03-05 | 大连理工大学 | 一种基于偏振光、光流矢量、双目视觉传感器的仿昆虫视觉组合导航方法 |
CN114216454A (zh) * | 2021-10-27 | 2022-03-22 | 湖北航天飞行器研究所 | 一种gps拒止环境下基于异源图像匹配的无人机自主导航定位方法 |
CN117739972A (zh) * | 2024-02-18 | 2024-03-22 | 中国民用航空飞行学院 | 一种无全球卫星定位***的无人机进近阶段定位方法 |
Non-Patent Citations (2)
Title |
---|
卢鸿谦;尹航;黄显林;: "偏振光/地磁/GPS/SINS组合导航方法", 宇航学报, no. 04, 31 July 2007 (2007-07-31) * |
林君;刁庶;张洋;张扬;赵静;: "地球物理矢量场磁测技术的研究进展", 科学通报, no. 23, 20 August 2017 (2017-08-20) * |
Also Published As
Publication number | Publication date |
---|---|
CN118050814B (zh) | 2024-07-12 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
EP1019862B1 (en) | Method and apparatus for generating navigation data | |
US9478034B1 (en) | Geoposition determination by starlight refraction measurement | |
CN110108984B (zh) | 电力巡线激光雷达***多传感器的空间关系同步方法 | |
GREJNER‐BRZEZINSKA | Direct exterior orientation of airborne imagery with GPS/INS system: Performance analysis | |
US10969493B2 (en) | Data processing device, data processing method, and data processing program | |
Delaune et al. | Range-visual-inertial odometry: Scale observability without excitation | |
Miller et al. | 3-D site mapping with the CMU autonomous helicopter | |
CN111426320B (zh) | 一种基于图像匹配/惯导/里程计的车辆自主导航方法 | |
CN105928515B (zh) | 一种无人机导航*** | |
CN104374388A (zh) | 一种基于偏振光传感器的航姿测定方法 | |
CN110296688A (zh) | 一种基于被动地理定位技术的侦测一体倾斜航测吊舱 | |
Hill et al. | Ground-to-air flow visualization using Solar Calcium-K line Background-Oriented Schlieren | |
Ye et al. | A review of small UAV navigation system based on multi-source sensor fusion | |
Li et al. | Unmanned aerial vehicle position estimation augmentation using optical flow sensor | |
CN105928519B (zh) | 基于ins惯性导航与gps导航以及磁力计的导航算法 | |
Kong et al. | A ground-based multi-sensor system for autonomous landing of a fixed wing UAV | |
CN112415535B (zh) | 导航***和导航方法 | |
CN118050814B (zh) | 低空光学定姿三分量磁测***及测定磁场矢量数据的方法 | |
Brommer et al. | INSANE: Cross-domain UAV data sets with increased number of sensors for developing advanced and novel estimators | |
Aaron et al. | Performance Analysis of Terrain Relative Navigation Using Blue Origin New Shepard Suborbital Flight Telemetry | |
Karimi et al. | Deep fusion of a skewed redundant magnetic and inertial sensor for heading state estimation in a saturated indoor environment | |
CN114994622A (zh) | 一种基于pos的机载雷达***静态校靶方法 | |
Yuan et al. | A robust multi-state constraint optimization-based orientation estimation system for Satcom-on-the-move | |
CN112485811B (zh) | 一种直升机外吊挂物运动测量方法 | |
Spielvogel et al. | Online 3-Axis Magnetometer Hard-Iron and Soft-Iron Bias and Angular Velocity Sensor Bias Estimation Using Angular Velocity Sensors for Improved Dynamic Heading Accuracy |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant |