GPS從入門到放棄（十七） — 對流層延時

對流層概念

對流層（Troposphere）地球大氣層靠近地面的一層。它同時是地球大氣層裡密度最高的一層，它蘊含了整個大氣層約75%的品質，以及幾乎所有的水蒸氣及氣溶膠。

• 對流層的下界與地面相接，上界高度随地理緯度和季節而變化，在低緯度地區平均高度為 17 ~ 18 公裡，在中緯度地區平均為 10 ~ 12 公裡，在高次元地區平均為 8 ~ 9 公裡，并且夏季高于冬季。
• 衛星導航領域所說的對流層與嚴格的大氣科學中的對流層定義有所差別。衛星導航領域中需要考慮對流層延時，是以把從地球表面到電離層之間的區域都認為是對流層，即從0km到約50km的大氣層區域。
• 雲、霧、雨、雪等衆多天氣現象都發生在對流層。
• 對流層中主要的傳播方式或效應有：大氣折射、波導傳播、對流層散射、多徑傳播、大氣吸收，以及水汽凝結體和其他大氣微粒的吸收和散射。

對流層對衛星導航的影響

因為衛星都處于大氣層外，衛星信号傳播到地球上必須要經過大氣層中的對流層，于是必然要受到對流層的影響。

對流層是非彌散性媒體，即媒體的介電常數與頻率無關，是以在對流層中，不同頻率的電磁波具有相同的傳播速度。是以在電離層延時一文中用到的利用雙頻接收機消除電離層延時的方法對對流層就無法使用了。這使得我們在實際中一般采用數學模型來估計對流層延時。當然，對于有基站的情況，也可以通過差分的方式來消除對流層延時（可參考文章差分GPS）。

對流層模型

可能是由于對流層模型太多，準确度并不是那麼高，在GPS的接口說明文檔中并沒有給出官方的參考對流層模型。這裡我們以 RTKLIB 中用到的 Saastamoinen 模型說一下。

标準大氣層模型可以表示為：

p = 1013.25 × ( 1 − 2.2557 × 1 0 − 5 h ) 5.2568 p = 1013.25 \times (1-2.2557\times 10^{-5} h)^{5.2568} p=1013.25×(1−2.2557×10−5h)5.2568

t = 15.0 − 6.5 × 1 0 − 3 h + 273.15 t = 15.0 - 6.5 \times 10^{-3} h +273.15 t=15.0−6.5×10−3h+273.15

e = 6.108 × exp ⁡ [ 17.15 T − 4684.0 T − 38.45 ] × h r e l e = 6.108 \times \exp\left[\frac{17.15T-4684.0}{T-38.45}\right] \times h_{rel} e=6.108×exp[T−38.4517.15T−4684.0​]×hrel​

其中 h h h 是海拔高度， h r e l h_{rel} hrel​ 是相對濕度， p p p 為大氣壓力， e e e 為大氣中的水汽壓力， T T T 為大氣溫度。

一般對流層模型中，将對流層天頂總延遲（ZTD：Zenith Total Delay）分為對流層靜力學延遲（ZHD：Zenith Hydrostatic Delay）和對流層濕延遲（ZWD：Zenith Wet Delay）。靜力學延遲約占總延遲量的90%，可以通過實測氣壓和氣溫來計算，而濕延遲影響因素較多，不太容易估算。

根據 Saastamoinen 模型，對流層延時 T r T_r Tr​ 可以表示為靜力學延遲 T h T_h Th​ 與 濕延遲 T w T_w Tw​ 的和：

T h = 0.0022768   p 1.0 − 0.00266 cos ⁡ ( 2 ϕ ) − 0.00028   h × 1 0 − 3 × 1 cos ⁡ z T_h = \frac{0.0022768\:p}{1.0-0.00266\cos(2\phi)-0.00028\:h\times 10^{-3}} \times \frac{1}{\cos z} Th​=1.0−0.00266cos(2ϕ)−0.00028h×10−30.0022768p​×cosz1​

T w = 0.0022768 ( 1255 T + 0.05 ) e × 1 cos ⁡ z T_w = 0.0022768\left(\frac{1255}{T}+0.05\right)e \times \frac{1}{\cos z} Tw​=0.0022768(T1255​+0.05)e×cosz1​

T r = T h + T w T_r = T_h + T_w Tr​=Th​+Tw​

其中 z = π / 2 − e l z = \pi/2-el z=π/2−el 為天頂角， e l el el 為衛星仰角。

下面摘錄了 RTKLIB 中用 Saastamoinen 模型計算對流層延時的函數代碼，可以跟以上公式進行一一對應。需要注意的是計算時的機關換算，azel[0], azel[1]分别是方位角和仰角，pos[0], pos[1]分别是接收機緯度和經度，它們的機關都是弧度，pos[2]為接收機高度，機關為米。

/* 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);
    trpw=0.002277*(1255.0/temp+0.05)*e/cos(z);
    return trph+trpw;
}