发明内容

为解决上述技术问题，本发明提供了一种联合GNSS与测风激光雷达的对流层温湿廓线反演方法，该方法在目前成熟的气象要素传感器、GNSS接收机和测风激光雷达硬件设备基础上进行算法层面的功能扩展，实现不依赖于历史探空数据的单站对流层大气温湿廓线反演获取，设备安装维护成本较低。因该方法不依赖历史探空数据，安装维护成本低，不仅适用于陆基观测平台，也适用于海基观测平台。

为达到上述目的，本发明的技术方案如下：

一种联合GNSS与测风激光雷达的对流层温湿廓线反演方法，包括如下步骤：

步骤一，获取基础观测量，包括：利用气象要素传感器测量GNSS天线位置处气象参数，包括气温、气压和比湿，利用GNSS接收机进行实时精密单点定位，并估计天顶总延迟，利用测风激光雷达测量分层高度上的垂直风速；

步骤二，提取中间变量，包括：根据GNSS天线位置处的气温和气压，利用压高方程计算气压廓线；根据天顶总延迟，利用映射函数计算斜路径总延迟；利用斜路径总延迟的理论值和观测值，结合数值仿真和参数优化方法反演湿折射率廓线；利用分层高度上的垂直风速计算Brunt-Väisälä频率，并利用Brunt-Väisälä频率反演虚位温廓线；

步骤三，求解最终输出变量，利用步骤二获得的气压廓线、湿折射率廓线和虚位温廓线，构造以各层高度处气温和比湿为未知量、各层高度处湿折射率和虚位温为因变量的二元非线性方程组，求解所述二元非线性方程组，剔除复数解和非物理含义解，得到各层高度处气温和比湿，形成对流层温湿廓线。

上述方案中，步骤一中，GNSS天线位置处气象参数包括气温T ₀ 、气压P ₀ 和比湿q ₀ ；实时精密单点定位包括GNSS天线位置处的经度lon、纬度lat和高程h ₀ ，由GNSS定位过程获得；天顶总延迟ZTD，由GNSS定位过程中通过参数估计获得；分层高度上的垂直风速表示为w _i (t)，w _i 表示第i层高度处的垂直风速，为一个时间序列，t代表时间。

上述方案中，步骤二中，计算气压廓线的具体方法如下：

对流层大气高度h _i 处的气压P _i 由下式进行计算：

（1）

其中，T ₀ 、P ₀ 分别为GNSS天线位置处的气温和气压。

上述方案中，步骤二中，计算斜路径总延迟的具体方法如下：

（1）由Saastamoninen模型计算天顶干延迟ZHD：

（2）

其中，lat为GNSS天线位置处的纬度，P ₀ 为GNSS天线位置处的气压；

（2）提取天顶湿延迟ZWD：

（3）

其中，ZTD为天顶总延迟；

（3）计算斜路径干延迟观测值 和斜路径湿延迟观测值 ：

由映射函数将ZHD、ZWD分别映射到各卫星斜路径方向，得到：

（4）

其中， 为卫星的仰角，为干项映射函数， 为湿项映射函数；

（4）计算各颗卫星的斜路径总延迟观测值 ：

（5）。

上述方案中，步骤二中，反演湿折射率廓线的具体方法如下：

首先，构造候选的折射率廓线N _i ，利用构造的折射率廓线N _i 进行射线追踪或抛物方 程仿真，得到各颗卫星的斜路径总延迟的理论值，记为 。

然后，计算斜路径总延迟的理论值和观测值之间的均方根误差，选择最小的均方 根误差对应的折射率廓线为最优解，即

（6）

其中，为各颗卫星的斜路径总延迟的观测值， 为各颗卫星的斜路径总延迟的理论值和观测值之间的均方根误差函数。

最后，计算湿折射率廓线：

折射率廓线的最优解 由干项 和湿项 组成，其中，干项 由Hopfield模型进行近似：

（7）

其中，h _i 是第i层的高度，h ₀ 为GNSS天线位置处的高度，T ₀ 为GNSS天线位置处的气温；

湿项由下式获得：

（8）

不同层高度处的湿项 形成湿折射率廓线。

上述方案中，步骤二中，计算Brunt-Väisälä频率的具体方法如下：

对测风激光雷达测得高度h _i 处的垂直风速w _i (t)进行傅里叶变换，频谱最大幅值对应的频率即为高度h _i 处的Brunt-Väisälä频率N _BV,i ：

