大气延迟

对流层延迟

在GNSS领域，对流层延迟是定位误差来源之一。卫星导航定位中的对流层延迟通常是泛指电磁波信号在通过高度为50km以下的未被电离的中性大气层时所产生的信号延迟。在研究信号延迟的过程中，我们不再将该大气层细分为对流层和平流层（如大气科学中那样），也不再顾及两者之间性质上的差别。

由于80%的延迟发生在对流层，所以我们将发生在该中性大气层中的信号延迟通称为对流层延迟。 对流层是大气层较低的部分，对于直到高达15GHz的频率来说它是非色散的。在这种介质中，与L1和L2上GPS载波和信号信息（PRN码和导航数据）相关联的相速和群速，相当于*空间传播被同等地延迟了。这种延迟随对流层折射率而变，其折射率取决于当地的温度、压力和相对湿度。

如果不补偿，这种延迟的等效距离可从卫星在天顶和用户在海平面上的2.4m左右到卫星在约5°仰角上的25m左右。

电离层延迟

电离层概念

电离层（Ionosphere）是地球大气的一个电离区域。它是受到太阳高能辐射以及宇宙线的激励而电离的大气高层。50千米以上的整个地球大气层都处于部分电离或完全电离的状态，电离层是部分电离的大气区域，完全电离的大气区域称磁层。

• 电离层的范围从离地面约50公里开始一直伸展到约1000公里高度的地球高层大气空域。
• 电离层的主要特性由电子密度、电子温度、碰撞频率、离子密度、离子温度和离子成分等空间分布的基本参数来表示。
• 电离层的研究对象主要是电子密度随高度的分布。电子密度（或称电子浓度）是指单位体积的*电子数，随高度的变化与各高度上大气成分、大气密度以及太阳辐射通量等因素有关。电离层内任一点上的电子密度，决定于上述*电子的产生、消失和迁移三种效应。在不同区域，三者的相对作用和各自的具体作用方式也大有差异。
• 电离层中存在相当多的*电子和离子，能使无线电波改变传播速度，发生折射、反射和散射，产生极化面的旋转并受到不同程度的吸收。

电离层延迟误差

电离层延迟误差（ionospheric delay error）亦称电离层折射误差(ionospheric refraction error)，是由电离层效应引起的观测值误差。

与其他电磁波一样，当导航卫星信号通过电离层时，导航卫星信号的路径也要发生弯曲，传播速度也会发生变化。电离层延迟对距离测量的影响，在天顶方向(仰角0度)最大可达50m，在近地平方向时(仰角20度)可达150m，因此必须仔细地加以改正，否则将严重影响伪距观测精度。

电离层折射可用三种方法来消除：

一是利用电离层模型加以改正；

二是利用同步观测值求差，这种方法对于短基线的效果尤为明显；

三是利用双频观测值，即通过不同频率的观测值组合来对电离层的延迟进行改正。

SPP中模型改正方法

声明

/*-----------------------------解码模块;-----------------------*/

/*星历结构体;*/
struct  GPSEPHEM
{
	unsigned long PRN;		//Satallite PRN number;

	GPSTIME TOC;
	GPSTIME TOE;
	
	/*第一数据块(第一子帧）;*/
	unsigned long week;		//GPS reference week number;
	double URN;				//User Range Accurancy variance,m2;
	unsigned long health;	//Health status;
	double tgd;				//Estimated group delay difference,seconds;
	unsigned long idoc;		//Issue of data clock;
	double a_f0;			//Clock aging paramenter,seconds(s);
	double a_f1;			//Clock aging paramenter(s/s);
	double a_f2;			//Clock aging paramenter(s/s/s);

	/*第二数据块（第二、三子帧）;*/
	//星历参考时刻toe时刻的轨道根数;
	double M0;				//Mean anomaly of reference time,radians;
	double w;				//Argument of perigee,radians;近地点角距;
	double I0;				//Inclination angle at reference time,radian;
	double w0;				//Right ascension,radians;
	double A;				//Semi-major axis,metres;半长轴;
	double ecc;				//Eccentricity,dimensionless-quantity defined for a conic section;
	double tow;				//Time stamp of subframe;子帧的时间标记;

