WRF进阶：撰写Vtable文件添加气象场数据源：以JRA55积雪数据为例

最新推荐文章于 2024-09-04 12:41:25 发布

什么都不会的张同学

最新推荐文章于 2024-09-04 12:41:25 发布

阅读量1.4k

点赞数 1

分类专栏：模式学习文章标签： linux bash WRF

本文链接：https://blog.csdn.net/weixin_43750300/article/details/130069280

版权

模式学习专栏收录该内容

20 篇文章 74 订阅

订阅专栏

在WRF的WPS前处理中，ungirb.exe可以将下载的气象场grib文件解码并转为中间格式，最终通过metgrid.exe与包含域信息的geo_em文件写入met_em文件中。
将气象数据转为metgrid.exe可以读取的中间格式有两种方法，一种为使用NCL或者Fortran直接写出：该方法我在博客NCL撰写WPS intermediate file添加海冰场中有所介绍，这里介绍另一种方法：撰写Vtable文件，再使用ungrib.exe对其解码，得到中间文件。

Vtable文件

ungrib,exe是通过链接Vtable文件，知道需要解码的GRIB文件格式、变量信息，对其解码，也就是对应我们的ln -sf ungrib/Variable_Tables/Vtable.** Vtable 这段代码，不同的数据来源需要链接不同的Vtable才能正常解码。
虽然WPS本身提供了众多数据集的Vtables，但是ungrib不可能以GRIB格式预测每种可能的气象数据来源。当ungrib.exe处理新的数据源时，用户可以从头开始创建新的Vtable，也可以使用现有的Vtable作为示例。
关于Vtable的价绍官方文档以及翻译：创建和编辑Vtable文件有很好的价绍，简单而言，用户创建编辑Vtable需要了解一下信息：

GRIB code：气象场的GRIB码（GRIB code，在该数据集集中是一个独一无二的数字，类似于ID。
Level信息
数据本身在metgrid中的描述和名称

想要了解这些，可以借助ungrib下的g1print.exe和g2print.exe工具并结合官方数据来源给出的信息。

数据下载

JRA55的数据可以去JMA下，也可以去NACR下，不过JMA实在是太慢了，而NACR会主动提供批量的py和csh脚本下载。
不知为何我的电脑在Windows下用python返回的数据请求都是404，于是放在服务器上用wget下了，速度也很快的，这里下的是6小时地面再分析中的雪水当量数据，官方描述为：
在这里插入图片描述
其实JRA55还有一个雪深再分析，但是每天只有18：00，这就不太适合用ungrib了，不过WRF里雪深雪水当量原本就是能互相换算的，自己把变量名对应好就行。
下面是NACR提供的脚本，把pwd和opt换成你的邮箱密码就行：

#!/usr/bin/env csh

set pswd = $1
if(x$pswd == x && `env | grep RDAPSWD` != '') then
  set pswd = $RDAPSWD
endif
if(x$pswd == x) then
  echo Usage: $0 YourPassword
  exit 1
endif
set v = `wget -V |grep 'GNU Wget ' | cut -d ' ' -f 3`
set a = `echo $v | cut -d '.' -f 1`
set b = `echo $v | cut -d '.' -f 2`
if(100 * $a + $b > 109) then
  set opt = 'wget --no-check-certificate'
else
  set opt = 'wget'
endif

set opt1 = '-O Authentication.log --save-cookies auth.rda_ucar_edu --post-data'
set opt2 = "email=youremail&passwd=$pswd&action=login"
$opt $opt1="$opt2" https://rda.ucar.edu/cgi-bin/login
if ( $status == 6 ) then
  echo 'Please check that your password is correct.'
  echo "Usage: $0 YourPassword"
  exit 1
