Project

General

Profile

Forums » Open discussion »

remapcon on lonlat gridtype

Added by esmail ghaemi almost 5 years ago

Hi!
I'm using CDO in my research and it is amazing!
Right now I'm trying to upscale my data to 1 km.
I tried to create grid file based on my original dataset (which has 200m resolution) with cdo griddes. For this, I just changed grid size to 1000m and grid number to 330 in the produced grid file.

but it just returns nan values.

I also tried without changing grid file (with the original grid file but the results are the same.

I use remapcon,grid.txt sample.nc upscale.nc to create new file.

I've attached grid.txt (the original one) and also a short period of my dataset.

Thank you very much.

Download all files
grid.txt (275 Bytes) grid.txt original gridtext
sample.nc (188 KB) sample.nc (part of) original data

Replies (8)

RE: remapcon on lonlat gridtype - Added by Karin Meier-Fleischer almost 5 years ago

Hi Esmail,

first, there are variables (time_bnds, UTM33N) missing in the file sample.nc which will return warnings using 'cdo sinfon'. They are set as attributes for time and Precip.

Delete the attributes:

ncatted -O -a grid_mapping,Precip,d,, -a bounds,time,d,, sample.nc sample_2.nc

I'm not pretty sure if the following is correct but I found the additional grid mapping settings for UTM33N at https://spatialreference.org/ref/epsg/wgs-84-utm-zone-33n/prettywkt/

my_gridfile.txt:

gridtype  = lonlat
gridsize  = 330
datatype  = float
xsize     = 22
ysize     = 15
xname     = X
xlongname = "X coordinate in UTM33N" 
xunits    = "meter" 
yname     = Y
ylongname = "Y coordinate in UTM33N" 
yunits    = "meter" 
xfirst    = 557600
xinc      = 1000
yfirst    = 5190600
yinc      = 1000
grid_mapping = crs
grid_mapping_name = transverse_mercator
latitude_of_origin = 0
central_meridian   = 15
scale_factor       = 0.9996
false_easting      = 500000
false_northing     = 0

Then do the remapping:

cdo -remapnn,my_gridfile.txt sample_2.nc out.nc

cdo sinfon out.nc

   File format : NetCDF4
    -1 : Institut Source   T Steptype Levels Num    Points Num Dtype : Parameter name
     1 : unknown  unknown  v instant       1   1       330   1  F32  : Precip        
   Grid coordinates :
     1 : projection               : points=330 (22x15)
                          mapping : transverse_mercator
                                X : 557600 to 578600 by 1000 meter
                                Y : 5190600 to 5204600 by 1000 meter
   Vertical coordinates :
     1 : surface                  : levels=1
   Time coordinate :  25 steps
     RefTime =  2007-01-01 00:00:00  Units = minutes  Calendar = gregorian
  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss  YYYY-MM-DD hh:mm:ss
  2016-08-16 20:30:00  2016-08-16 20:35:00  2016-08-16 20:40:00  2016-08-16 20:45:00
  2016-08-16 20:50:00  2016-08-16 20:55:00  2016-08-16 21:00:00  2016-08-16 21:05:00
  2016-08-16 21:10:00  2016-08-16 21:15:00  2016-08-16 21:20:00  2016-08-16 21:25:00
  2016-08-16 21:30:00  2016-08-16 21:35:00  2016-08-16 21:40:00  2016-08-16 21:45:00
  2016-08-16 21:50:00  2016-08-16 21:55:00  2016-08-16 22:00:00  2016-08-16 22:05:00
  2016-08-16 22:10:00  2016-08-16 22:15:00  2016-08-16 22:20:00  2016-08-16 22:25:00
  2016-08-16 22:30:00
cdo    sinfon: Processed 1 variable over 25 timesteps [0.01s 16MB].

-Karin

RE: remapcon on lonlat gridtype - Added by esmail ghaemi almost 5 years ago

Dear Karin
Thank you very much for your quick reply.
Actually I want to use conservative remapping, not bilinear.
In order to check, I did exactly the method that you described above(also used remapbil), but I got only NaN values in the output file.
I'm using cdo version 1.9.8.
Thanks again.
Esmail

RE: remapcon on lonlat gridtype - Added by Karin Meier-Fleischer almost 5 years ago

You can't use remapcon because your coordinate units are meter. See remapcon documentation: 'This module contains operators for a first order conservative remapping of fields between grids in spherical coordinates.'

Use remapnn.

cdo info out.nc

    -1 :       Date     Time   Level Gridsize    Miss :     Minimum        Mean     Maximum : Parameter ID
     1 : 2016-08-16 20:30:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     2 : 2016-08-16 20:35:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     3 : 2016-08-16 20:40:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     4 : 2016-08-16 20:45:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     5 : 2016-08-16 20:50:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     6 : 2016-08-16 20:55:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     7 : 2016-08-16 21:00:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
     8 : 2016-08-16 21:05:00       0      330       0 :     0.19336     0.19336     0.19336 : -1            
     9 : 2016-08-16 21:10:00       0      330       0 :     0.39258     0.39258     0.39258 : -1            
    10 : 2016-08-16 21:15:00       0      330       0 :      8.9443      8.9443      8.9443 : -1            
    11 : 2016-08-16 21:20:00       0      330       0 :      7.2119      7.2119      7.2119 : -1            
    12 : 2016-08-16 21:25:00       0      330       0 :      1.2998      1.2998      1.2998 : -1            
    13 : 2016-08-16 21:30:00       0      330       0 :     0.61523     0.61523     0.61523 : -1            
    14 : 2016-08-16 21:35:00       0      330       0 :     0.37012     0.37012     0.37012 : -1            
    15 : 2016-08-16 21:40:00       0      330       0 :     0.26855     0.26855     0.26855 : -1            
    16 : 2016-08-16 21:45:00       0      330       0 :     0.22559     0.22559     0.22559 : -1            
    17 : 2016-08-16 21:50:00       0      330       0 :    0.060547    0.060547    0.060547 : -1            
    18 : 2016-08-16 21:55:00       0      330       0 :    0.090820    0.090820    0.090820 : -1            
    19 : 2016-08-16 22:00:00       0      330       0 :     0.22754     0.22754     0.22754 : -1            
    20 : 2016-08-16 22:05:00       0      330       0 :     0.11719     0.11719     0.11719 : -1            
    21 : 2016-08-16 22:10:00       0      330       0 :    0.086914    0.086914    0.086914 : -1            
    22 : 2016-08-16 22:15:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
    23 : 2016-08-16 22:20:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
    24 : 2016-08-16 22:25:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
    25 : 2016-08-16 22:30:00       0      330       0 :      0.0000      0.0000      0.0000 : -1            
cdo    info: Processed 1 variable over 25 timesteps [0.01s 16MB].

RE: remapcon on lonlat gridtype - Added by esmail ghaemi almost 5 years ago

Thanks again!

That's strange, because I did the exact procedure but I get nan values.

About the method, I thought maybe I can do something with nco to justify the grid type and then do the conservative remapping. Probably I should think about other ways.

Cheers.:)

