测距三角高程的误差分析

测距三角高程的误差分析

1928年，自温宁马尼兹提出了计算垂直面波动的公式以来，许多人研究了使用重力数据计算垂直面波动的问题。在20世纪50年代末，人们创造了一个不断变化的环境，为重力场和垂直面的确定创造了新的局面。土壤力学测量系统的实施有效地促进了垂直面偏离的确定。然而，由于缺乏详细的重力数据，山地运输变得困难，上述工作的发展也变得困难。基于以上原因,本文提出了用地形代替重力数据的方案,并对有关计算方法和精度等问题作了讨论。1垂线偏差的非线性变化在测距三角高程中垂线偏差的非线性改正(Δu)公式:Δu=(uA+u02-um)SABρΔu=(uA+u02−um)SABρ(1)式中uA、uB为A、B两点的垂线偏差,SAB为AB的距离,ρ=206265″,um为AB线路上垂线偏差的积分平均值,即:um=1SAB∫udsum=1SAB∫uds(2)如果在A、B二点上作过天文大地测量或在其周围作过详细的重力测量,则uA、uB可求,又在沿两点间的路线上作过上述测量,则um也可求,然而以上工作特别是在山区除了用于试验之外是不可能做到的。为弥补这一不足,本文将在后面讨论这一问题。现在仅就Δu值作一基本估计,实际上,我们在中已估计该值为1″,它对高差的影响达10mm之多,从山区的实际情况看,这一误差必须予以考虑(不像有些作者认为的那样是可以忽略的),例如在阿尔卑斯山,有的作者用实际数据估算,在相距0.5、1.0km的两点,估计其误差分别达到1.8″和2.8″。这从下例中也可以看出,在某山区的相邻的四个点上,A、B相距仅1.6km,而垂线偏差之差达6″,B、C二点相距不到1.1km,其差值10″,若以C、D两点而言,在南北方向上的垂线偏差分别为-7″、-11.8″,C点大致在CD线上,依直线内插,该值为-8.9″,而实测值为+8.4″,两者竟相差达17.3″。可见垂线偏差的非线性变化却是很大的。图2为第二个实例,在一剖面上的A、B、C三点均有天文大地坐标,由此可计算出在AB方向上的垂线偏差,其值如下,用直线内插得到C点的值为4.6″,而实测值为6.6″,其间差值为2″,此例也说明,在山区垂线偏差并不呈线性变化。为估算此项的影响,今设垂线偏差随距离呈二次变化,这时,在距A点的i处有:ui=uA+aSAi+bS2Ai2Ai(3)对于B点:uB=uA+aSAB+bSLABLAB(4)式中a、b为常数,将式(3)代入式(2),有:um=uA+a2SAB+b3S2ABum=uA+a2SAB+b3S2AB(5)这时(uA+uB2-um)=16bS2AB(uA+uB2−um)=16bS2AB(6)将上图中的数值代入式(3)、(4)、(5)、(6)后可得:uA+uB2-um=1.445″uA+uB2−um=1.445′′由此引起的测距高程误差达9.7mm,若顾及更高阶的影响,该值有可能还要更大,这正说明了在山区的测距三角高程中必须考虑垂线偏差的非线性影响。2区内的争议系数由于垂线偏差与计算点到扰动物质处的距离的平方成反比,故在计算点附近影响很大,这从温宁·迈乃兹函数表达式中也可以看出,尽管如此,这里我们仍然顾及远区域的影响,这时:ξ=ξc+ξN+ξF(7)式中ξF为远区结果,ξN为近区结果,ξc为中央区结果,由于现今用球函数表示的扰动位模型已经完整到360阶,并正向720阶过度,其分辨率已达到0.5°或更小,然而由模型计算的垂线偏差与实际的比较,其平均误差大达6.51″,为此我们作如下设计,在距计算点的0.5°以外用模型计算,在0.5°以内用实际资料(地形、重力)来计算,这时可用平面公式计算重力垂线偏差,直角坐标的公式:ξ=12πγ∑Δgi∫∫xdxdy(x2+y2)32=12πγ∑Δgifiξ=12πγ∑Δgi∫∫xdxdy(x2+y2)32=12πγ∑Δgifi(8)η的表达式与上式基本相同,只要将x,y互换即得,fi为内的值。如果按棋盘式布点,计算点在其中心,点距为2λ,各方块产生误差的总和为δξ=12πγ(∑iδ2gif2i)12(9)若重力点的分布及单元面积的划分是扇形环带状,圆周的半径等于(2n±1)λ,其中n是区号,每个区内的点数等于8n。该区带来的误差可以表示成公式:δξ=12πγln2n+12n-1(∑cos2αm(Δα)2)12δg(10)这时Δα=2π8n∑cos2αm=4n(11)则:δξ=14γn-1/2δgln2n+12n-1(12)如果重力的代表误差以±5×10ms-2估计,若取5个带则由此引起总的误差为±0.32″,而第一带的误差已达0.29″～0.295″,由此可以看出,无论按哪一方案布点,1区(1λ～3λ)所引起的误差为整个误差之大部,可见近区的影响是非常大的,但应注意的是这只有在均匀布点时才能按此法计算。当然,如此进行重力点的布测,人们很难作到,而这样选取高程点,(如有大比例地形图)则是比较容易的,以上是关于近区即ξN的计算及其误差的讨论。对于中央区(ξc)的计算及其误差可采用以下几种方式:(1)在文献中,设重力异常按下述方式变化:Δg=Δg0+Δgοxx+οΔgοgy(13)式中Δg0为计算点的重力异常,οΔgοx,οΔgοy分别为重力异常在x,y方向梯度,而垂线偏差及其误差:ξc=-0.105οΔgοxr0(14)δξC=±0.105mgxr0(15)式中mgx为重力异常水平梯度的误差,r0为距计算点的半径,设mgx=±4mgal/km(1gal=1cm/s2),对于r0=1km,则δξc=±0.42″。(2)由文献有:ξc=ρ4πγ[4ln√2+1]14(-Δg1+Δg2+Δg3-Δgc)=ρ4πγ[4ln√2+1]Δˉg(16)式中Δg1,Δg2,Δg3,Δgc分别为以原点为中心四个象限(方块)内的重力异常,Δˉg为以上代数和的均值,则上式的均方根差:δξc=ρ4πγ×3.77δΔˉg(17)设δΔˉg=±6×10-5ms-2,则δξc=±0.35″关于重力异常的误差δg一般用代表误差来代替,它与所代表地区的面积有关,人们常用亨特(Hunter)公式来表示,亦即:δg=C(√x+√y)(18)式中C为代表误差系数,x,y为该面积上的边长,若x=y,则:δg=2C√x(19)C值的大小既与所在地区的重力异常的复杂程度有关,又与所采用的何种重力异常有关,重力异常的复杂程度又和地形的粗糙度有关。对此我国已作了不少试验,于空间异常而言,高差变化较小的平原(高原),C值在0.6以下,高差变化小于300m的丘陵,C值在0.6～1.0之间,高差变化在300～700m的小山区,C值在1.0～1.7之间,高差变化在700～1200m的大山区,C值在1.7～3.0之间,高差变化在1200m以上的特大山区,C值在3.0～3.9之间。对于布格异常而言,C值比空间异常要小得多,例如在鄂西北的大山区,C值分别为2.41、0.68,而对法那异常而言,C值仅为0.39。由以上可以看出,在推估空间异常时为什么要采用布格异常和法那异常作间接内插,亦即所谓移去-恢复的方法。由于测距三角高程大多是在山区进行,为讨论方便,今以大山区布格异常的C值代入上式,得到以x为边长的任一方块内的代表误差:δg=1.36√x(20)对于x为1km和16km的δg分别为1.4和5.4×10-5ms-2,如果在计算点附近以0.5km的间距布点(实测或内插),半径2km～55km以1km的间距布点,其外由球函数计算,则垂线偏差δξ(δη)的总误差约为0.4″。3垂线偏差的计算以上讨论了用重力数据计算垂线偏差及其误差的原理和方法,但由于山区重力点的测量相当困难,以致分布不均匀且比较稀疏,经过推估或间接内插,目前各地基本上可以组成10′×10′乃至5′×5′为单元的重力异常,估计由此求得垂线偏差的精度也只能达到4″。尽管如此,由于近区重力资料的缺乏,还无法像上节中所述的那样来计算垂线偏差。为此我们设计:10′到30′的近区仍用重力资料来计算,而在10′以内则用地形数据来计算,问题是在这样范围内地形数据能否取代重力数据?据许多人研究,研究的区域越小,其代表性愈强,如在1°×1°范围内,空间异常与高程的相关性可以达到0.9左右,某大山区相关系数为0.96,而鄂西北大山区也在0.9以上,在国外也有不少人做了类似研究,得到许多相似结果。我们认为在半径为10′的范围几乎完全可以用地形数据来代替。事实上,空间重力异常的高频部分,绝大部分为地形起伏所制约,这种地形起伏正是扰动物质的场源。根据这一特点,我国在珠穆朗玛峰垂线偏差的计算中已成功地得到应用。由地形质量计算垂线偏差的途径有二:其一是用地形数据直接计算,现列出柱坐标和直角坐标的表示式,对于计算点(原点)任一单元质量引起的效应有:ξ=fδg∫∫∫xdxdydz(x2+y2+z2)3/2=fδg∫∫∫rcosAdAdrdz(r2+z2)3/2(21)式中f为引力常数,δ为岩石密度,g为重力加速度。对于柱坐标,某扇形的ξnk有:ξΤnk=-0.00773δδ0(sinAn-sinAn-1)hlnrk+1+√r2k+1+h2rk+√r2k+h2(22)而总的结果:ξr=n∑n=1k∑k=1ξnk(23)式中h=H-H0,δ为地壳密度,δ0为地球平均密度,r为向径,A为方位角,n为扇形号,k为环带号。如取rk+1/rk=1.426,又设mh=±15m,代入相应的误差公式,对k环有:mξTk≤±0.08″(24)若取k=9,并顾及中央区的误差(设之为0.10″),则总的误差:mξT=∑ξTK=0.25″(25)假设由地形代替重力使误差增大0.5倍,同时顾及近区10′～30′以及远区的误差,则垂线偏差(ξ)最终误差不会超过0.55″。至于由直角坐标表示的公式可以用以下两种形式表示,一种是:ξΤ=fyδ{(ysinh-1z(x2+y2)12+zsinh-1y(x2+y2)12-xtg-1yzx(x2+y2+z2)12)}|x2x1|y2y1|z2z1(26)式中x2,x1,y1,y2,z1,z2,为对任一平行六面体积分的上、下限,其余的计算类同,最后求和可得总的结果。另一种形式:ξ=fδ(z2-z1)g|lnx+ry2x+ry1|x2x1(27)式中r=(x2+y2+z2m)12‚zm=z1+z22‚ry2=(x2+y2+z2m)12‚ry1=(x2+y21+z2m)12。其二是用间接法进行计算,所谓间接法是先求出地形质量的重力效应,然后代入温宁·迈乃兹公式计算,这和前述的用重力资料计算垂线偏差的方法相同,此处不再重复。在中央区域可分为相等的八个三角形,对于平面标为(0,0),(0,1),(1,1)第一个三角形三点上的高差h00,h01和h11,且h00=0,设三角域内高差呈线性变化,即h=a0+a1x+a1y,则可以求出:a1=h01-h00,a2=h11-h01,这时该柱体内的重力效应为Δc1=fρa(ln(1+√2)-1√1+a22ln1+a22+a1a2+√1+a22√2+(a1+a2)2a1a2+√1+a21√1+a22(28)式中f为引力数,a为等腰直角三角形直边边长。其它三角域上的重力效应可类似求出。该文中还列出了非中央区域的计算公式。此外还有不少文献,列出了快速富氏变换的计算公式,此处不再赘述。4重力改良项由于在测距三角高程计算出的正常高和正高,故还需作一项改正,以往是用到重力资料,如今亦可用高程资料来代替。对于正常高改正一般用下式表示:ΔΗAB=(ΗB-ΗA)-∑idh(29)式中HA、HB为A、B点的上的正常高,∑idh为水准测量或三角高程测量所得的高差。文献根据重力异常与高程呈线性及∫gdh+∫hdg=gBHB-gAHA的关系,有:ΔHAB=ε+Δ(30)式中ε为正常水准面不平行改正,Δ为重力异常改正。ε=-0.00154ΔφHmsin2φ(31)Δ=0.0000571(H2B-H2A)=0.00102(HBΔgBB-HAΔgBA)-0.00102Hm(ΔgBB-ΔgBA)(32)式中Ηm=12(ΗB+ΗA)‚ΔgBA、ΔgBB为A、B点上的布格异常,Δφ为纬度差以角分表示,由于上式中后二项很小,可省略,这时ΔHAB=-0.00154ΔφHmsin2φm+0.0000286HmΔH(33)在上式中ΔH为A、B两点的高差,φm为平均纬度,可用下式估计δ(重力改正项)的误差:m=√2Ηm,正高改正∶ΟCAB=(ΗΔgB-hA)-∑Δhi=1g0[∑(gi-g0)Δhi+hA(ˉgA-g0)-hB(ˉgB-g0)](34)式中g0为重力近似值,gi沿水准或测距三角高程点上的重力值,ˉgA、ˉgB为在A、B两点处沿垂线方向到大地水准面上的重力平均值,hA、hB为A、B两点的正高,计算可以近似值代替。仍根据重力与高程基本是线性相关以及∫gdh+∫hdg=gBhB-gAhA的恒等关系,可以推得:ΟC=0.114×10-3ˉhΔΗ-0.83×10-3sin2φˉhS(35)式中OC以mm为单位,S以km为单位。对于山区,若ˉh=1500m,Δh=500m,S=25km,Ψ=45°时的正高改正达85mm。关于正高改正的误差为σoc,其平方为σ2oc=1g20[0.025ˉh2σ2h+nˉh2σ