endif
set opt1 = "-N --load-cookies auth.rda_ucar_edu"
set opt2 = "$opt $opt1 https://rda.ucar.edu"
set filelist= (\ /data/ds628.0/anl_land125/2022/anl_land125.225_soilw.2022120100_2022123118  \
)
while($#filelist > 0)
  set syscmd = "$opt2$filelist[1]"
  echo "$syscmd ..."
  $syscmd
  shift filelist
end

中间的filelist你可以自己构建，也可以去NACR选的时候自带提供……无所谓。

查看GRIB信息

由于JRA55都是GRIB1格式，直接用g1print,exe打印以下就是了：./g1print.exe /JRA55/land/anl_land125.065_snwe.2020090100_2020093018
返回：

Fortran Unit = 0
UNIX File descriptor: 3

----------------------------------------------------
 rec GRIB GRIB  Lvl  Lvl  Lvl         Time      Fcst
 Num Code name  Code one  two                   hour
----------------------------------------------------
   1  65            1    0    0  2020-09-01_00:00 + 00
   2  65            1    0    0  2020-09-01_06:00 + 00
   3  65            1    0    0  2020-09-01_12:00 + 00
   4  65            1    0    0  2020-09-01_18:00 + 00
   5  65            1    0    0  2020-09-02_00:00 + 00
   6  65            1    0    0  2020-09-02_06:00 + 00
   7  65            1    0    0  2020-09-02_12:00 + 00
   8  65            1    0    0  2020-09-02_18:00 + 00
   9  65            1    0    0  2020-09-03_00:00 + 00
  10  65            1    0    0  2020-09-03_06:00 + 00
  11  65            1    0    0  2020-09-03_12:00 + 00
  12  65            1    0    0  2020-09-03_18:00 + 00
  13  65            1    0    0  2020-09-04_00:00 + 00
  14  65            1    0    0  2020-09-04_06:00 + 00
  15  65            1    0    0  2020-09-04_12:00 + 00
  16  65            1    0    0  2020-09-04_18:00 + 00
  17  65            1    0    0  2020-09-05_00:00 + 00
  18  65            1    0    0  2020-09-05_06:00 + 00
  19  65            1    0    0  2020-09-05_12:00 + 00
  20  65            1    0    0  2020-09-05_18:00 + 00
  21  65            1    0    0  2020-09-06_00:00 + 00
  22  65            1    0    0  2020-09-06_06:00 + 00
  23  65            1    0    0  2020-09-06_12:00 + 00
  24  65            1    0    0  2020-09-06_18:00 + 00
  25  65            1    0    0  2020-09-07_00:00 + 00
  26  65            1    0    0  2020-09-07_06:00 + 00
  27  65            1    0    0  2020-09-07_12:00 + 00
  28  65            1    0    0  2020-09-07_18:00 + 00
  29  65            1    0    0  2020-09-08_00:00 + 00
  30  65            1    0    0  2020-09-08_06:00 + 00
  31  65            1    0    0  2020-09-08_12:00 + 00
  32  65            1    0    0  2020-09-08_18:00 + 00
  33  65            1    0    0  2020-09-09_00:00 + 00
  34  65            1    0    0  2020-09-09_06:00 + 00
  35  65            1    0    0  2020-09-09_12:00 + 00
  36  65            1    0    0  2020-09-09_18:00 + 00
  37  65            1    0    0  2020-09-10_00:00 + 00
  38  65            1    0    0  2020-09-10_06:00 + 00
  39  65            1    0    0  2020-09-10_12:00 + 00
  40  65            1    0    0  2020-09-10_18:00 + 00
  41  65            1    0    0  2020-09-11_00:00 + 00
  42  65            1    0    0  2020-09-11_06:00 + 00
  43  65            1    0    0  2020-09-11_12:00 + 00
  44  65            1    0    0  2020-09-11_18:00 + 00
  45  65            1    0    0  2020-09-12_00:00 + 00
  46  65            1    0    0  2020-09-12_06:00 + 00
  47  65            1    0    0  2020-09-12_12:00 + 00
  48  65            1    0    0  2020-09-12_18:00 + 00
  49  65            1    0    0  2020-09-13_00:00 + 00
  50  65            1    0    0  2020-09-13_06:00 + 00
  51  65            1    0    0  2020-09-13_12:00 + 00
  52  65            1    0    0  2020-09-13_18:00 + 00
  53  65            1    0    0  2020-09-14_00:00 + 00
  54  65            1    0    0  2020-09-14_06:00 + 00
  55  65            1    0    0  2020-09-14_12:00 + 00
  56  65            1    0    0  2020-09-14_18:00 + 00
  57  65            1    0    0  2020-09-15_00:00 + 00
  58  65            1    0    0  2020-09-15_06:00 + 00
  59  65            1    0    0  2020-09-15_12:00 + 00
  60  65            1    0    0  2020-09-15_18:00 + 00
  61  65            1    0    0  2020-09-16_00:00 + 00
  62  65            1    0    0  2020-09-16_06:00 + 00
  63  65            1    0    0  2020-09-16_12:00 + 00
  64  65            1    0    0  2020-09-16_18:00 + 00
  65  65            1    0    0  2020-09-17_00:00 + 00
  66  65            1    0    0  2020-09-17_06:00 + 00
  67  65            1    0    0  2020-09-17_12:00 + 00
  68  65            1    0    0  2020-09-17_18:00 + 00
  69  65            1    0    0  2020-09-18_00:00 + 00
  70  65            1    0    0  2020-09-18_06:00 + 00
  71  65            1    0    0  2020-09-18_12:00 + 00
  72  65            1    0    0  2020-09-18_18:00 + 00
  73  65            1    0    0  2020-09-19_00:00 + 00
  74  65            1    0    0  2020-09-19_06:00 + 00
  75  65            1    0    0  2020-09-19_12:00 + 00
  76  65            1    0    0  2020-09-19_18:00 + 00
  77  65            1    0    0  2020-09-20_00:00 + 00
  78  65            1    0    0  2020-09-20_06:00 + 00
  79  65            1    0    0  2020-09-20_12:00 + 00
  80  65            1    0    0  2020-09-20_18:00 + 00
  81  65            1    0    0  2020-09-21_00:00 + 00
  82  65            1    0    0  2020-09-21_06:00 + 00
  83  65            1    0    0  2020-09-21_12:00 + 00
  84  65            1    0    0  2020-09-21_18:00 + 00
  85  65            1    0    0  2020-09-22_00:00 + 00
  86  65            1    0    0  2020-09-22_06:00 + 00
  87  65            1    0    0  2020-09-22_12:00 + 00
  88  65            1    0    0  2020-09-22_18:00 + 00
  89  65            1    0    0  2020-09-23_00:00 + 00
  90  65            1    0    0  2020-09-23_06:00 + 00
  91  65            1    0    0  2020-09-23_12:00 + 00
  92  65            1    0    0  2020-09-23_18:00 + 00
  93  65            1    0    0  2020-09-24_00:00 + 00
  94  65            1    0    0  2020-09-24_06:00 + 00
  95  65            1    0    0  2020-09-24_12:00 + 00
  96  65            1    0    0  2020-09-24_18:00 + 00
  97  65            1    0    0  2020-09-25_00:00 + 00
  98  65            1    0    0  2020-09-25_06:00 + 00
  99  65            1    0    0  2020-09-25_12:00 + 00
 100  65            1    0    0  2020-09-25_18:00 + 00
 101  65            1    0    0  2020-09-26_00:00 + 00
 102  65            1    0    0  2020-09-26_06:00 + 00
 103  65            1    0    0  2020-09-26_12:00 + 00
 104  65            1    0    0  2020-09-26_18:00 + 00
 105  65            1    0    0  2020-09-27_00:00 + 00
 106  65            1    0    0  2020-09-27_06:00 + 00
 107  65            1    0    0  2020-09-27_12:00 + 00
 108  65            1    0    0  2020-09-27_18:00 + 00
 109  65            1    0    0  2020-09-28_00:00 + 00
 110  65            1    0    0  2020-09-28_06:00 + 00
 111  65            1    0    0  2020-09-28_12:00 + 00
 112  65            1    0    0  2020-09-28_18:00 + 00
 113  65            1    0    0  2020-09-29_00:00 + 00
 114  65            1    0    0  2020-09-29_06:00 + 00
 115  65            1    0    0  2020-09-29_12:00 + 00
 116  65            1    0    0  2020-09-29_18:00 + 00
 117  65            1    0    0  2020-09-30_00:00 + 00
 118  65            1    0    0  2020-09-30_06:00 + 00
 119  65            1    0    0  2020-09-30_12:00 + 00
 120  65            1    0    0  2020-09-30_18:00 + 00

***** End-Of-File on C unit   3

这样GRIB code和LEV code都有啦，开始写Vtable吧！

创建VTABLE

直接写一个肯定很麻烦，我这里是从ECMWF直接cp再改的，这无所谓：

GRIB | Level| Level| Level| metgrid  |  metgrid | metgrid                                  |
Code | Code |   1  |   2  | Name     |  Units   | Description                              |
-----+------+------+------+----------+----------+------------------------------------------+
  65 |  1   |   0  |      | SNOW     | kg m-2   |Water Equivalent of Accumulated Snow Depth|
  66 |  1   |   0  |      | SNOWH    | m        | Physical Snow Depth                      
-----+------+------+------+----------+----------+------------------------------------------+
#

格式检查

此时链接我们创建的Vtable下载的积雪数据，便能成功解码，你可以通过util下的d_intermediate.exe工具检验一下你的中间格式是什么样子。

================================================
FIELD = SNOW
UNITS = kg m-2 DESCRIPTION = Water Equivalent of Accumulated Snow Depth
DATE = 2020-01-26_00:00:00 FCST = 0.000000
SOURCE = unknown model and orig center
LEVEL = 200100.000000
I,J DIMS = 288, 145
IPROJ = 0  PROJECTION = LAT LON
  REF_X, REF_Y = 1.000000, 1.000000
  REF_LAT, REF_LON = 90.000008, 0.000000
  DLAT, DLON = -1.250000, 1.250000
  EARTH_RADIUS = 6367.470215
DATA(1,1)=-1000000015047466219876688855040.000000

================================================
FIELD = SNOWH
UNITS = m DESCRIPTION = Physical Snow Depth
DATE = 2020-01-26_00:00:00 FCST = 0.000000
SOURCE = unknown model and orig center
LEVEL = 200100.000000
I,J DIMS = 288, 145
IPROJ = 0  PROJECTION = LAT LON
  REF_X, REF_Y = 1.000000, 1.000000
  REF_LAT, REF_LON = 90.000008, 0.000000
  DLAT, DLON = -1.250000, 1.250000
  EARTH_RADIUS = 6367.470215
DATA(1,1)=-4999999721059844884160053248.000000

修改METGRID.TBL

METGRID.TBL文件是一个文本文件，它定义了要用网格插值的每个气象场的参数。换句话说，只有参数在TBL文件中被定义，它才能被metgrid.exe文件读取插值。
我们应当注意，尽管在METGRID.TBL存在着SNOW这个变量，但它的一些属性并不一定与我们的SNOW对应，此时我们要修改使得其被正确读取。
比如这里我们就需要加入missing_value=-1.E20这一行。
在这里插入图片描述
想要知道更详细的信息，可以用NCL的filedump看一看：ncl_filedump JRA55/land/anl_land125.065_snwe.2020010100_2020013118.grib

file global attributes:
   dimensions:
      initial_time0_hours = 124
      g0_lat_1 = 145
      g0_lon_2 = 288
   variables:
      float SnWe_GDS0_SFC ( initial_time0_hours, g0_lat_1, g0_lon_2 )
         sub_center :   241
         center :       Japanese Meteorological Agency - Tokyo (RSMC)
         long_name :    Water equivalent of accumulated snow depth
         units :        kg/m^2
         _FillValue :   1e+20
         level_indicator :      1
         gds_grid_type :        0
         parameter_table_version :      200
         parameter_number :     65
         forecast_time :        0
         forecast_time_units :  hours