Geostationary file support
Added by Miles Sowden over 2 years ago
Is geostationary projection fully supported in CDO (or NCO for that matter)?
At this point I’m trying to achieve a rescaling into ½ or double grid spacing as per https://stackoverflow.com/questions/55746419/change-grid-size-of-a-netcdf-file. However other operators such as "cdo chiftx,1" also fail. (sidenote I can edit the grid file to gridtype=generic to "fix" the shiftx problem)
I do not want to reproject into latitude/longitude (given as variables) as my files are huge (1.8 Gb every ten minutes). But wish to work with the files as X & Y indices for now.
For this example I cropped a larger file 5500 x 5500 with the command below and attached the file for reference
ncks -h -d x,2700,2800 -d y,4200,4300 AHI_SENSOR_2000.nc AHI.nc
The larger grid description is given in AHI.txt (attached). (obtained with cdo griddes AHI_SENSOR_2000.nc > AHI.txt)
Reading the CF conventions http://cfconventions.org/Data/cf-conventions/cf-conventions-1.9/cf-conventions.html#_geostationary_projection
I think I need to add
scanningMode = 64
grid_mapping = geostationary
grid_mapping_name = geostationary
longitude_of_projection_origin = 140.7
latitude_of_projection_origin = 0
false_easting = 0.
false_northing = 0.
perspective_point_height = 35785863.
sweep_angle_axis="x"
fixed_angle_axis="y"
I think the problem is that both CDO and NCO do not yet fully support the regular grid obtained from geostationary satellites.
Thanks
Replies (4)
RE: Geostationary file support - Added by Karin Meier-Fleischer over 2 years ago
Hi Miles,
I am not an expert for geostationary data, but I would guess that the following grid description file should work.
gridtype = projection gridsize = 10201 xsize = 101 ysize = 101 xname = x xlongname = "projection_x_coordinate" xunits = "m" yname = y ylongname = "projection_y_coordinate" yunits = "m" xfirst = -99000 xinc = 2000 yfirst = -2901000 yinc = -2000 # grid_mapping = geostationary grid_mapping_name = geostationary longitude_of_projection_origin = 140.7 latitude_of_projection_origin = 0 perspective_point_height = 35785863 semi_major_axis = 6378137 semi_minor_axis = 6356752.3 sweep_angle_axis = y fixed_angle_axis = x false_northing = 0 false_easting = 0 scanningMode = 64 GeoTransform = -5500000. 2000. 0. 5500000. 0. -2000. proj_params = "+proj=geos +lon_0=140.7 +h=35785863 +x_0=0 +y_0=0 +a=6378137 +b=6356752.3 +units=m +sweep=y +no_defs"
-Karin
RE: Geostationary file support - Added by Miles Sowden over 2 years ago
Thanks Karin,
One step forward. Repeating myself for clarity
1: Get the grid information with " cdo griddes in.nc > grid.txt "
2: Amend the grid information as above (note my file was missing the grid mapping information)
3: Update the grid information with " cdo setgrid,grid.txt in.nc out.nc " (note did this on both 1 and 2 km input grids)
Step above worked on the cropped file but I got an error as below on the larger file
cdo remapnn,Generic_1000.nc Generic_2000.nc Remap_1000.nc
Error (gridDefBoundsGeneric) : Allocation of 3872000000 bytes failed. [ line 2985 file grid.c ]
4: Reproject to my double size grid with " cdo remapnn,Crop_1000.nc Crop_2000.nc Remap_1000.nc "
5: Compare the remapped to origional with " cdo sub Remap_1000.nc Crop_1000.nc Delta_1000.nc "
6: Now compare the change on extrapolation to the change between neighbouring pixels with
" cdo sub -shifty,1 -shiftx,1 Crop_1000.nc Crop_1000.nc Near_1000.nc "
This fails with cdo(2) shiftx: Unsupported grid type: projection
I can then fix this by changing the grid type to generic, do the shift function, change back to projection.
I think the grid functions don't yet know how to fully handle geostationary projection mappings.
RE: Geostationary file support - Added by Miles Sowden over 2 years ago
sorry error message after 3 should actually be after 4
RE: Geostationary file support - Added by Karin Meier-Fleischer over 2 years ago
First, without the data file it is hard to reproduce the problem, but the error message tells you that there is a memory issue when you use remapnn. If you have CDO compiled with OpenMP you can use the -P option https://code.mpimet.mpg.de/projects/cdo/embedded/cdo.pdf#subsection.1.2.4
cdo -P 8 remapnn,Crop_1000.nc Crop_2000.nc Remap_1000.nc
The operators shiftx and shifty can be used with rectilinear and curvilinear grids only.