一种基于再分析数据的天文台址水汽评估方法

1.本发明涉及大气水汽评估技术领域，更具体的说是涉及一种基于再分析数据的天文台址水汽评估方法。


背景技术：

2.水汽是大气的重要组成成分，在天气预报、气候变化研究、辐射传输等大气研究中具有重要的地位，也是衡量一个射电天文台址重要的性能指标之一。在不同的信号接收频段，大气对射电信号的水汽吸收效应各异(参见图1)。大气中的水汽不仅引起信号的衰减，水汽的起伏所引起的大气噪声辐射的起伏同样会成为望远镜灵敏度的限制；水汽含量在空间分布上的小尺度不均匀性和变化还会引起到达天线口面的无线电波光程的变化，这对天线观测也会引起天线指向的变化和天线增益的降低。较低的大气含量可有效降低射电信号延迟，提高信号透射率。因此研究和探索台址水汽变化规律对射电天文观测具有重要的实践价值。
3.然而，水汽在大气中分布极不均匀，变化十分迅速，常规地面资料很难准确反映水汽输送和分布状况。因此，用新技术新手段精确测定大气中的水汽含量，是当今气象学乃至天文选址水汽探测过程中所面临的重要任务。目前，探测对流层大气水汽含量的方法及优缺点如下：
4.(1)无线电探空
5.无线电探空(radiosonde，sonde)是目前气象业务中探测大气水汽的最常用手段，其探测原理是通过施放气球，收集高空大气的温度、气压、湿度等观测数据，但是无线电探空施放气球和探空仪的成本很高，目前建设的探空站数量也有限，在海洋和极地，水汽资料缺乏，即使在陆地，高原和无人区测站也非常稀疏，且探空站分布不均匀，空间分辨率极差，且目前探空站的观测仅进行早晚两次，时间分辨率太低，不足以分辨水汽的时空变化，对中小尺度天气的监测特别不利。因此，无线电探空己经不能满足水汽时间空间多变性的要求。
6.(2)卫星观测
7.卫星观测是气象业务中水汽观测的一种辅助手段，卫星探测相较探空观测，时间和空间分辨率均有极大的优势，但反演的精度还有待进一步提高，卫星探测的红外技术，在晴空区才能反映水汽含量的垂直分布，得不到云下的水汽信息，而搭载的分光计和微波辐射计探测的垂直分辨率又很有限，反演精度还不高。
8.(3)微波辐射计(wvr)(主要手段)
9.微波辐射计(water vapor radiometer，wvr)是目前探测大气水汽最为精确的一种设备，但微波辐射计在有浓云时候穿透能力会下降，雨天更是会产生较大的误差，不能满足全天候观测的要求，而星载的微波辐射计受背景温度的影响较大，在陆地上的精度比海面上的差。同时，微波辐射计的价钱也比较昂贵，且需要标定，限制了它的普及应用。
10.(4)雷达探测
11.新一代多普勒天气雷达能探测垂直积分液态水含量，而raman激光雷达可探测大
气中水汽的分布，但是，雷达的探测成本高，且探测范围有限，难以应用于大范围的观测。
12.(5)飞机探测
13.由飞机携带水汽观测仪器，也能完成水汽的探侧，但是飞机探测的时间和范围都很有限，且只能按照规定航线，故只能用于特定任务、科研、人工增雨等的特殊观测，且成本也极高。
14.(6)地面湿度计观测
15.地面湿度计观测是常规气象业务中观测地面大气水汽的一种方法，其测量是2m近地层空气中的水汽状况，不能反映整层大气的水汽含量和变化，因此对预报降水和剧变天气的作用极其有限。
16.因此如何克服采用水汽辐射计或者太阳光度计测量大气水汽含量时，测量结果均受到天气影响的缺陷，是本领域技术人员亟需解决的问题。


技术实现要素：

