本帖最后由 丈量大地 于 2012-6-27 20:48 编辑
介绍一种利用Excel表格解算七参数的方法 目前计算机的操作系统下都安装有Office2000---2003，利用Office下的Excel表格和表格中提供的多种矩阵运算函数，可以解算在GPS RTK作业中经常用到的由地心坐标系转换到国家坐标系所使用的坐标转换参数或称作七参数。利用EXCEL表格解算七参数的方法非常方便和实用且便于掌握。下面把解算过程中用到的有关矩阵运算的函数列出如下： 1. ＝Mmult(数组1，数组2) ――――――－两矩阵乘积函数 2. ＝Transpose(数组)――――――――――矩阵转置函数 3. ＝Minverse(数组)―――――――――― 矩阵求逆函数 4. ＝数组1 ± 数组2――――――――――矩阵加减函数 这里的数组其含义是：例如我们把一个矩阵存放到第3行～20行，第A列～G列，那么在进行矩阵运算时数组应输入A3：G20。 下面给出解算七参数的数学模型如下： Xi Xi Xi 0 -Z Y Εx/ρ ⊿X Yi ＝ Yi ＋ Yi K×10- 6 ＋ Z 0 -X ΕY/ρ ＋ ⊿Y Zi Zi Zi -Y X 0 ΕZ/ρ ⊿Z Ⅱ Ⅰ Ⅰ Ⅰ Xi Xi Li ＝ Yi － Yi （2） Zi Zi Ⅱ Ⅰ 角标Ⅱ表示地心空间大地直角坐标 角标Ⅰ表示参心空间大地直角坐标 Li ＝BiY （3） Y ＝[ ⊿X ⊿Y ⊿Z K Εx ΕY ΕZ ] （4） Vi ＝BiY － Li ( i = 1，… t ) t为公共点个数 （5） V ＝BY － L 所有误差方程 （6） Y ＝( BTB)- 1 BTL （7） 上述求解过程认为两种坐标系统的坐标是等权的，或不知道它们的精度信息的情况下所采用的解算公式。 若已知了公共点上两套坐标的精度信息，就可以确定其相应的权矩阵，来确定转换参数。若设
σ2xi 对称 σ2xi 对称 ΣⅡi ＝ σYX σ2Yi ΣⅠi ＝ σYX σ2Yi (8) σZX σZY σ2Zi σZX σZY σ2Zi ΣLi ＝ΣⅡi ＋ ΣⅠi (9) PLi ＝(ΣLi)- 1 (10) 转换参数的最小二乘解为： Y ＝ ( BTPB)- 1 BTPL (11) 精度评定： σ20 =( V TPV) / (3t－7) (12) σ2Y ＝σ20 ( BTPB)- 1 (13) 下面介绍利用Excel表格计算的过程如下： 1. 首先确定误差方程的系数矩阵B，B矩阵由下式确定：
1 0 0 Xi×10-6 0 －Zi/ρ Yi/ρ Bi ＝ 0 1 0 Yi×10-6 Zi/ρ 0 －Xi/ρ (14) 0 0 1 Zi×10-6 －Yi/ρ Xi/ρ 0 i ＝1…t, ρ＝206265 组成常数项矩阵L XⅡi －XⅠi Li ＝ YⅡi －YⅠi i＝1…t ZⅡi －ZⅠi 如果t＝3时，则B矩阵为9行×7列，L为9行×1列矩阵。 2. 求B矩阵的转置矩阵BT 利用=Transpose(数组)转置矩阵函数进行转置运算，首先在表格中确定将要转置的这个矩阵的行数和列数并把其范围确定下来。然后在输入栏中输入转置矩阵函数，数组输入B矩阵的所在单元格的行数和列数例如，=TANSPOSE(A3：I20)， 然后按CTRL+SHIFT+ENTER键，所得结果就显示在所划定的范围内。 3. 法方程组成（BTB）或（BTPB） 利用矩阵乘积函数=Mmult（数组1，数组2）组成法方程，首先在表格中确定乘积矩阵的行数和列数并把其范围确定下来。然后在输入栏中输入矩阵乘积函数，例如： = Mmult（A24：R30，A3：G20），然后按CTRL+SHIFT+ENTER键，所得结果就显示在所划定的范围内。 4. 法方程常数项矩阵组成（BTL）或（BTPL） 操作步骤与步骤3相同。 5. 法方程解算（矩阵求逆） （BTB）-1或（BTPB）-1 利用矩阵求逆函数=Minverse（数组）解算法方程，具体操作同步骤3， 例如： = Minverse（A33：G39）。 6. 求转换参数Y=（BTB）-1（BTL）或（BTPB）-1（BTPL） 操作步骤与步骤3相同。 精度评定： 7. 解算改正数矩阵V ＝BY － L 利用矩阵加法或减法函数进行计算，例如： =B55：B72－D55：D72 8. 求V的转置矩阵VT 计算步骤与步骤2相同 9. 计算VTV或VTPV 操作步骤与步骤3相同。利用公式（12）计算出单位权中误差σ20 10. 利用公式（13），计算转换参数的方差协方差矩阵，其中主对角线上的元素 分别为所解算的七个参数的精度 | 
雷达卡
发表于 2012-6-27 20:45
提升卡
沉默卡
喧嚣卡
变色卡
抢沙发
显身卡
发表于 2013-6-20 09:04