（9）

其中，为傅里叶变换。

上述方案中，步骤二中，反演虚位温廓线的具体方法如下：

（1）计算GNSS天线位置处的虚位温θ _v,0 ：

（10）

其中，P _k 为1000hPa，T ₀ 、P ₀ 和q ₀ 分别为GNSS天线位置处的气温、气压和比湿；

（2）反演第i层高度处的虚位温θ _v,i ：

由Brunt-Väisälä频率N _BV,i 与虚位温θ _v,i 之间的关系，有：

（11）

其中，g为重力加速度，h _i 是第i层的高度，h _i-1 是第i-1层的高度，h ₀ 为GNSS天线位置处的高度，θ _v,i-1 为第i-1层高度处的虚位温；

根据(11)式，由θ _v,0 依次计算得到不同层高度h _i 处的虚位温θ _v,i ，形成虚位温廓线。

上述方案中，步骤三中，构造二元非线性方程组的具体方法如下：

由于各层高度h _i 处的虚位温θ _v,i 、气温T _i 、比湿q _i 、气压P _i 和水汽压e _i 之间存在如下关系：

(12)

其中，P _k 为1000hPa，比湿q _i 单位为kg/kg；

由(12)式可知，虚位温θ _v,i 是气温T _i 、气压P _i 和比湿q _i 的函数，由于气压P _i 在步骤二 中已经获得，此时将虚位温θ _v,i 简化为气温T _i 和比湿q _i 的二元函数，记为 ：

(13)；

湿折射率 是气温T _i 、气压P _i 和比湿q _i 的函数：

(14)

由于气压P _i 在步骤二中已经获得，此时将湿折射率 简化为气温T _i 和比湿q _i 的二元函数，记为：

(15)；

联立式(13)和(15)，得到以气温T _i 和比湿q _i 为自变量，虚位温θ _v,i 和湿折射率 为因变量的二元非线性方程组：

(16)。

进一步的技术方案中，步骤三中，求解对流层温湿廓线的具体方法如下：

采用数值方法求解式（16）所列二元非线性方程组，对于多解情况，首先剔除复数解，其次依据气温T _i 垂直分布趋势和比湿q _i 分布规律剔除非物理含义解，得到最终输出的不同高度处的气温T _i 和比湿q _i ，形成对流层温湿廓线。

通过上述技术方案，本发明提供的一种联合GNSS与测风激光雷达的对流层温湿廓线反演方法具有如下有益效果：

本发明将GNSS对流层延迟和测风激光雷达分层高度上的垂直风速两种基础观测量相结合，在GNSS接收机和测风激光雷达原有的定位和测风功能基础上，扩展出对流层大气温湿廓线反演探测的新功能。考虑到目前的陆基和海基气象观测平台，通常安装有气象要素传感器、GNSS接收机和测风激光雷达中的一种或多种，本发明能够在气象观测平台已有仪器设备基础上，通过补充缺少的气象要素传感器、GNSS接收机或测风激光雷达，并增加数据处理单元和算法软件，较为便捷地实现单站对流层大气温湿廓线反演探测功能。该方式设备安装和维护的成本低，能够为大规模组网观测，尤其是海洋浮标网系统观测提供技术支撑。

本发明通过实时精密单点定位方式实现GNSS对流层延迟的实时获取，通过测风激光雷达实时测量分层高度上的垂直风速，并结合气象要素传感器实时测量的气象要素对上述两种基础观测量进行处理。在实际的应用场景中，根据对流层温湿廓线数据使用者对数据实时性的需求配置具备相应计算能力的数据处理单元，并根据数据使用者对精度的需求设置数据处理过程中参数优化算法的取值范围、约束条件、收敛标准等，实现对流层温湿廓线的实时或近实时输出。相对于传统探空探测方式，具有更高的时空分辨率和显著的灵活性。

本发明在对GNSS对流层延迟观测量进行数据处理的过程中采用参数优化方式获取对流层大气湿折射率廓线，结合对测风激光雷达分层高度上的垂直风速观测量进行数据处理并获取的虚位温廓线，构造二元非线性方程组并进行数值求解，得到最终输出的温湿廓线。该方式不依赖于历史探空数据,相对于其他建立在历史探空数据基础上的训练或回归建模方法，例如微波辐射计普遍采用的神经网络建模反演工作模式，该方法具有更高的普适性，例如对于缺乏历史探空资料的海基观测平台，则同样适用。