RE: remapcon on lonlat gridtype - Added by Karin Meier-Fleischer almost 5 years ago

Which CDO version (cdo -V) do you use?

RE: remapcon on lonlat gridtype - Added by esmail ghaemi almost 5 years ago

I'm using cdo 1.9.8, but I also tested it with the 1.9.6 version.
Just for clarification, I get nan values for remapbil.
But for remapnn and remapdis, I get just one value in all pixels (for each time-step).
Thanks a lot.

RE: remapcon on lonlat gridtype - Added by esmail ghaemi over 4 years ago

So after a while, I finally found out how to do it. I'm gonna write it here.
In the first step, both datasets should be in the same coordinate system. This can be done easily by using pyproj.
In the next step (if it is needed), coordinates are defined using CDO:
ncatted -a units,lon,o,c,'degreeE' -a units,lat,o,c,'degreeN' des_1.nc
Then, the corners of the original data and also destination data should be calculated. Below I explain it for one dataset. You should also do the two following steps to get the corners of the original dataset.
I did it using an NCL code:
begin

_ dfile = addfile("des_1.nc","r")
dstGridName = "des_corners.nc" 

Opt                = True
      Opt@ForceOverWrite = True
      Opt@PrintTimings   = True

;;---If we have 2D lat/lon arrays.
; Mask2D is necessary if dst lat/lon has missing values.
Opt@Mask2D = where(.not.ismissing(dfile->lat).and.\
.not.ismissing(dfile->lon),1,0)
curvilinear_to_SCRIP(dstGridName,dfile->lat,dfile->lon,Opt)

;;---Clean up
delete(Opt)

end_

After that, corners should be added to the des_1.nc using CDO:
cdo -f nc -setgrid,des_corners.nc des_1.nc des_1_with_corners.nc

Now using cdo griddes des_1_with_corners.nc>destination_grid.txt we can generate txt file for the next step, which is regredding using conservative method:
CDO remapcon,destination_grid.txt original_data_with_corners.nc new_conserv.nc

Cheers,:)

    (1-8/8)