17.有鉴于此，本发明提供了一种基于再分析数据的天文台址水汽评估方法，克服了上述缺陷。
18.为了实现上述目的，本发明提供如下技术方案：
19.一种基于再分析数据的天文台址水汽评估方法，具体步骤为：
20.步骤1、获取全球再分析数据；
21.步骤2、根据全球再分析数据提取待评估站点所在网格的四个格点的气象参数提取，并利用气象参数计算各个格点的大气水汽含量；
22.步骤3、将各个格点的大气水汽含量利用双线性差值的方法进行优化，获得待评估站点的大气水汽含量。
23.可选的，气象参数包括分层气压数据、分层比湿数据、分层高度数据、时间信息。
24.可选的，步骤2中待评估站点的区域大气水汽含量计算的步骤为：
25.步骤2-1、获取待评估站点的地理坐标；
26.步骤2-2、根据待评估站点的地理坐标确定其所在网格，网格由全球再分析数据形成；
27.步骤2-3、提取所在网格上各个格点的气象参数；
28.步骤2-4、根据各个格点的气象参数计算各个格点的大气水汽含量。
29.可选的，步骤2-4的具体步骤为：
30.根据各个格点的分层高度数据确定地面海拔高度；
31.根据各个格点的分层气压数据确定当前地面气压值；
32.根据各个格点的分层比湿数据确定当前地面比湿值；
33.依据地面海拔高度、地面气压值、地面比湿值获得各个格点的大气水汽含量；
34.各个格点的大气水汽含量组成待评估站点的区域大气水汽含量。
35.可选的，确定地面海拔高度过程包括分层海拔高度的计算，其计算公式为：
36.h＝z/g；
37.式中，z为重力势；g为重力加速度。
38.可选的，各个格点的大气水汽含量为：
[0039][0040]
式中，g为重力加速度；p为对流层分层气压；ps为当前地面气压值，q为当前地面比湿值。
[0041]
可选的，步骤3中待评估站点的大气水汽含量的具体步骤为：将各个格点的大气水汽含量计算值利用双线性差值法插回到待评估站点，获得待评估站点的大气水汽含量。
[0042]
可选的，双线性差值法包括在x方向和y方向上分别进行一次线型差值。
[0043]
可选的，x方向上的线性插值为：
[0044][0045][0046]
式中，q
11
、q
12
、q
21
、q
22
均为格点；x，y为待评估站点的坐标；x1、x2、y1、y1为格点的坐标；r1和r2均为在x方向线性插值后得到的坐标点。
[0047]
可选的，y方向上的线性插值为：
[0048][0049]
式中，y为待评估站点的纵坐标；y1、y2为格点的纵坐标；r1和r2均为在x方向线性插值后得到的坐标点。
[0050]
经由上述的技术方案可知，与现有技术相比，本发明公开了一种基于再分析数据的天文台址水汽评估方法，通过国内外气象组织获取气象数据，不依赖任何观测设备，因此不受测量设备和天气的影响，能够实时对天文台址水汽评估，提高了评估的准确性。
附图说明
[0051]
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据提供的附图获得其他的附图。
[0052]
图1为在不同的信号接收频段大气对射电信号的水汽吸收效果示意图；
[0053]
图2为本发明的方法流程图；
[0054]
图3为本发明中ecmwf数据分层示意图；
[0055]
图4为本发明中ecmwf数据网格示意图；
[0056]
图5为本发明中奇台观测站水汽精度比较结果示意图；
[0057]
图6为本发明中不同站点随月份水汽变化结果示意图；
[0058]
图7为本发明中不同季节水汽分析示意图。
具体实施方式
[0059]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他
实施例，都属于本发明保护的范围。
[0060]
实施例1
[0061]
本发明实施例公开了一种基于再分析数据的天文台址水汽评估方法，通过提取欧洲中期天气预报中心(简称ecmwf)的再分析资料中的有效气象参数，利用水汽与气象参数的相应关系计算所选站点的水汽参量，用来进行野外台址水汽评估，同时引用线性相关函数进一步提高所选站点水汽的精确及合理性，大大提高了野外站点水汽评估效率，具体步骤，如图2所示，具体为：
[0062]
步骤1、获取全球再分析数据；
[0063]
在需要数据时，登录ecmwf网站，在产品页面选择需要提取的气象参数下载全球的再分析数据(在此以2015年1月1日00时，海拔h0为120m，所选参考点p(x,y)为例)，其参数包括：产品类型：reanalysis；变量名称：geopotential(重力势z)、specific humidity(比湿q)、pressure level(分层气压p)、year(年)、month(月)、day(日)、time(小时)等；格式类型：netcdf格式；选好参数后，点击下载，完成原始数据下载。
[0064]
其中，ecmwf可以为用户提供全球空间分辨率为0.25
°×
0.25
°
和时间分辨率为1h的完备气象数据，如图3所示，它把对流层大气分成37层，每层包含气温、湿度、气压等气象要素，其数据格式分别为netcdf4和grib。
[0065]
步骤2、根据全球再分析数据提取待评估站点所在网格的四个格点的气象参数提取，并利用气象参数计算各个格点的大气水汽含量，具体步骤为：
[0066]
步骤2-1、获取待评估站点的地理坐标；
[0067]
步骤2-2、根据待评估站点的地理坐标确定其所在网格，网格由全球再分析数据形成；
[0068]
步骤2-3、提取所在网格上各个格点的气象参数；
[0069]
步骤2-4、根据各个格点的气象参数计算各个格点的大气水汽含量。
[0070]
其具体为：
[0071]
获得原始数据后依据所选区域的地理坐标(经度和纬度)和时间信息，提取相应的气象数据，具体参数如下：分层气压数据pi；分层比湿数据分层高度数据时间信息t。
[0072]
由于所下载的原始数据为全球的再分析气象据，按照0.25
°×
0.25
°
的格点分布(即再分析资料的气象信息分布在图4网格中黑色圆点上)，要想获取所选站点的气象参数数据，首先确定所选站点的经纬度坐标，然而在提取所选站点的气象参数时，由于其候选站点的坐标几乎无法与再分析数据中的格点位置重合，即：所选站点坐标处无法提取再分析气象参数。
[0073]
在实际操作中，任意站点坐标必然位于四个格点围成的方形网格内(此处假设所选站点为p(x,y)位于以q
11
、q
12
、q
21
、q
22
四个格点围成的四边形区域内)，首先对四个格点的气象参数进行提取，再进一步计算出各个格点的水汽含量。
[0074]
步骤2-4-1、获取分层高度信息(h)，确定地面海拔高度。
[0075]
首先，根据分层海拔高度h与重力势z之间的关系，获取分层海拔高度h，其计算公式为：
[0076]
h＝z/g；
[0077]
其中g为重力加速度，取值9.8m/s2。
[0078]
其次，比较分层海拔高度h和本地海拔高度h0数据差异，选择37层h高度值大于h0为有效海拔数据。
[0079]
其中，重力势z的提取方法即读取geopotential(重力势z)数据：利用ncl语言首先对z数据进行提取。geopotential(重力势)参数中包含了scale_factor和add_offset。
[0080]
因此z的提取方式如下：
[0081]
z＝z*scale_factor add_offset；
[0082]
其参数说明如下：
[0083][0084]
步骤2-4-2、获取分层气压，确定当前地面气压值。
[0085]
气压参数由于没有修正量，直接利用ncl语言获取37层气压分量。仿照有效海拔数据的获取方法，找出对应的有效气压数据。
[0086]
其参数说明如下：
[0087]
int32 level(level)
[0088]
units:millibars
[0089]
long_name:pressure_level
[0090]
步骤2-4-3、获取分层比湿q，确定当前地面比湿值
[0091]
q参数中亦包含了scale_factor和add_offset，因此同z提取方法一样，得到分层比湿数据：
[0092]
q＝q*scale_factor add_offset；
[0093]
仿照有效海拔数据的获取方法，找出对应的有效比湿q数据。
[0094]
其参数说明如下：
[0095][0096]
步骤2-4-4、各个格点的pwv计算：
[0097][0098]
其中，g代表重力加速度9.8m/s2；p代表对流层分层气压；ps代表地面气压，q代表分层比湿。
[0099]
步骤3、将区域大气水汽含量计算值利用双线性差值的方法进行优化，获得待评估站点的大气水汽含量。
[0100]
把再分析资料计算的各个格点的pwv值利用双线性插值的方法插回到相应的站点，从而提高再分析资料的实用性和可靠性。
[0101]
双线性插值，又称为双线性内插。在数学上，双线性插值是有两个变量的插值函数的线性插值扩展，其核心思想是在两个方向分别进行一次线性插值。想得到未知函数f在点p＝(x,y)的值，假设我们已知函数f在q
11
＝(x1,y1)、q
12
＝(x1,y2),q
21
＝(x2,y1)以及q
22
＝(x2,y2)四个点的值，首先在x方向进行线性插值，得到：
[0102][0103][0104]
然后在y方向进行线性插值，得到：
[0105][0106]
这样就可以得到候选站址的数值，且线性插值的大小和顺序无关。
[0107]
本实施例中还包括数据有效性验证方法，具体为：
[0108]
数据准确度是分析结果与真值的相符程度。准确度通常用误差来表示，误差越小，表示分析结果的准确度越高。选择利用绝对误差和相对误差来评估pwv结果的有效性。
[0109]
假定仪器测定值pwv为参考值t0，ecmwf解算值pwv为待测值t1；绝对误差(ea)：
[0110]
ea＝t
1-t0；
[0111]
相对误差(er)：
[0112]er
＝(ea/t1)
×
100％；
[0113]
本实施例公开了实际的验证效果：
[0114]
1.新疆天文台奇台观测站水汽有效性验证
[0115]
所选站点坐标：纬度(43.60
°
)，经度(89.68
°
)；
[0116]
时间：2017-2018；
[0117]
海拔：1765m；
[0118]
参考设备：太阳辐射计；
[0119]
经计算该站点水汽精度如下：
[0120]
如图5所示，月均值绝对误差约2mm，相对误差小于20％。
[0121]
2.新疆天文台喀什慕士塔格光学观测站水汽有效性验证
[0122]
所选站点坐标：纬度(38.33)，经度(74.90
°
)；
[0123]
时间：2018年；
[0124]
海拔：4526m；
[0125]
参考设备：水汽辐射计；
[0126]
经计算该站点水汽精度如下：
[0127]
平均绝对误差约1mm，相对误差小于19％(见表1)。
[0128]
表1水汽辐射计与ecwmf数据偏差(2018)
[0129]
[0130][0131]
实施例2
[0132]
新疆克拉玛依地区野外站点水汽评估
[0133]
站点5坐标：经度(85.077)，纬度(46.030)，海拔600m；
[0134]
站点6坐标：经度(84.553)，纬度(45.784)，海拔1050m；
[0135]
站点9坐标：经度(83.795)，纬度(45.106)，海拔500m；
[0136]
时间：2011年-2020年；
[0137]
由于每个站点的数据均不在全球再分析数据格点上，对三个站点的水汽分别进行双线性插值计算，得到各自pwv值。通过图6和图7可以看出不同站点pwv随时间变化关系：6号站点pwv最小，9号站点pwv最大。
[0138]
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。
[0139]
对所公开的实施例的上述说明，使本领域专业技术人员能够实现或使用本发明。对这些实施例的多种修改对本领域的专业技术人员来说将是显而易见的，本文中所定义的一般原理可以在不脱离本发明的精神或范围的情况下，在其它实施例中实现。因此，本发明将不会被限制于本文所示的这些实施例，而是要符合与本文所公开的原理和新颖特点相一致的最宽的范围。