等位面上逆Stokes_Hotine_Meinesz运算数值积分Fortran代码逆Stokes_Hotine_Vening Meinesz运算 二维码
260
发表时间:2024-12-19 10:40网址:http://www.zcyphygeodesy.com ![]() [计算目标] 由等位面大地高格网(m)及其面上残差高程异常(m)或垂线偏差向量(S,W,″)格网,按逆Stokes、Hotine或Vening-Meinesz运算数值积分方法,计算等位面上其他类型残差扰动重力场元。 扰动重力场元逆运算积分属Stokes边值问题,要求被积高程异常或垂线偏差向量(S,W)位于等位面,通常用于海洋卫星测高重力场反演计算。 为实现有限半径积分,通常需采用参考重力场移去恢复法,先移去等位边界面上的模型空间异常/扰动重力,再积分得到计算点处的残差高程异常,最后恢复计算点处的模型高程异常。 等位面可采用参考重力场模型(不大于360阶)构造,在高度不大于10千米的近地空间,可用等正(常)高面大地高格网表示。 [测试入口程序] InvStokesHotineMeinesznumintg.f90 运算类型参数knd-=0逆Stokes运算,=1逆Hotine运算,=2逆Vening-Meinesz运算。 输入参数dr-积分半径(m)。 输入参数calcpntfl-计算点文件名。头文件一行,记录格式:点号/点名 经度(度小数) 纬度(度小数)......。 输入参数dwmhgrdfl-等位面大地高格网文件名。 输入参数gravgrdfl-等位面上残差高程异常(knd=0或1)或垂线偏差向量(knd=2)格网文件名。 [主要调用模块] (1)扰动重力场元逆运算数值积分算法模块 InversenumIntegral(calcpntfl,dwmhgrdfl,gravgrdfl,knd,dr) 输出文件reslt.txt。文件记录格式: 当knd=0时,在输入文件记录的基础上,增加一列由等位面大地高格网内插得到的计算点大地高,和一列该点的残差空间异常积分值。 当knd=1时,在输入文件记录的基础上,增加一列由等位面大地高格网内插得到的计算点大地高,和一列该点的残差扰动重力积分值。 当knd=3时,在输入文件记录的基础上,增加一列由等位面大地高格网内插得到的计算点大地高,和3列该点的残差高程异常、残差扰动重力和残差空间异常积分值。 (2)扰动场元径向梯度数值积分算法模块 real*8 function RaGradientBLH(BLH,gra,dwm,nlat,nlon,hd,dr,GRS) 输入BLH(3)-计算点的经纬度(度小数)和大地高(m)。 输入dwm(nlat,nlon)-等位边界面大地高格网,用于精确计算积分面元与计算点的积分距离。 输入gra(nlat,nlon)-等位面上残差扰动场元格网。 输入dr, hd(6)-积分半径(m)和格网规格参数(最小最大经度,最小最大纬度,经度间隔,纬度间隔)。 输入GRS(6)-gm, ae, j2, omega, 1/f, 缺省值。 返回-残差扰动场元径向梯度积分值。单位等于残差扰动场元gra单位/m,方向为向径方向(由地心指向地球外部)。 (3)逆Vening-Meinesz运算数值积分算法模块 AnVMdeflBLH(BLH,kai,nta,dwm,nlat,nlon,hd,dr,rst,GRS) 输入kai(nlat,nlon)-等位面上残差垂线偏差南向分量格网(弧度)。 输入nta(nlat,nlon)-等位面上残差垂线偏差西向分量格网(弧度)。 返回rst(3)-计算点的残差高程异常(m)、残差空间异常(m/s2)和残差扰动重力(m/s2)。 (4)场元一阶水平梯度向量估计模块 HGradientOne(BLH,rga,hgt,nlat,nlon,hd,tx,GRS) 输入rga(nlat,nlon)-场元数值格网。 输入hgt(nlat,nlon)-场元所在面的大地高格网。 返回参数tx(2)-计算点BLH处场元rga的1阶水平梯度向量估值,tx(1)北向,tx(2)东向,单位等于场元gra单位/m。 (5)正常重力场元计算模块 normdjn(GRS,djn);GNormalfd(BLH,NFD,GRS) 返回NFD(5)-正常重力位,正常重力,正常重力梯度,正常重力线方向,正常梯度方向。 (6)勒让德函数及其导数计算模块 LegPn_dt2(pn,dp1,dp2,n,t) 计算勒让德函数Pn(t)及其对ψ一、二阶导数t=cosψ。 (7)大型法方程解算模块库 EqHestens(BPB,xx,nn,BPL);EqJordon(BPB,xx,nn,BPL) EqCholesky(BPB,XX,nn,BPL);EqueSVD(BPB,XX,nn,BPL) RidgeEstimate(BPB,xx,nn,BPL);…… (8)大地坐标形式变换包 BLH_RLAT(GRS, BLH, RLAT);BLH_XYZ(GRS, BLH, XYZ) RLAT_BLH(GRS, RLAT, BLH) (9)其他辅助模块 CGrdPntD2(lon,lat,dt,row,col,hd);PickRecord(line, kln, rec, nn) [编译连接] Fortran固定格式代码,fortran编译器,外部连接库mkl_lapack95_ilp64.lib。 [算法公式] PAGravf4.5说明书 1.4.1大地测量数据文件格式约定 7.9.3扰动重力场元逆运算积分公式 7.1(4)低阶勒让德函数及其一、二阶导数算法 附件rar压缩包:visual studio_intel fortran 集成环境测试项目、DOS可执行测试程序和全部测试输入输出数据。 |