对流层折射一般泛指非电离层大气对电磁波的折射,其中非电离大气包括对流层和平流层
1.RTKLIB改正模型
1.1 saas模型
RTKLIB中使用的第一个改正模型为Saastamoinen
模型,其公式如下:
对应manualP149
对流层模型
标准大气可以表示为:
其中p为总气压(hPa), T为空气的绝对温度(K), h为MSL(平均海平面)以上的大地高度,e为水汽分压(hPa), hrel为相对湿度。
对流层天顶总延迟(ZTD:Zenith Total Delay)分为对流层静力学延迟(ZHD:Zenith Hydrostatic Delay)和对流层湿延迟(ZWD:Zenith Wet Delay)
对流层延时 Tr 可以表示为静力学延迟 Th 与 湿延迟 Tw 的和:
T
h
=
0.0022768
p
1.0
−
0.00266
cos
(
2
ϕ
)
−
0.00028
h
×
1
0
−
3
×
1
cos
z
T
w
=
0.0022768
(
1255
T
+
0.05
)
e
×
1
cos
z
T
r
=
T
h
+
T
w
\begin{gathered} T_h=\frac{0.0022768 p}{1.0-0.00266 \cos (2 \phi)-0.00028 h \times 10^{-3}} \times \frac{1}{\cos z} \\ T_w=0.0022768\left(\frac{1255}{T}+0.05\right) e \times \frac{1}{\cos z} \\ T_r=T_h+T_w \end{gathered}
Th=1.0−0.00266cos(2ϕ)−0.00028h×10−30.0022768p×cosz1Tw=0.0022768(T1255+0.05)e×cosz1Tr=Th+Tw
对流层延迟Trs用Saastamoinen模型表示,p、T、e由标准大气导出。其中 z=π/2−ELrs 为天顶角,ELrs 为卫星仰角
对应代码如下:
/* troposphere model -----------------------------------------------------------
* compute tropospheric delay by standard atmosphere and saastamoinen model
* args : gtime_t time I time
* double *pos I receiver position {lat,lon,h} (rad,m)
* double *azel I azimuth/elevation angle {az,el} (rad)
* double humi I relative humidity
* return : tropospheric delay (m)
*-----------------------------------------------------------------------------*/
extern double tropmodel(gtime_t time, const double *pos, const double *azel,
double humi)
{
const double temp0=15.0; /* temparature at sea level */
double hgt,pres,temp,e,z,trph,trpw;
if (pos[2]<-100.0||1E4<pos[2]||azel[1]<=0) return 0.0;
/* standard atmosphere */
hgt=pos[2]<0.0?0.0:pos[2];
pres=1013.25*pow(1.0-2.2557E-5*hgt,5.2568);
temp=temp0-6.5E-3*hgt+273.16;
e=6.108*humi*exp((17.15*temp-4684.0)/(temp-38.45));
/* saastamoninen model */
z=PI/2.0-azel[1];
trph=0.0022768*pres/(1.0-0.00266*cos(2.0*pos[0])-0.00028*hgt/1E3)/cos(z);//静力学延迟 Th
trpw=0.002277*(1255.0/temp+0.05)*e/cos(z);//湿延迟 Tw
return trph+trpw;
}
其误差的标准差为0.3
1.2 sbas模型
RTKLIB中使用的第二个改正模型为SBAS
模型,其定义如下:
对应manualP149
SBAS对流层模型
如果将处理选项ʺ对流层校正ʺ设置为ʺSBASʺ,则应用SBAS接收机规范中定义的SBAS对流层模型。该模型通常称为ʺMOPS模型ʺ。详情请参考[8]A.4.2.4。
对应代码如下:
/* tropospheric delay correction -----------------------------------------------
* compute sbas tropospheric delay correction (mops model)
* args : gtime_t time I time
* double *pos I receiver position {lat,lon,height} (rad/m)
* double *azel I satellite azimuth/elavation (rad)
* double *var O variance of troposphric error (m^2)
* return : slant tropospheric delay (m)
*-----------------------------------------------------------------------------*/
extern double sbstropcorr(gtime_t time, const double *pos, const double *azel,
double *var)
{
const double k1=77.604,k2=382000.0,rd=287.054,gm=9.784,g=9.80665;
static double pos_[3]={0},zh=0.0,zw=0.0;
int i;
double c,met[10],sinel=sin(azel[1]),h=pos[2],m;
trace(4,"sbstropcorr: pos=%.3f %.3f azel=%.3f %.3f\n",pos[0]*R2D,pos[1]*R2D,
azel[0]*R2D,azel[1]*R2D);
if (pos[2]<-100.0||10000.0<pos[2]||azel[1]<=0) {
*var=0.0;
return 0.0;
}
if (zh==0.0||fabs(pos[0]-pos_[0])>1E-7||fabs(pos[1]-pos_[1])>1E-7||
fabs(pos[2]-pos_[2])>1.0) {
getmet(pos[0]*R2D,met);
c=cos(2.0*PI*(time2doy(time)-(pos[0]>=0.0?28.0:211.0))/365.25);
for (i=0;i<5;i++) met[i]-=met[i+5]*c;
zh=1E-6*k1*rd*met[0]/gm;
zw=1E-6*k2*rd/(gm*(met[4]+1.0)-met[3]*rd)*met[2]/met[1];
zh*=pow(1.0-met[3]*h/met[1],g/(rd*met[3]));
zw*=pow(1.0-met[3]*h/met[1],(met[4]+1.0)*g/(rd*met[3])-1.0);
for (i=0;i<3;i++) pos_[i]=pos[i];
}
m=1.001/sqrt(0.002001+sinel*sinel);
*var=0.12*0.12*m*m;
return (zh+zw)*m;
}
1.3 precise 模型
RTKLIB中使用的第三个改正模型为精密对流层模型,其定义如下:
对应manualP149
精密对流层模型
如果将处理选项ʺ对流层校正ʺ设置为ʺEstimate ZTDʺ或ʺEstimate ZTD+Gradʺ,则应用更精确的对流层模型,严格映射函数为:
RTKLIB采用Saastamoinen模型(E.5.4)给出了天顶角z=0和相对湿度hrel=0的对流层天顶静液延迟。对于映射函数,RTKLIB默认使用NMF (Niell映射函数)[70]。在参数估计过程中,天顶总延迟ZT,r和梯度参数GN,r,GH,r(在ʺEstimate ZTD+Gradʺ的情况下)作为未知参数进行估计。对于映射函数,RTKLIB可以通过设置编译器选项- DIERS_MODEL since ver来使用GMF[71]。2.4.2。
对应代码如下:
/* precise tropospheric model ------------------------------------------------*/
static double trop_model_prec(gtime_t time, const double *pos,
const double *azel, const double *x, double *dtdx,
double *var)
{
const double zazel[]={0.0,PI/2.0};
double zhd,m_h,m_w,cotz,grad_n,grad_e;
/* zenith hydrostatic delay */
zhd=tropmodel(time,pos,zazel,0.0);
/* mapping function */
m_h=tropmapf(time,pos,azel,&m_w);
if (azel[1]>0.0) {
/* m_w=m_0+m_0*cot(el)*(Gn*cos(az)+Ge*sin(az)): ref [6] */
cotz=1.0/tan(azel[1]);
grad_n=m_w*cotz*cos(azel[0]);
grad_e=m_w*cotz*sin(azel[0]);
m_w+=grad_n*x[1]+grad_e*x[2];
dtdx[1]=grad_n*(x[0]-zhd);
dtdx[2]=grad_e*(x[0]-zhd);
}
dtdx[0]=m_w;
*var=SQR(0.01);
return m_h*zhd+m_w*(x[0]-zhd);
}
其中调用了tropmapf
函数,而它实际上就是nmf
投影函数
信号传播路径上的对流层延迟STD与测站天顶方向的对流层延迟ZTD间有:STD=m*ZTD;m为映射函数
其误差的标准差为0.01
2. SPP与PPP中对应计算
2.1 SPP
SPP中对应函数pntpos.c-rescode
中对对流层进行选择模型改正
/* tropospheric correction */
if (!tropcorr(time, nav, pos, azel + i * 2, opt->tropopt, &dtrp, &vtrp))
计算对流层延时代码
extern int tropcorr(gtime_t time, const nav_t *nav, const double *pos,
const double *azel, int tropopt, double *trp, double *var)
{
trace(4,"tropcorr: time=%s opt=%d pos=%.3f %.3f azel=%.3f %.3f\n",
time_str(time,3),tropopt,pos[0]*R2D,pos[1]*R2D,azel[0]*R2D,
azel[1]*R2D);
/* Saastamoinen model */
if (tropopt==TROPOPT_SAAS||tropopt==TROPOPT_EST||tropopt==TROPOPT_ESTG) {
*trp=tropmodel(time,pos,azel,REL_HUMI);
*var=SQR(ERR_SAAS/(sin(azel[1])+0.1));
return 1;
}
/* SBAS (MOPS) troposphere model */
if (tropopt==TROPOPT_SBAS) {
*trp=sbstropcorr(time,pos,azel,var);//标准差为0.12
return 1;
}
/* no correction */
*trp=0.0;
*var=tropopt==TROPOPT_OFF?SQR(ERR_TROP):0.0;//ERR_TROP=3
return 1;
}
根据配置对流层选择对应模型分别是tropmodel(saas)和sbstropcorr
2.2 PPP
PPP中对应函数ppp.c-ppp_Res
中对对流层进行选择模型改正
/* tropospheric and ionospheric model */
if (!model_trop(obs[i].time,pos,azel+i*2,opt,x,dtdx,nav,&dtrp,&vart)||
!model_iono(obs[i].time,pos,azel+i*2,opt,sat,x,nav,&dion,&vari)) {
continue;
}
计算对流层延时代码
static int model_trop(gtime_t time, const double *pos, const double *azel,
const prcopt_t *opt, const double *x, double *dtdx,
const nav_t *nav, double *dtrp, double *var)
{
double trp[3]={0};
if (opt->tropopt==TROPOPT_SAAS) {
*dtrp=tropmodel(time,pos,azel,REL_HUMI);
*var=SQR(ERR_SAAS);//ERR_SAAS=0.3 标准差为0.3
return 1;
}
if (opt->tropopt==TROPOPT_SBAS) {
*dtrp=sbstropcorr(time,pos,azel,var);//标准差为0.12
return 1;
}
if (opt->tropopt==TROPOPT_EST||opt->tropopt==TROPOPT_ESTG) {
matcpy(trp,x+IT(opt),opt->tropopt==TROPOPT_EST?1:3,1);
*dtrp=trop_model_prec(time,pos,azel,trp,dtdx,var);//标准差为0.01
return 1;
}
return 0;
}