Wrong conservative remapping weights
Added by Harald Rybka almost 4 years ago
Hello,
I have tried to perform a conservative remapping procedure of RADAR data (1km resolution) onto a curvilinear model grid (approx. 3km resolution) but it seems like the calculated weights are wrong. I think this is because of the (probably) unknown/unsupported grid type of the RADAR observations.
Here are the following steps I have performed (and some additional information):
cdo griddes CCLM-0275_grid_bnds.nc > CCLM-0275_griddes.txt # get grid description file of target grid cdo gencon,CCLM-0275_griddes.txt RADKLIM_nomiss.nc gencon_wgt.nc # calculate remapping weights cdo remap,CCLM-0275_griddes.txt,gencon_wgt.nc RADKLIM_nomiss.nc RADKLIM_remapcon.nc # perform conservative remapping
The result looks like this (comparing original - 1km res with 3km output res):
During the generation of remapping weights many errors occur including map weights above 1 and below 0 (see image).
cdo gencon: Map weight < 0! grid1idx=982738 grid2idx=140606 nlink=1752516 wts=-3.32649e+15 cdo gencon: Map weight > 1! grid1idx=982739 grid2idx=140192 nlink=1752517 wts=8.96177e+14 [...] cdo gencon: Error: sum of wts for map1 42616 0.199219 1 [...]
Here are some outputs of the verifygrid and ncdump -h of the input and output grids:
>cdo -verifygrid RADKLIM_nomiss.nc cdo verifygrid: Grid consists of 990000 cells (type: curvilinear), of which cdo verifygrid: 990000 cells have 4 vertices cdo verifygrid: 61996 cells are not unique cdo verifygrid: 990000 cells have their vertices arranged in a clockwise order cdo verifygrid: Processed 2 variables ( 0.94s 182MB ) >cdo -verifygrid CCLM-0275_grid_bnds.nc cdo verifygrid: Grid consists of 180525 cells (type: curvilinear), of which cdo verifygrid: 180525 cells have 4 vertices cdo verifygrid: Processed 1 variable ( 0.19s 45MB ) >ncdump -h RADKLIM_nomiss.nc netcdf RADKLIM_nomiss { dimensions: time = UNLIMITED ; // (1 currently) x = 900 ; y = 1100 ; vertices = 4 ; variables: double time(time) ; time:standard_name = "time" ; time:long_name = "Time" ; time:units = "hours since 2001-01-01 00:50:00.0" ; time:calendar = "standard" ; time:axis = "T" ; double lon(y, x) ; lon:standard_name = "longitude" ; lon:long_name = "Longitude" ; lon:units = "degrees_east" ; lon:_CoordinateAxisType = "Lon" ; lon:bounds = "lon_bnds" ; double lon_bnds(y, x, vertices) ; double lat(y, x) ; lat:standard_name = "latitude" ; lat:long_name = "Latitude" ; lat:units = "degrees_north" ; lat:_CoordinateAxisType = "Lat" ; lat:bounds = "lat_bnds" ; double lat_bnds(y, x, vertices) ; double x(x) ; x:standard_name = "projection_x_coordinate" ; x:long_name = "RADOLAN Grid x coordinate of projection" ; x:units = "m" ; x:axis = "X" ; double y(y) ; y:standard_name = "projection_y_coordinate" ; y:long_name = "RADOLAN Grid y coordinate of projection" ; y:units = "m" ; y:axis = "Y" ; int crs ; crs:long_name = "RADOLAN Grid" ; crs:grid_mapping_name = "polar_stereographic" ; crs:latitude_of_projection_origin = 90. ; crs:straight_vertical_longitude_from_pole = 10. ; crs:semi_major_axis = 6370040. ; crs:semi_minor_axis = 6370040. ; crs:earth_radius = 6370040. ; crs:standard_parallel = 60. ; crs:false_easting = 0. ; crs:false_northing = 0. ; crs:inverse_flattening = 0. ; crs:scale_factor_at_projection_origin = 0.9330127019 ; crs:units = "m" ; crs:proj4 = "+proj=stere +lat_0=90 +lat_ts=90 +lon_0=10 +k_0=0.933012701892 +x_0=0 +y_0=0 +a=6370040 +b=6370040 +units=km +no_defs" ; float RR(time, y, x) ; RR:standard_name = "rainfall_amount" ; RR:long_name = "hourly rainfall " ; RR:units = "kg m-2" ; RR:grid_mapping = "crs" ; RR:coordinates = "lat lon" ; RR:_FillValue = 999.f ; RR:missing_value = 999.f ; >ncdump -h CCLM-0275_grid_bnds netcdf CCLM-0275_grid_bnds { dimensions: rlon = 415 ; rlat = 435 ; nv = 4 ; variables: double rlon(rlon) ; rlon:long_name = "longitude in rotated pole grid" ; rlon:standard_name = "grid_longitude" ; rlon:units = "degrees" ; rlon:axis = "X" ; double rlat(rlat) ; rlat:long_name = "latitude in rotated pole grid" ; rlat:standard_name = "grid_latitude" ; rlat:units = "degrees" ; rlat:axis = "Y" ; double lon(rlat, rlon) ; lon:long_name = "longitude" ; lon:standard_name = "longitude" ; lon:units = "degrees_east" ; lon:bounds = "lon_vertices" ; double lat(rlat, rlon) ; lat:long_name = "latitude" ; lat:standard_name = "latitude" ; lat:units = "degrees_north" ; lat:bounds = "lat_vertices" ; double lon_vertices(rlat, rlon, nv) ; lon_vertices:units = "degrees_east" ; double lat_vertices(rlat, rlon, nv) ; lat_vertices:units = "degrees_north" ; char rotated_pole ; rotated_pole:long_name = "coordinates of the rotated North Pole" ; rotated_pole:grid_mapping_name = "rotated_latitude_longitude" ; rotated_pole:grid_north_pole_latitude = 39.25 ; rotated_pole:grid_north_pole_longitude = -162. ; byte dummy(rlat, rlon) ; dummy:long_name = "dummy_variable" ; dummy:grid_mapping = "rotated_pole" ; dummy:coordinates = "lat lon" ;
Does anyone have a suggestion on how to solve this issue?
-- Harald
RADKLIM_ori.png (76.5 KB) RADKLIM_ori.png | Original RADAR data | ||
RADKLIM_remapped.png (133 KB) RADKLIM_remapped.png | RADAR data remapped |
Replies (1)
RE: Wrong conservative remapping weights - Added by Harald Rybka almost 4 years ago
I found the error....
The cell bounds have been defined clockwise not counterclockwise!