集成gnss和ccd天顶筒的高精度垂线偏差快速测量方法
【技术领域】
[0001] 本发明涉及一种垂线偏差的天文大地测量方法，尤其涉及一种集成GNSS和C⑶天 顶筒的高精度垂线偏差快速测量方法。
【背景技术】
[0002] 地面上一点的重力向量g和相应椭球面上的法线向量η之间的夹角定义为该点的 垂线偏差。垂线偏差表征重力的方向，是大地测量归算、地球重力场和大地水准面模型精 化、研宄地球形状的必要地球物理量，在大地测量和空间技术等方面有着重要用途。
[0003] 垂线偏差是地球内部质量分布的表现，可用于反演大气反常折射研宄，监测地球 内部质量迀移和能量积累，监测地震、火山爆发等自然灾害，也可用于探测地下质量异常， 进行资源勘探。在军事领域，垂线偏差是导弹、火箭发射等必需的基础数据。
[0004] 现有技术中，常用的垂线偏差观测方法，主要包括传统天文大地测量方法、重力测 量方法、天文重力测量方法以及GPS水准测量方法。但是，上述这些方法普遍存在基础数据 采集过程复杂，或者根本就无法有效获取、后期数据处理过程繁琐复杂、费时费力等缺点或 不足：
[0005] 传统天文大地测量方法，通常采用Τ4、J05等天文经炜仪，配以一套守时、授时系 统。其天文炜度的观测常用Telgte方法，天文经度的观测常用双星等高法并同时收录无线 电时号。这种方法获得的垂线偏差精度可达到0.3"。
[0006] 但是，对于这种传统天文大地测量方法，通常，其为了使观测精度足够高，所使用 的天文经炜仪的体积都很大、十分笨重；且仪器的操作对观测者经验要求过高、观测时间很 长、效率极低。而且，其常用的守时、授时系统，不仅操作复杂，精度较低，且容易出错。
[0007] 重力测量方法，主要借助于大地水准面和地球椭球面上的重力异常，利用Stokes 或VeningMeinesz公式计算垂线偏差。
[0008] 然而，这种方法的前提是假定大地水准面之外没有扰动物质，且其前提是，需要对 全球重力异常数据已知，这十分困难、且根本就不现实。因此，该方法还没有得到独立的应 用，并且该方法的垂线偏差的测量精度也偏低。
[0009] 天文重力测量方法，其实质是结合天文大地测量方法和重力测量方法来确定垂线 偏差。该方法首先需要一定数量的天文大地点的垂线偏差值，然后利用重力测量方法加密 重力异常数据，采用最小二乘配置或多面函数方法确定其他点的垂线偏差值。
[0010] 这种方法的主要缺点在于，求出的垂线偏差精度偏低。
[0011] GPS水准测量方法，需要同时进行GPS测量和水准测量，获取基线长、大地方位角 和高程异常之差，就可以求出垂线偏差。
[0012] 然而，这种方法的使用方面有条件限制，比如，对于某些地形起伏较大的区域进行 精密水准测量的难度太大，基线也不能过长，且由于基线长和大地方位角的影响，导致其求 出的垂线偏差精度偏低。
【发明内容】
[0013] 本发明的目的是，提供一种具有普适特性，且高精度垂线偏差的数据采集过程省 时省力、数据的处理方法相对简单的一种集成GNSS和CCD天顶筒的高精度垂线偏差快速测 量方法。
[0014] 本发明为实现上述目的所采用的技术方案是，一种集成GNSS和CCD天顶筒的高精 度垂线偏差快速测量方法，其特征在于，包括以下步骤：
[0015] 第一步，在测站点上架设CCD天顶筒，使其光轴指向天顶方向，对天顶区域连续在 0°和180°方向进行拍摄，获取CCD恒星图像，并由高精度电子水平仪获取天顶筒的倾斜 数据，由GNSS获取观测时间信息和测站点的大地坐标；
[0016] 第二步，对CXD天顶筒拍摄的FITS格式的C⑶恒星图像进行读取，再进行背景噪 声消除；
[0017] 第三步，对噪声消除后的CCD恒星图像进行恒星星象区域的自动搜索；
[0018] 第四步，运用改进的修正二维矩方法，计算上述自动搜索出的恒星星象区域的能 量中心的影像坐标；
[0019] 第五步，对适用星表进行处理，通过GNSS获取的时间信息以及测站点的大地坐 标，结合EGM2008计算测站点粗略天文坐标。根据结果对星表中与CCD恒星图像相对应的 区域进行截取，并计算曝光时刻天顶区域的恒星赤道坐标与视星等信息，并将恒星赤道坐 标投影到切平面上转换为平面坐标；
[0020] 第六步，对C⑶影像中的星象与星表投影到平面上的恒星进行识别、匹配；
[0021] 第七步，依据匹配成功的星对，利用最小二乘法计算摄影测量投影变换公式中的 坐标转换系数，并迭代计算测站点的初始赤道坐标，再转换为测站点的初始天文坐标；
[0022] 第八步，将高精度电子水平仪获取的天顶筒倾斜改正数据和极移改正数据转换至 子午方向和卯酉方向上，确定出测站点的精确天文坐标；
[0023] 第九步，结合GNSS测量数据解算出的测站点大地坐标，利用Helmert公式，即可计 算出测站点的垂线偏差。
[0024] 上述技术方案直接带来的技术效果是，垂线偏差快速测量方法基于GNSS和CCD天 顶筒，在满足测量精度和分辨率需要的基础上，利用改进的数据采集及处理方法，降低了垂 线偏差获取过程所消耗的时间，提高了测量效率，能够满足高的测量精度和无地域限制等 要求。
[0025] 优选为，上述背景噪声消除所采用的方法为3X3中值滤波法，具体步骤如下：
[0026] 步骤一，确定一个以某个像素为中心点的3X3邻域；
[0027] 步骤二，将该邻域中各个像素的灰度值（X1, X2, X3, ...，Xn)按从小到大进行排序：
[0028] X1^ X2^x3... xn；
[0029] 步骤三，按公式计算出各个像素的灰度值的中值Y :
[0030]
[0031] 步骤四，以中值Y代替原来中心点的灰度值，进行中值滤波；
[0032] 步骤五，确定一个以另一个像素为中心点的3X3邻域，重复上述步骤二至步骤 四，直至全部像素滤波完成。
[0033] 该优选技术方案直接带来的技术效果是，去噪前后恒星识别率提高了约60%，增 加了观测量，大大削弱了随机噪声对于CCD恒星图像处理的影响。垂线偏差子午分量和卯 酉分量的标准差相较于去噪前分别平均减小了 〇. 016"和0. 021"，内符精度有所提高。
[0034] 3X3中值滤波是一种基于排序统计理论的能去除噪声的非线性处理方法，基本原 理是把数字图像中的灰度值用该点的3X3邻域中各点灰度值的中值代替。
[0035] 这里的邻域通常被称为窗口，当窗口在图像中上下左右进行移动后，利用中值滤 波可以很好地对图像进行去噪处理。
[0036] 进一步优选，上述恒星星象区域的自动搜索是采用区域生长算法进行的。
[0037] 该优选技术方案直接带来的技术效果是，采用区域生长算法进行C⑶恒星图像中 的星象自动搜索，其方法简单，所需花费的时间较短。
[0038] 原因在于，区域生长算法的基本原理为设定一个最小灰度阈值，如果CXD图像中 某区域由灰度值相近的若干相邻像素组成，且该区域内各像素灰度值均大于最小灰度阈 值，则判定该区域为一颗恒星星象区域。
[0039] 进一步优选，上述改进的修正二维矩方法，是在以往的修正二维矩方法上，运用迭 代计算方法对阈值进行改进，其步骤如下：
[0040] 第1步，提取C⑶恒星图像中最大灰度值gmax和最小灰度值g min，令
[0041]
[0042] 第2步，根据阈值Tk把C⑶恒星图像区分为星象区域和背景区域，其中，星象区域 位置（Xi，yP的灰度值gk+1 (Xi，以> 〇,而背景区域的灰度值为〇,如公式所示：
[0043]
[0044]贝U，按公式（4)求得星象区域的平均灰度值：
[0045]
[0046] 4
，将式⑷代入，可得：
[0047]
[0048] 第3步，设置一个正数μ，如果I Tk-Tk+1 I < μ，则终止迭代计算过程，按公式（6) 求得最优阈值：
[0049] T = Tk+B (6)
[0050] 如果I Tk-Tk+11彡μ，则令Tk= T k+1，继续进行第2步的迭代计算。
[0052]
[0051] 第4步，按公式（7)和式（8)，求出(XD恒星图像中星象区域中心坐标：
[0053]
[0054] 该优选技术方案直接带来的技术效果是，在以往二维修正矩算法的基础上，用迭 代法改进阈值T，使得星象中心定位精度得到较大改善，定位精度标准差控制在0. 1-0. 15 像素，对垂线偏差测量精度的影响控制在〇. Γ内。进一步优选，上述适用星表的处理是通 过调用AURIGA中的NOVAS-F程序，并结合GNSS获取的时间信息，对观测历元观测到的恒星 位置与适用星表历元恒星平位置存在差异的各影响要素进行改正，以获取星表中恒星的视 位置；
[0055] 然后，对星表中与CCD恒星星像相对应的区域进行截取，截取过程分为两步：
[0056] 第一步，在集成CCD和GPS数字天顶筒进行观测之前，首先对星表中本观测时段可 能通过测站点天顶区域的恒星信息进行截取，星表截取范围以将观测时可能出现在天顶区 域视场内的恒星全部覆盖为原则，其测站点赤道坐标（α2, δ2)由如下公式（9)求得：
[0057] δ 2= Φ 2
[0058] Q2=A 2+GAST (9)
[0059] 上式（9)中，（Φ2, A2)为测站点的初始天文坐标，该初始天文坐标是采用GNSS大 地测量结果和EGM2008重力场模型计算出测站点的垂线偏差，再推得测站点的天文坐标， 取近似值得到的；
[0060] 所述星表截取范围大于数字天顶筒视场；
[0061] 在截取范围确定后，运用恒星视位置计算模型对观测时段内可能出现在该范围内 的恒星信息进行截取，并将这些信息进行存储；
[0062] 第二步，在进行恒星匹配工作前，需对以上星表区域进行进一步截取，使该区域能 够完全与CCD恒星图像区域对应，其截取条件如下：
[0063] amin< a a max
[0064] δω?η< δ δ max (10)
[0065] 上式（10)中：
[0066] ( α ρ δ J为天顶区域视场内恒星的赤道坐标；
[0067] amin、αΜΧ、δ_、δ_取决于测站点的赤道坐标（α 2, δ2)和数字天顶筒的视场 Fi，其计算公式如下：
[0068] Qmin= a 2-Fi/2
[0069] Qmax= a 2+Fi/2 (11)
[0070] Smin= δ 2-Fi/2
[0071] 5 max= 5 2+Fi/2
[0072] 所述将恒星赤道坐标投影到平面上转换为平面坐标是按投影变换公式（1