dtcenter / MET

Model Evaluation Tools
https://dtcenter.org/community-code/model-evaluation-tools-met
Apache License 2.0
78 stars 24 forks source link

Enhance MET's grid code to support rotated projections (rotated lat/lon in particular) #1029

Closed dwfncar closed 6 years ago

dwfncar commented 6 years ago

This request was made by Bob Craig at the Air Force. He sent us a sample GALWEM data file over Afghanistan that's on a rotate lat/lon projection:


rotated lat-lon grid:(840 x 1020) units 1e-06 input WE:SN output WE:SN res 48
    lat -6.884500 to 6.871999 by 0.013500
    lon 354.330000 to 5.656499 by 0.013500 #points=856800
    south pole lat=-58.250000 lon=66.599999 angle of rot=0.000000


This task is to enhance the MET library code to handle this projection in particular and rotated projections more generally. This sample data file can be found on dakota in:
   /d3/projects/MET/MET_development/jira/MET-1029


Here's the RT MET-Help ticket for this request:
   https://rt.rap.ucar.edu/rt/Ticket/Display.html?id=86320


Some issues to consider:


Need to identify a job number to which this development work could be charged. Currently assigned for met-8.0... but need to pick a realistic target release.

[MET-1029] created by johnhg

dwfncar commented 6 years ago

On 9/7/18, I tested out this code and there is more work to be done...


I've attached a single GRIB2 record (rotated_ll.grib2) for testing. This is the first record in the Afghanistan data file that Bob Craig sent to us:
PS.557WW_SC.U_DI.C_GP.GALWEM-RD_GR.C1P5KM_AR.AFGHANISTAN_DD.20180101_CY.00_FH.012_DF.GR2


(1) I ran IDV version 5.5 and made the attached plot: absv_rotated_ll.png
   This shows the pattern of the data and location over Pakistan and Afghanistan.


(2) I ran the development version of MET's plot_data_plane tool to make the attached plot: absv_rotated_ll_MET.png
   > trunk/met/bin/plot_data_plane rotated_ll.grib2 absv_rotated_ll_MET.ps 'name="ABSV"; level="P1013";'
   > convert -rotate 90 -background white -flatten absv_rotated_ll_MET.ps absv_rotated_ll_MET.png


(3) MET is reading the data in the same order as IDV. The pattern looks correct, but the location on the earth is wrong. MET is plotting this over the horn of Africa. To demonstrate, I ran the following regrid_data_plane and plot_data_plane commands to regrid to a lat/lon grid over the west edge of Africa. The resulting image (absv_AFRICA_MET.png) is attached.
   > trunk/met/bin/regrid_data_plane rotated_ll.grib2 'latlon 300 300 -15 -15 0.1 0.1' absv_AFRICA_MET.nc -field 'name="ABSV"; level="P1013";'
   > trunk/met/bin/plot_data_plane absv_AFRICA_MET.nc absv_AFRICA_MET.ps 'name="ABSV_P1013"; level="(,)";'
   > convert -rotate 90 -background white -flatten absv_AFRICA_MET.ps absv_AFRICA_MET.png


So the MET tools are reading the rotated lat/lon data from GRIB2 files correct, but not placing them on the correct spot on earth.

by johnhg

dwfncar commented 6 years ago

Randy made great progress and fixed MET so that it can parse the rotated lat/lon info from a GRIB2 file. John did some testing and found that we also need to enhance MET to read the rotated lat/lon info from MET's NetCDF files:


bin/pcp_combine -add rotated_ll.grib2 'name="ABSV"; level="P1013";' rotated_ll.nc
bin/plot_data_plane rotated_ll.nc rotated_ll.ps 'name="ABSV_P1013"; level="(,)";'


ERROR : read_netcdf_grid_v3() -> Projection "Rotated LatLon" not a currently supported type ("LatLon", "Mercator", "Lambert Conformal", "Polar Stereographic"). by johnhg

dwfncar commented 6 years ago

Randy added support for MET to read rotated lat/lon from its own NetCDF output files.


More testing... run wgrib2 to create some more rotated lat/lon data. Use this beta version of wgrib2 on theia:
   /scratch4/NCEPDEV/cpc/noscrub/Wesley.Ebisuzaki/grib2/wgrib2/wgrib2


Or download it and compile from here:
   http://ftp.cpc.ncep.noaa.gov/wd51we/tmp/wgrib2.tgz


Follow new_grid instructions from here:
   http://www.cpc.ncep.noaa.gov/products/wesley/wgrib2/new_grid.html
-new_grid rot-ll:sp_lon:sp_lat:sp_rot lon0:nlon:dlon lat0:nlat:dlat outfile


And run test to put GFS data onto existing Afghanistan domain:
wgrib2 gfs_TMP_Z2.grib2 -new_grid rot-ll:66.599999:-58.25:0.0 -5.67:1020:0.0135 -6.8845:840:0.0135 gfs_TMP_Z2_RotLL.grib2


Confirm that the grid matches:
wgrib2 -V gfs_TMP_Z2_RotLL.grib2
1:0:vt=2012040912:2 m above ground:12 hour fcst:TMP Temperature [K]:

    ndata=856800:undef=0:mean=300.916:min=266.084:max=315.884

    grid_template=1:winds(N/S):

        rotated lat-lon grid:(1020 x 840) units 1e-06 input WE:SN output WE:SN res 48

        lat -6.884500 to 4.442000 by 0.013500

        lon 354.330000 to 8.086500 by 0.013500 #points=856800

        south pole lat=-58.250000 lon=66.599999 angle of rot=0.000000


Looks good! Try defining a rotation of 30 degrees:
wgrib2 gfs_TMP_Z2.grib2 -new_grid rot-ll:66.599999:-58.25:30.0 -5.67:1020:0.0135 -6.8845:840:0.0135 gfs_TMP_Z2_RotLL_30deg.grib2


And plot it:
plot_data_plane gfs_TMP_Z2_RotLL_30deg.grib2 gfs_TMP_Z2_RotLL_30deg_MET.ps 'name="TMP"; level="Z2";'


Tried plotting in IDV but IDV errors out when trying to load the data.


Next tried ncl_convert2nc but that doesn't support a non-zero rotation:
ncl_convert2nc gfs_TMP_Z2_RotLL_30deg.grib2
warning:GdsRLLGrid: Nonzero rotation angle not yet supported for Grid template 1; no coordinates will be returned


Next I tried Panoply 4.9.4 with a similar result. It can plot the data with a 0 degree rotation but not the 30 degree rotation. by johnhg

dwfncar commented 6 years ago

Randy finished this for the met-8.0 for GRIB2, MET NetCDF, and the python interface.


May still need to be added for GRIB1 and CF-Compliant NetCDF. by johnhg