	//9个轨道摄动参数;
	double delta_N;			//Mean motion difference,radians/second;
	double V_I;				//Rate of inclination angle,radians/second;
	double V_w0;			//Rate of Right ascension,radians;
	double cuc;				//Argument of latitude(amplitude of cosine,radians);
	double cus;				//Argument of latitude(amplitude of sine,radians);
	double crc;				//Orbit radius(amplitude of cosine,metres);
	double crs;				//Orbit radius(amplitude of sine,metres);
	double cic;				//Inclination(amplitude of cosine,radian);
	double cis;				//Inclination(amplitude of sine,radian);

	//包含的参数;
	unsigned long IODE1;
	unsigned long IODE2;
	unsigned long z_week;	//Z count week number;
	double toe;				//Reference time for ephemeris ,seconds;
	double toc;				//SV clcok correction term,seconds;
	int AS;					//0=FALSE;1=TRUE;
	double N;				//Corrected mean motion,radians/seconds;	

};

/*Ion与UTC解码结构体;*/
/*Message_ID==8 IONUTC Ionospheric and UTC Data*/
struct IONOSPHERIC
{
	double a0;	//Alpha parameter constant term;
	double a1;	//Alpha parameter 1st order term;
	double a2;	//Alpha parameter 2nd order term;
	double a3;	//Alpha parameter 3rd order term;
	double b0;	//Beta parameter constant term;
	double b1;	//Beta parameter 1st order term;
	double b2;	//Beta parameter 2nd order term;
	double b3;	//Beta parameter 3rd order term;
	IONOSPHERIC()
	{
		a0 = a1 = a2 = a3 = b0 = b1 = b2 = b3 = 0;
	}
};

struct UTC
{
	unsigned long utc_wn;	//UTC reference week number;
	unsigned long tot;		//Reference time of UTC parameters;
	double A0;				//UTC constant term of polynomial;
	double A1;				//UTC 1st order term of polynomial;
	unsigned long wn_lsf;	//Future week number;
	unsigned long dn;		//Day number(1=Sunday and 7=Saturday);
	long deltat_ls;			//Delta time due to leap second;
	long deltat_lsf;		//Future delta time due to leap second;
	unsigned long deltat_utc;//Time difference;
};

/*接收机观测值结构体;*/
/*Message_ID==43 RANGE Satellite Range Information*/
struct RANGE
{
	long obsnum;			//Number of oberservations with information to follow;
	unsigned short PRN_slo;	//GPS:1-32,SBAS:120-138,GLONASS:38-61;
	unsigned short glofreq;	//GLONASS Frequency+7;
	double psr;				//Pseudorange measurement(m);
	float psr_std;			//Pseudorange measurement standard deviation(m);
	double adr;				//Carrier phase,in cycles(accumulated Doppler range);
	float adr_std;			//Estimated carrier Doppler frequency(Hz);
	float dopp;				//Instantaneous carrier Doppler frequency;
	float C_No;				//Carrier to nosie density radio;C/No=10[log10(S/No)](dB-Hz);
	float locktime;			//#of second of continuous tracking(no cycle slipping);
	unsigned long ch_tr_status;	//Tracking status;

};

/*存储接收机自己解算的位置结构体;*/
/*Message_ID==496 PDPPOS PDP filter position*/
struct PDPPOS
{
	int sol_status;		//Solution status;
	int pos_type;		//Position type ;
	double lat;			//Latitude;
	double lon;			//Longitude;
	double hgt;			//Height above mean sea level;海拔高 ellipsoid(m);
	float undulation;	//Undulation-the relationship between the geoid and the WGS84 ellipsoid;
	int datum_id;		//Datum ID number;
	float lat_delta;	//Latitude standard deviation;
	float lon_delta;	//Longitude standard deviation;
	float hgt_delta;	//Height standard deviation;
	float diff_age;		//Differential age in seconds;
	float sol_age;		//Solution age in seconds;
	unsigned char sats;			//Number of satellite vehicles tracked;
	unsigned char sats_soln;	//Number of satellites in the solution;

};

/*每个历元的观测值数据;*/
struct EPOCHOBS
{
	GPSTIME	Time;					//GPS time;
	short SATnum[CHANNELNUM];		//which SAT	
	long Obsnums;					//Number of observation with information to follow;
	RANGE RangeObs[CHANNELNUM];		//Number of Channels;
	unsigned long ReceiverStatus;	//Receiver Status;
	EPOCHOBS()
	{
		Obsnums = 0;
		ReceiverStatus = -1;
	}
};

