《区域性地表形变分析中稳健坐标基准的确定.doc》由会员分享，可在线阅读，更多相关《区域性地表形变分析中稳健坐标基准的确定.doc（12页珍藏版）》请在三一文库上搜索。
1、精品论文大全区域性地表形变分析中稳健坐标基准的确定施一民1,2 罗彦 11.上海同济大学 测量系, 上海市(200092)2.现代工程测量国家测绘局重点实验室, 上海市 (200092)E-mail:摘 要: 为更客观合理地反映实际位移，提出了按位移分量的一次范数为最小的抗差估计方 法来确定相对稳定位移，与通常采用的拟稳平差不同，它毋须基于对点位形变情况的先验假 设。为有效地鉴定分水平位置和垂直方向的相对稳定点，采用了新大地坐标下的点位位移表 述。利用上海地区 GPS 监测网的两期实测数据按该方法进行了计算，初步证明了该方法应用 于形变分析的可行性和有效性.关键词地壳形变, 水平位移,稳健坐标
2、基准, 新型大地坐标系1 引言对GPS形变监测网点在垂直及水平方向上的点位形变表述, 目前大多是基于三维大地坐 标。通常是以若干IGS站在某参考历元的全球框架(如ITRF-2000)坐标为基准(强制固定或强 约束)来求解观测站在全球框架中的坐标。而因各期网观测历元并不相同，还须按站速度归 算至同一参考历元。为了能反映局部区域内所关注的地表形变,在各期的求解中通常要固定 一点或若干点的坐标。但固定点的选择一般是按照先验认识，常与实际情况不一致。若按拟 稳平差处理1,也同样需对拟稳点作先验假定并在平差后进行统计检验、迭代更换，难免也 有主观随意性。要使由两期网坐标之差异能更好地反映研究区域内所发生
3、的形变,能否不再仅依赖于先 验认识来寻求能统一各期网平差的坐标基准，使由两期网坐标之差异能客观反映其间相对位 置的变动。本文将利用抗差估计来确定局部区域中基于ITRF的稳定可靠的坐标基准，以客观反映本区 域所关心的显著的地表形变。在此采用新型大地坐标系将平面点位和水平位移统一表示在 椭球面上来进行分析。2 按三维点位变动量一次范数最小求解的原理和方法形变监测网中，为了更客观反映监测网点的位移量，应是基于若干相对稳定的基准点给 出基准数据。若是基准点确实不存在位移，则可将这些基准数据作为网平差中的固定值。然 而在变形体中实际上往往并不存在绝对意义上的“不动点”，我们所要寻求的只是相对意义 基金项
4、目：国家自然科学基金资助项目(40471114)作者简介：施一民 (1942-), 男 ,浙江宁波人, 教授，博士生导师E-mail:5上的“稳定点”，但这是在数据处理中难以可靠地断定的，往往只能是凭经验或是对以往形变的粗略了解。通常就基于Moore-penrose广义逆(伪逆)来求解，实质上它是由观测值在最小二乘意义下 的最佳相对网形通过平移、旋转乃至伸缩最优拟合于基准数据。在此所说的基准数据，指的 是在某一确定的坐标系中表述点位的一组数据（可包括点位坐标、方位角、边长），据此即 可经过平差定位实现该坐标系。这是从实现控制网定位这一角度来定义基准数据，而不论基 准数据经平差后变动与否2，3。
5、按伪逆求解是将首期所有网点坐标全部作为基准数据，所得 的坐标待定参数，其实就是两期监测网点之间拟合后的余差，亦即首期所有网点坐标与后期 网拟合定位后的最佳相对网形之间的坐标差异。按此平差所得的坐标改正数未必就能客观地 反映实际发生的形变量。当个别点上存在较大位移时，并不能准确地予以反映，往往被平滑 到其余网点上。当然亦可采用拟稳平差：仅将一部分网点（常称为拟稳点）的首期坐标作为 拟合中的基准数据，另一部分网点（常称为非拟稳点）的首期坐标则不作为基准数据, 此时 的二次范数最小的条件对应于仅在拟稳点上实现最优拟合，使得在拟稳点上坐标余差的平方 和为最小。但这又需对拟稳点作先验假定并在平差后进行统
6、计检验、迭代更换，亦有一定的 主观随意性，从而可能影响形变分析的正确性。为此 80 年代起开始采用了抗差估计（稳健估计）来限制观测值中粗差的影响和确定稳定 基准4,5。由于形变分析的目的就在于发现异常位移，因此在应用抗差估计法时，就不能如 对观测值粗差那样用抗差估计中的选权迭代法来予以剔除。本文是按位移分量的一次范数为最小的估计方法来确定相对稳定点6 ，这是由于按一次 范数最小具有较好的抗差性，不会掩盖形变点上的异常位移。但由于一次范数不是可微函数， 不能得出目标函数的一阶导数方程（法方程），可采用直接搜索目标函数极小值的单纯形法 求解。72.1 三维直角坐标形式下位移向量的间接观测方程式对于
7、三维空间直角坐标系中的基准点的确定，可按在三个直角坐标轴方向上各点位移绝 对值和为最小来求解各点的位移量，从而来确定一种相对稳定的坐标基准。设两期含有 t 个 网点的 GPS 监测网经过无约束网平差已获得各基线向量平差值，在其中选取互为独立且遍及全部 t 个网点的 t 1 条由首期 到第 K 期两期基 线向量平 差值之差 (K 1 , K1 , K 1 )ij ij ij ij ij1 1 1K K KX ij X ijYij YijZ ij Z ij，在此 (X , Y , Z )和 (X , Y, Z ij ) 分别表示第 1 期和第 K 期由第 i 至第 j 网点的三维坐标差平差值,且
8、i, j 1, t ，可用维数为3(t-1)1 的列向量 L 总括表示之。由基线向量与三维空间直角坐标之间的线性关系，即可得出由两期基线向量平差值之差估计位移向量的间接观测方程式AX = LK AX 1 = L(1)式(1)中 X 1 表示首期网点平差后三维空间直角坐标，因此解向量 X 即为由首期到第 K 期的三维位移量，表示为 X= (X K , Y K , Z K , X K , Y K , Z K .X K ,.Y K , Z K )T ，系1112 2 2t t t数阵 A 表示了 t-1 个基线向量与其两端点三维坐标的线性关系，可见维数为 3(t-1)3t 的系数阵 A 中的各元素必
9、为 1 或零。2.2 按改进的单纯形法的求解原理和方法为了以位移量绝对值总和最小为目标函数按照改进的单纯形法来求解，首先按照单纯形 法对变量有非负性的限制，将式(1)改写为AX + AX = L式(2)中 X + , X 0 ，(2)位移向量一次范数为最小为目标函数则为minCX + + CX + (t 1)C X (3)式(3)中， X 是为了获取初始基本容许解而引入的由 3t (1 t t 1) 个非负变量组成的人工向量，目标函数的价值系数 C, C 分别为 3t1 和 3 t 1 的行向量，其所有的分量都为 1。按照单纯形法的求解，将全部变量分成大于零的基变量所组成的 X B 和等于零的
10、非基变量所组成的 X N 两类向量,系数矩阵和价值系数行向量也按此划分，因此可以将由式(2)和 (3)组成的线性规划问题表述为以下等价形式minC B X B+ C N X NAB X B + AN X N = L(4)X B 0, , X N = 0基变量的系数矩阵 AB 的 3(t-1)1 个列向量必定须线性无关，由此可得到基本容许解X= A1 L A1 A X= A1 L = L相应的目标函数为BBX N = 0BNNB(5)B BBZ = C A 1 L = C L(6).在逐次换基迭代以求最优解的过程中,每得到一个基本容许解之后都要判断它是否为最优 解，为此须对每个非基变量计算出其判
11、别数j B Bjj = C C A1 A= C j A j(7)式(7) 中,j 遍及每个非基变量。A j 为此时非基变量所对应系数矩阵第 j 列，称 为对于 基 AB 的单纯形乘子。当所有的判别数 j 0 ，则表明己获得使目标函数最小的最优解，否则须由kj = min| j 0rk( i = 1,2,.,3(t 1)(10)由此确定 xr 为出基变量， AB 的第 r 列 A 为换入基。以 Ak 取代 Ar 得到新的基向量。每次换基都使得 AB1 改变一列，也就是相当于以 a为主元对整个系数矩阵和常数列向量 L进行消元，由此而得到新的基本容许解。再次求出判别数进行逐次换基迭代,直至所有的判别
12、数 j 0 而获得最优解。通过上述解算得出的结果中,在 t 个网点的 3t 个坐标分量中,作为基变量的坐标分量的个数为 3(t -1),仅有 3 个坐标分量作为非基变量而等于零，由此可将这 3 个非基变量作为固定 数据而作独立网平差。由于这 3 个零位移是在 3 个坐标轴的方向上，未必是在同一个点上， 仅一个坐标轴方向的位移为零，该点的另两个方向的位移却可能很大，因此不能据此判定即 是相对稳定点。当然可以比较这几个含有零位移分量的点位上其他的位移分量，考察各点总 位移的大小，但毕竟能否定出相对稳定点并无确切的把握。为弥补这个缺点，可将平面位移 表示在椭球面的新大地坐标系来表述平面位移。3 利用
13、新大地坐标系来寻求相对稳定点地形变分析常需要分解为垂直方向及水平位置两个方面来进行研究，在某些地区如城市， 垂直形变和水平形变往往因不同原因造成而存在有某种不一致性，由三维地心直角坐标并不 能直接表述垂直及水平方向的点位形变。而在上节中，在空间直角坐标系下寻求相对稳定点 也确有不便。在此转而采用新大地坐标系下来寻求相对稳定点，文献8提出,椭球面上一点 i沿其南北方问和东西方问的水平位移量可用新大地坐标系9中的纵向和横向坐标差分量( s, n s)，并证得它与基于大地坐标系所表示的 (M B , Ncos B L ) 是完全一致Li Bi Biiiiii的。为了得出新大地坐标形式下位移向量的间接
14、观测方程式,可先将三维直角坐标下的位移量分解表示成沿经线方向、纬线方向及大地高方向的位移量，即有X i Bi Yi = Ri Li (11)Z i H i 式(11)中，X i ,Yi ,Z i 为监测点 i 在空间直角坐标下的位移量，而Bi ,Li ,H i 则为点 i在大地坐标系中的水平和垂直方向的位移量。11 (M+ H )sinB cosL (N+ H ) cosB sinLcosB cosL i i i ii i i ii i Ri = (Mi + Hi )sinBi sinLi(Ni + Hi ) cosBi cosLicosBi sinLi (Mi + Hi ) cosBi0si
15、nBi 在上式中的 Li , Bi 分别为点 i 的经纬度，M i , N i 分别为点 i 的子午线曲率半径和卯酉圈曲率半径。由文献8所得出的微分关系式LidsLids= M i dBi , dBi =M(12)ids= N i cos Bi dLi , dLidsB niBinii =N i cos Bi(13)B0 为所设定的新大地坐标系原点的大地纬度,对应于原点的卯酉圈曲率半径为N0， ni 为点i的长度归化因子，求解公式可参看文献9。将式(12)和(13)代入式(11) 中，即得： (M i + H i ) sin Bi cos Li ( N i + H i ) sin Liicos
16、 Bi cos Li X M iN isLi (M i+ H i) sin Bisin Li( N i+ H i) cos Lii Yi = cos Bi sin Li ni sB Z M iN ii (M i + H i ) cos BiM i0sin Bi H isLii= F. ni sB (14)H i由式(14)，可将以空间直角坐标系中的位移量表示成新大地坐标系下水平和垂直方向的位移 量，以向量 Y 表示某一期对第 1 期的水平和垂直位移量并作为未知量,同样可列出观测方程BY = L(15)式(15)中右端的常数项 L 仍为平差后的两期 t1 个基线向量的三维坐标差值，而系数 矩阵
17、B 的构成应顾及式(14), 即有 F000 1AB = AF = 0 L 0F200 LF00LLt (16)式(16) 中的 F 阵为以 Fi 为各子块的分块对角阵，可见矩阵 A 的各元素已不再是仅为 1, 0 。如此遂将全部 GPS 网点的空间直角坐标系下位移量的观测方程(1)转换成新大地坐标下位移量的观测方程(15)。所得的位移量已分解成经线方向，纬线方向和大地高方向。按位移绝对值和为最小用上述的单纯形法来求解各点的水平及垂直位移量，从而来确定一种相对稳定的 坐标基准。4 实测上海 GPS 形变监测网的计算验证本文取用了共有 35 个点的上海地区 GPS 监测网于 02 年 12 月和
18、 03 年 6 月获得的两期 观测值分别经无约束平差的所得的各点的 WGS84 三维直坐标和三维大地坐标,由此即可得出34 条两期基线向量平差值之差(在此是选取 No.4 至其余各点的基线向量)。按三个直角坐标 轴方向上各点的位移绝对值和为最小,用单纯形法解得的各点的位移量如表 1 所示。再选定 监测网的中心附近一点(No.4)作为坐标原点建立起新大地坐标系，通过新大地坐标与大地坐 标的转换关系式，将各点的大地经纬度转换为新大地纵横坐标。利用式(14)和(16)可算得观 测方程式(15)中的系数阵 B。在新大地坐标系下对式(15) 用单纯形法进行位移计算，结果如表 2 所示。在此,若sL 为正
19、,表示向北位移,为负数则表示为向南位移；n sB 为正表示向东位 移，为负则表示向西位移；H 为正是上升，为负则表示沉降。表 1在空间直角坐标系中按位移分量一次范数最小所确定的位移分量No Xi（m） Yi（m） Zi（m）No Xi（m） Yi（m） Zi（m）10.00370.0030-0.0031180-0.0082-0.009720.0171-0.00250.0054190.0189-0.0078-0.01803-0.0072-0.0062-0.010520-0.00130.02050.00574-0.004000.0003210.01240.01860.006550.0305-0.0
20、060-0.033122-0.00570.03220.004360.01450.01740.010723-0.00700.0125-0.01237-0.0043-0.01190.002624-0.0072-0.00590.004880.01360.0003-0.013925-0.0177-0.0072-0.015190.0272-0.0280-0.029226-0.02200.03780.0264100.0020-0.0093-0.007927-0.0042-0.0034-0.021211-0.00680.0178028-0.0052-0.00970.0064120.00910.0189-0.
21、002229-0.0049-0.00650.0102130.00250.03600.006630-0.0004-0.0394-0.016514-0.02030.03920.006931-0.00860.01610.025215-0.00050.02770.0065320.0183-0.0210-0.0091160.01260.0184-0.0149330.0076-0.00380.0100170.00010.0178-0.0116340.0261-0.01880.0029i由表 1 可知，Xi,Yi ,Z i的 3 个零位移分量分别分布在 No.18, 4，11，这 3 个点位上，仅凭 1 个
22、位移分量为零, 当然不足以判断该点即为相对稳定点，还须参考其上另外两个坐标分量。因在这 3 个点中, No. 4 的 3 个位移分量相对较小，为-0.0040, 0, 0.0003，亦可将该 点取为固定点，即取该点的 3 个位移分量均为零,再经平移得出各点的位移向量。表 2 在新大地坐标系中按位移分量一次范数最小确定的位移分量NosL i(m)sB i(m)H i(m)NosL i(m)sB i(m)H i(m)1-0.00600.00780180.00370.04030.01952-0.0252-0.0003-0.012919-0.00180.0465-0.00983-0.00330.003
23、7-0.0225200.00180.01780.0086400.0004-0.0196210.0047-0.01530.006550.0251-0.0169-0.0394220.007800.00376-0.0254-0.0165-0.0274230.00740.01490.00687-0.01320.0188-0.032724-0.00480.0246-0.00188-0.01820.0422-0.0532250.00710.0086-0.04269-0.00150.0836-0.008626-0.00410.0055-0.0145100.01330.03220.0182270.01400.
24、0023-0.001411-0.0261-0.02870.016728-0.01590.0029-0.017812-0.0322-0.03580.010329-0.00360.0044-0.001413-0.0259-0.01900.012930-0.01010.0082-0.0479140.00180.03510.016131-0.00540.00230.000115-0.0013-0.02000.009232-0.00550.0029-0.0462160.0052-0.04640.015633-0.02590.0086-0.0182170.00710.03530.015534-0.0037
25、0.0124-0.0649L由表 5.4 的数据结果可知sii, ni s BH i的 3 个零位移分量分别在 No.4,22,1 点上,No.4点的水平位移量均很微小,但该点垂直位移却较大。垂直零位移在点 No.1 上No.31 亦接近于零，因此宜分别取 No.4 和 No.1 为水平位移和垂直位移的基准点，由此来确定相对稳定 基准。在此基准下解求各点的水平和垂直位移向量。为此可利用第 2 期 GPS 监测网 GPS 基线 向量的平差值，作为固定值的基准数据则取首期网平差后 No.4 点的大地经纬度及 No.1 点的 大地高，由此求出各点的大地经纬度和大地高，再将经纬度换算为新大地坐标，来与
26、第 1 期的新大地坐标相比较，从而得出各点的垂直位移与水平位移。更简便的做法是，可足够精确地认为表 2 所得的垂直位移H i 即为垂直位移正确值，从而得出第 2 期各点的大地高。由此不难得出第 2 期中各点的大地经纬度以至水平位移分量sin B jcos L jsin B jsin L jcos B j X X 0 B M + H4 jY4 jM + HM + H jjjjj 0 =j Y Y(17) L j 4 jsin L jcos L j4 jj4 jZ0 4 j 4 j4 j 0 ( N j + H j ) cos B j( N j + H j ) cos B j 4 jZ 4 j j
27、jj式(17)中, (XY4 jZ 4 j) 为水平基准点 No.4 至 No.j 第 2 期的基线向量平差值,(BLH )为 No.j 第 1 期平差后坐标, (X 00 Z 0 ) 为由第 1 期平差后坐标所得的基线向量。由此即可求得在此基准下第 2 期平差后大地经纬度变动量，再根据式 (12)和(13)即可求得新大地坐标系下的各点水平位移分量。将数据列于表 3 中表 3 固定 NO.4 为水平基准点求解新大地坐标系下的水平位移分量NosL(m)ni sB(m)NosL(m)ni sB(m)10.00710.008118-0.00190.008520.0112-0.0005190.0030
28、0.000630.00580.009920-0.00420.003240.00060.000221-0.0004-0.00755-0.0073-0.010222-0.00190.000960.0154-0.008723-0.00790.012270.00820.0029240.00660.022080.00170.001025-0.01090.024090.00490.0043260.01030.011710-0.00030.006227-0.00780.0181110.0043-0.0023280.01030.0135120.0023-0.0048290.01450.0136130.0005
29、-0.0080300.01320.034014-0.01850.0099310.02200.0117150.0032-0.0011320.01660.008416-0.0074-0.0075330.01900.004317-0.00770.0034340.02810.0003i i i i由各点的位移即可绘制出水平位移矢量(图 1)与地表沉降图(图 2)。由于新大地坐标是一种曲面坐标系，在平面图上,点位坐标可近似地采用 xi= ni s Bi, yi = s L 。（m）图 1 GPS 监测网点在新大地坐标系下的水平位移矢量图i20mm（m）图.1 中，箭头的指向表示第二期对第一期的点位的移动
30、方向，箭头的长度用来按比例尺表示位移的大小；而在图.2 中，向上指的箭头表示监测点垂直上升，向下指的箭头则表示监测 点垂直下降，以箭头长度表示沉降量的大小。图 2 GPS 监测网点在新大地坐标系下的沉降图（m）20mm（m）5 结语和建议以位移量一次范数最小为原则解算各监测点的位移量，可直接基于两期观测值各自确定 的最小二意义下的最优相对网形来确定相对稳定的位移分量，由于它毋须事先确定拟稳点， 可以避免由于假设失真而不能如实反映变形。水平位移与垂直位移之间虽有一定的联系，也有一定的独定性，尤其是在城市地区，同 时就 GPS 测量而言，垂直位置与水平位置在精度上目前还有差距，因此分别选取水平和垂
31、直 位置的基准点的新尝试，也许更符合于实际变形情况，而位移量以新大地坐标来表述更有利于进行这方面的分析和研究。参考文献1 周江文. 监测网的拟稳平差M. 北京：测绘出版社,1987.2 施一民. 论秩亏自由网平差的性质及稳健基准的意义J.同济大学学报，1991，19(3)：279-288.3 施一民.论测量控制网定位的各种处理方法J.同济大学学报，2002（11）：1331-13364 Caspary,W.Chen,Y.Q.,Knig,R.Kongruenzuntersuchungen in Deformationsnetzen durch Minimierung derSumme der K
32、laffungsbetrge.Schriftenreihe der Hochschule der Bundeswehr,Heft 9,Mnchen 1983 5 杨元喜.相关观测抗差估计.抗差估计论文集M.北京：测绘出版社，19926 施一民 监测网按一次范数最小的两阶段解算法J测量与地球物理集刊 1988，9 .11-18 7 卢开澄.组合数学算法与分析M.北京：清华大学出版社，19838 施一民.采用新型大地坐标系进行地形变分析的探索J.中国科技论文在线,2006 年 6 月，9 施一民，朱紫阳，范业明.坐标参数为长度量的一种新型的大地坐标系J. 同济大学学报， 2005,33（11）10
33、 施一民,朱紫阳，范业明.新型大地坐标系与大地坐标系之间的转换J.中国科技论文在线，2005 年 11月Establishment of Robust Datum of Analysisof Regional Crustal DeformationSHI Yi-min1,2 LUO Yan11.Department of Surveying and Geoinformatics,Tongji University,Shanghai, 200092,China；2.Key Laboratory of Modern Engineering Surveying, State Bureau of Su
34、rveying and Mapping, Shanghai, 200092, ChinaAbstractIn order to indicate actual displacement objectively and reasonably, a robust estimation with minimum of first order norm of displacement component is offered to determine displacement of relative stability, and it does not have to make any assumpt
35、ion for deformation of points, being different from ordinary quasi-stable adjustment. New form of geodetic coordinates is adopted to express displacement for effectively judging relative stable points with horizontal and vertical orientation. In this paper, according to this method, data of two stages of GPS deformation monitoring net in Shanghai are calculated to preliminarily prove that this method is feasible and effective for strain analysis.Key words： Crustal Deformation, horizontal displacements, robust datum, the new form of geodetic coordinate system