? ?做自己的事WRF项目模拟时,需要用中国土地利用数据替换WRF嵌入旧的土地利用数据,需要更换土地利用。我使用的中国土地利用数据来自中国土地利用数据(1980-2015)。我下载的lucc以1980为例演示。
1、坐标转换
? ?(1)下载的lucc如下图所示,1980数据: ? ?(2)在ArcMap打开属性表可以看到,lucc1980年的坐标系是投影坐标系Krasovsky_1940_Albers,数据格式为GRID。这里需要将投影坐标系转换为地理坐标系WGS_84,并且将GRID转换为tiff格式。
2、重分类
? ?因为根据WRF对lucc目前只能输入要求USGS-24分类、IGBP-20分类或IGBP-21分类系统,我使用的中国土地利用数据有6个一级分类26个二级分类(详见lucc1980中的lucc分类系统)和USGS与IGBP分类系统不同。所以这里我需要将就lucc将分类系统转换为USGS-在这里,我参考了我师姐的转换标准,如下图所示:
3、tiff格式转二进制
? ?WRF所有静态地理数据都需要二进制格式,我们分类得很好lucc1980为tiff格式可参考用户手册转二进制方法和气象家园高精度地形数据的应用WRF中(ASTER,SRTM)中的方法。 ? ?将分类好的lucc1980拷贝到Ubuntu中的WRF模型下的WPS_GEOG文件夹下,如下图所示。Ubuntu将文件夹更改为命令行lucc1980。 其中的命令gdal_translate -of ENVI xxx1.tif xxx2.bil 中 xxx我们需要转移的文件,xxx22为我们转换的文件 在运行这步前需要安装GDAL(Geospatial Data Asbstraction Library),然而,试错后,气象家园的博主认为命令应该改为gdal_translate -of ENVI -co INTERLEAVE=BSQ xxx1.tif xxx2.bil,证明博主的方法是正确的。如果你的话Ubuntu没有安装GDAL如果需要安装,安装gdal)
test@test:~/Build_WRF/WPS_GEOG/lucc1980$ gdal_translate -of ENVI -co INTERLEAVE=BSQ lucc1980_ProjectRaster_reclassification.tif.tif data.bil ERROR 4: `lucc1980_ProjectRaster_reclassification.tif.tif' does not exist in the file system, and is not recognised as a supported dataset name. GDALOpen failed - 4 `lucc1980_ProjectRaster_reclassification.tif.tif' does not exist in the file system, and is not recognised as a supported dataset name. test@test:~/Build_WRF/WPS_GEOG/lucc1980$ gdal_translate -of ENVI -co INTERLEAVE=BSQ lucc1980_ProjectRaster_reclassification.tif data.bil Input file size is 6157, 3393 0Warning 6: Driver ENVI does not support PIXELTYPE creation option ...10...20...30...40...50...60...70...80...90...100 - done.
获得以下三份文件: 打开data.hdr可以看到以下基本信息
ENVI
description = {
data.bil}
samples = 6157
lines = 3393
bands = 1
header offset = 0
file type = ENVI Standard
data type = 1
interleave = bsq
byte order = 0
map info = {
Geographic Lat/Lon, 1, 1, 66.2922997247769, 54.9785484716575, 0.0117836266203959, 0.0117836266203959,WGS-84}
coordinate system string = {
GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]]}
band names = {
Band 1}
在lucc1980中新建一个空白文档并命名为index,这里可以参考官网文档Writing Static Data to the Geogrid Binary Format和博客WRF更换静态下垫面数据。根据data.hdr头文件和数据内容,index文件设置如下:
type=categorical
category_min=1
category_max=24
projection=regular_ll
dx=0.0117836266203959
dy=0.0117836266203959
known_x=1.0
known_y=1.0
known_lat=14.9967033487
known_lon=66.2922997248
wordsize=1
tile_x=6157
tile_y=3393
tile_z=1
units="category"
description="USGS 24-category land use categories"
mminlu="USGS"
missing_value=128
iswater=16
islake=-1
isice=24
isurban=1
row_order=top_bottom
这里需要注意的是:known_lat和known_lon是tiff文件左下角的经纬度,这需要在ArcMap中打开数据的属性,可以找到。missing_value是因为我们的tiff是一张矩形的影像,里面包括了中国土地利用,但是除此之外还有一些不是分类体系里面的像元,这些就是缺省值,现在你可能不知道,当你运行geogrid.exe时,会出现,因为我已经执行过一遍了,故知道缺省值,你到时只要遇到加上去就可以了。row_order这一项是因为读取需要从顶部往底部读取,如果不设置那么默认值为bottom_top,则读出来的影像是倒过来的。其他选项如何设置,自行查看官方文档。 接下来,将data.bil改为00001-06157.00001-03393,其中6157为lucc1980的列数,3393为lucc1980的行数,然后把其余文件删除即可,这样tiff转换二进制就成功了,如下图所示。
4、修改GEOGRID.TBL
在WPS文件夹中的geogrid文件夹中,找到GEOGRID.TBL文件,打开后找到name=LANDUSEF项,分别在每一项的后面面加上 landmask_water = lucc1980:16 # Calculate a landmask from this field、interp_option = lucc1980:nearest_neighbor、rel_path = lucc1980:lucc1980/,然后保存即可。
5、修改namelist.wps
其中最重要的是需要将 geog_data_res设置为 ‘usgs_30s+default’,‘lucc1980+default’,特别提醒,粗域和嵌套域的土地利用类型的分类数都必须是一样的,比如这里usgs_30s是24类,lucc1980也必须是24类,如果设置不同的话在执行real.exe会报错。我的namelist.wps设置如下所示:
&share
wrf_core = 'ARW',
max_dom = 2,
start_date = '2021-01-01_00:00:00','2021-01-01_00:00:00',
end_date = '2021-01-03_18:00:00','2021-01-01_00:00:00',
interval_seconds = 21600
io_form_geogrid = 2,
debug_level=0
/
&geogrid
parent_id = 1,1,
parent_grid_ratio = 1,3,
i_parent_start = 1,30,
j_parent_start = 1,30,
e_we = 100,121,
e_sn = 100,121,
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!! IMPORTANT NOTE !!!!!!!!!!!!!!!!!!!!!!!!!!!!
! The default datasets used to produce the MAXSNOALB and ALBEDO12M
! fields have changed in WPS v4.0. These fields are now interpolated
! from MODIS-based datasets.
!
! To match the output given by the default namelist.wps in WPS v3.9.1,
! the following setting for geog_data_res may be used:
!
! geog_data_res = 'maxsnowalb_ncep+albedo_ncep+default', 'maxsnowalb_ncep+albe$
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!! IMPORTANT NOTE !!!!!!!!!!!!!!!!!!!!!!!!!!!!
geog_data_res = 'usgs_30s+default','lucc1980+default',
dx = 27000,
dy = 27000,
map_proj = 'mercator',
ref_lat = 43.935972,
ref_lon = 87.528372,
truelat1 = 30.0,
truelat2 = 60.0,
stand_lon = 87,
geog_data_path = '/home/test/Build_WRF/WPS_GEOG/'
/
&ungrib
out_format = 'WPS',
prefix = 'FILE',
/
&metgrid
fg_name = 'FILE'
io_form_metgrid = 2,
/
5、运行程序
我这个案例是模拟的新疆区域,先查看模拟域位置对不对。
test@test:~/Build_WRF/WPS-4.1$ ncl util/plotgrids_new.ncl
然后运行geogrid.exe
test@test:~/Build_WRF/WPS-4.1$ ./geogrid.exe
运行完后,看到Successful completion of geogrid,说明运行成功了,并且在WPS中生成以下两个文件。 整个流程到这里就结束了,还请给批评指正!希望对大家有帮助!