/*存储所有数据，包括星历，观测值等的结构体;*/
struct RAW
{
	EPOCHOBS EpochObs;				//每个历元的观测数据;
	GPSEPHEM EphemObs[SATNUM];		//卫星星历数据;
	IONOSPHERIC IonObs;				//电离层数据;
	UTC	UTCObs;						//时间数据;
	PDPPOS PDPpos;
};

/*对流层气象参数结构体，包含在SATPOS结构体中;*/
struct TROPOSPHERE
{
	double T;	//测站上的干温;
	double p;	//测站上的气压;
	double RH;	//测站上的相对湿度;
	double H;	//测站高度;

	TROPOSPHERE()
	{
		T = 0.0;	//测站上的干温;
		p = 0.0;	//测站上的气压;
		RH = 0.0;	//测站上的相对湿度;
		H = 0.0;	//测站高度;
	}
};

/*存储接收机解算结果的结构体;*/
struct RECIVIERPOS
{
	GPSTIME Time;	//什么时刻的位置;

	XYZ XYZPOS;	//空间直角坐标系下的位置;
	BLH BLHPOS;	//大地坐标系下位置;
	double Velocity[3];	//接收机速度;
	double delta_t;	//接收机钟差;
	double Vdelta_t;//接收机钟速;

	TROPOSPHERE TropPara;

	double PDOP;
	double DOP;
	
	double sigmaP;
	double sigmaV;


	RECIVIERPOS()
	{
		DOP = 0;
		PDOP = 0;
		delta_t=0;	//接收机钟差;
		Vdelta_t=0;//接收机钟速;
		Velocity[0] = Velocity[1] = Velocity[2] = 0.0;
		sigmaV = sigmaP = 0;
	}
};

/*解算卫星星历过程中得到的卫星轨道参数;*/
struct SATORBPAR
{
	double n0;
	double n;
	double tk;
	double Mk;
	double Ek;
	double vk;
	double PHIk;

	double delta_uk;
	double delta_rk;
	double delta_ik;

	double uk;
	double rk;
	double ik;

	double xk1;
	double yk1;
	double wk;

	double xk;
	double yk;
	double zk;

	double e;
	double A;

};

/*存储卫星星历解算结果的结构体;*/
struct SATPOS
{
	short PRN;
	GPSTIME Time;
	XYZ XYZPOS;
	BLH BLHPOS;
	EAR O2S_EARPOS;
	double delta_Ion;
	double delta_Trop;
	double SignalTime;
	double Velocity[3];
	double delta_t;
	double Vdelta_t;
	SATPOS()
	{
		Velocity[0] = Velocity[1] = Velocity[2] =  0.0;
		delta_t = 0;
		Vdelta_t = 0;
		SignalTime = 0;
		delta_Ion = 0;
		delta_Trop = 0;
	}
};

实现

/*-----------------------电离层对流层误差解算结果;-------------------------*/
//模块功能：用以解算对流层与电离层对定位造成的误差干扰，并给出改正数;
//函数包含：;
//	Klobuchar:电离层改正函数;
//	Hopefiled:对流层改正函数;
//	公式来源：《GPS测量与数据处理》中关于电离层与对流层改正部分，王老师所给PPT;
/*------------------------------------------------------------------------*/
#include "Fuction.h"

/*-------------------------电离层改正模型;---------------------------------*/
//函数功能：计算电离层对定位的影响;
//函数输入：;
//	RawObs:接收机输出数据;
//	SatPos:卫星位置结构体;
//	RecivierPos:接收机位置结构体;
//	delta_Ion:电离层延迟;
//公式来源：《GPS测量与数据处理》中关于电离层与对流层改正部分，王老师所给PPT;
/*--------------------------------------------------------------------------*/
bool Klobuchar(RAW*RawObs, SATPOS *SatPos,RECIVIERPOS *RecivierPos, double *delta_Ion)
{
	ENU ENUPos;
	XYZ2BLH(RecivierPos->XYZPOS, &RecivierPos->BLHPOS);
	XYZ2ENU(RecivierPos->XYZPOS, SatPos->XYZPOS, &ENUPos);
	ENU2EAR(ENUPos, &SatPos->O2S_EARPOS);

	double psi, phi_I, lambda_I, phi_m, t = 0;
		double A_I,P_I,F_I,X_I=0;
	//计算测站点在地心的夹角;
	psi = 0.0137 / (SatPos->O2S_EARPOS.E/PI + 0.11) - 0.022;

	phi_I = RecivierPos->BLHPOS.B / PI + psi*cos(SatPos->O2S_EARPOS.A);
	if (phi_I > 0.461)phi_I = 0.461;
	if (phi_I < -0.461)phi_I =- 0.461;

	//升交点经度;
	lambda_I = RecivierPos->BLHPOS.L / PI + (psi*sin(SatPos->O2S_EARPOS.A)) / cos(phi_I*PI);

	//地磁维度;
	phi_m = phi_I + 0.064*cos((lambda_I - 1.617)*PI);

	//计算地方时;
	t = 43200 * lambda_I + RawObs->EpochObs.Time.SecOfWeek;
	t = t - (int)(t / 86400.0) * 86400;
	while (t<0 || t>86400)
	{
		if (t >= 86400)t = t - 86400;
		if (t < 0)t = t + 86400;
	}

	//振幅;
	A_I = RawObs->IonObs.a0 + RawObs->IonObs.a1*phi_m + RawObs->IonObs.a2*phi_m*phi_m +
		RawObs->IonObs.a3*phi_m*phi_m*phi_m;
	if (A_I < 0)A_I = 0.0;

	//周期;
	P_I = RawObs->IonObs.b0 + RawObs->IonObs.b1*phi_m + RawObs->IonObs.b2*phi_m*phi_m +
		RawObs->IonObs.b3*phi_m*phi_m*phi_m;
	if (P_I < 72000)P_I = 72000;

	//相位延迟;
	X_I = 2 * PI*(t - 50400) / P_I;

	//倾斜因子;
	F_I = 1.0 + 16.0*(0.53 - SatPos->O2S_EARPOS.E / PI)*(0.53 - SatPos->O2S_EARPOS.E / PI)
		*(0.53 - SatPos->O2S_EARPOS.E / PI);

	//计算电离层延迟;
	if (fabs(X_I) <= 1.57)
	{
		*delta_Ion = (5e-9 + A_I*cos(X_I))*F_I;
	}
	else
	{
		*delta_Ion = 5e-9*F_I;
	}
	*delta_Ion = *delta_Ion*WGS84_V_Light;

	return true;
}
/*-----------------------------Hopefield 模型;----------
函数功能：用于计算对流层改正参数;
函数输入：;
	SatPos:卫星位置结构体;
	ReciverPos:接收机位置结构体;
	Trop:接收机位置气象参考;
	delta_Trop:对流层改正;
	公式来源：《GPS测量与数据处理》中关于电离层与对流层改正部分，王老师所给PPT;
	----------------------------------------------------*/

bool Hopefield(SATPOS *SatPos, RECIVIERPOS ReciverPos, TROPOSPHERE Trop,double *delta_Trop)
{
	ENU Obs2Sat;
	BLH BLHObsPos;
	XYZ2ENU(ReciverPos.XYZPOS, SatPos->XYZPOS, &Obs2Sat);
	ENU2EAR(Obs2Sat, &SatPos->O2S_EARPOS);
	XYZ2BLH(ReciverPos.XYZPOS, &BLHObsPos);
	
	if (SatPos->O2S_EARPOS.E < 0){ *delta_Trop = 0; return false; }
	//标准气象元素;
	double H0, T0, p0, RH0 = 0;
	H0 = 0;
	T0 = 15 + 273.15;//单位是开尔文;
	p0 = 1013.25;
	RH0 = 0.5;

	double RH, p, T, e, hw, hd, Kw, Kd = 0;
	Trop.H = BLHObsPos.H;

	Trop.RH = RH0*exp(-0.0006396*(Trop.H - H0));

	Trop.p = p0*pow(1 - 0.0000226*(Trop.H - H0)*1.0, 5.225);
	Trop.T = T0 - 0.0065*(Trop.H - H0);
	e = Trop.RH*exp(-37.2465 + 0.213166*Trop.T - 0.000256908*Trop.T*Trop.T);

	hw = 11000.0;
	hd = 40136 + 148.72*(T0 - 273.16);
	Kw = (155.2e-7)*4810.0*e*(hw - Trop.H) / (Trop.T*Trop.T);
	Kd = (155.2e-7)*Trop.p*(hd - Trop.H) / Trop.T;
	
	double delta_d, delta_w = 0;
	delta_d = Kd / sin(sqrt(pow(SatPos->O2S_EARPOS.E, 2) + 6.25*D2R*D2R ));
	delta_w = Kw / sin(sqrt(pow(SatPos->O2S_EARPOS.E, 2) + 2.25*D2R*D2R));
	
	*delta_Trop = delta_w + delta_d;
	return true;

}