How to rotate longitudes and reorder data?

Help
2013-07-26
2016-12-16
  • Phil Bentley

    Phil Bentley - 2013-07-26

    Hi folks,

    I'm wondering if it's possible to use NCO to shift/rotate longitudes from -180/180 to 0/360 and reorder any dependent variables which reference the longitude dimension.

    I can, for example, use ncks to split the source file into eastern and western hemispheres, thus:

      $ ncks -d lon,0.0,180.0 foo.nc foo_e.nc
      $ ncks -d lon,-180.0,-0.01 foo.nc foo_w.nc
    

    And I can use ncap2 to shift the western hemisphere longitudes into the range 180 - 360, thus:

      $ ncap2 -s "lon=lon+360.0" foo_w.nc foo_w_360.nc
    

    But the final step of recombining the two hemispheres is eluding me.

    I'd be interested to learn if there's a smart(er) way of doing this with NCO.

    TIA

    Phil

     
  • Charlie Zender

    Charlie Zender - 2013-07-26

    Hi Phil,
    This is a great question because the full solution, which you are 2/3s of the way to, requires understanding the numerics and disk-storage conventions of gridded data. The best way to finish the task is to use ncpdq to reorder everything so lon is the record dimension, then use ncrcat to concatenate the west to the east file, then use ncpdq to restore the desired final ordering. There are some examples in the manual http://nco.sf.net/nco.html#ncpdq that show this technique of stitching files together. It automatically handles all variables in the file, so once you get the hang of it it's pretty easy. Should be a FAQ. cz

     
  • Christopher Lynnes

    Depending on your situation, you can also do it within ncap2 by just copying hyperslabs to a new variable and tweaking the corresponding longitudes.  Here's a snippet from one of my ncap2 scripts (note I am rotating a 3-d variable, hence the extra plev dim):

    defdim("lat", 180);
    defdim("lon", 360);
    defdim("plev", 17);
    ta(0:7,:,0:179) = temperature(0:7,:,180:359);
    ta(0:7,:,180:359) = temperature(0:7,:,0:179);
    lat = 0.;
    lon = 0.;
    lat(:) = LatitudeU271;
    lon(0:179) = LongitudeU272(180:359);
    lon(180:359) = LongitudeU272(0:179) + 360.;

     
  • henry Butowsky

    henry Butowsky - 2013-07-26

    HI Guys
    heres  my script :

    /****************************************************/
    lon=lon;
    where(lon >=180.0f) lon=lon-360.0f;
    ?***************************************************/

    /* after running above script need to apply dim changes
       use the command -
    ncks -msa_usr_rdr -d lon,450,899 -d lon,0,449 m1.nc m2.nc

    */

     
  • Charlie Zender

    Charlie Zender - 2013-07-26

    Ahhh, we should have better documentation on the -msa_usr_rdr switch, which, as Henry points out, is perfectly suited to this application and which I am remiss at popularizing. I would add that if you reverse the order of Henry's suggestion, and do ncks first and ncap2 last, then you can use coordinate specifications instead of indices in the ncks arguments. Then the command would apply to any grid.

    Chris' entirely ncap2 solution is also fine though, as he well knows, it must be explicitly coded to handle multiple variables

    Thanks everyone for your suggestions.
    cz.

     
  • Charlie Zender

    Charlie Zender - 2013-07-26

    To my chagrin, I realized that -msa_usr_rdr is one of NCO's secret switches
    (viewable with ncks -secret), mainly because there was no documentation.
    Rotating the longitude grid is a perfect use case!
    I've written a section that explains it, and the -msa switch is no longer
    "secret", as of the current snapshot. No more excuses for not using it.
    Long live MSA! Full example is here, and pasted below for completeness
    http://nco.sf.net/nco.html#msa

    The MSA user-order switch ‘-msa_usr_rdr’ (or ‘-msa_user_order’, both of which shorten to ‘-msa’) requests that the multislabs be output in the user-specified order from the command-line, rather than in the input-file on-disk storage order. This allows the user to perform complex data re-ordering in one operation that would otherwise require cumbersome steps of hyperslabbing, concatenating, and permuting. Consider the recent example of a user who needed to convert datasets stored with the longitude coordinate Lon ranging from [−180,180) to datasets that follow the [0,360) convention.

         % ncks -H -v Lon in.nc
         Lon=-180
         Lon=-90
         Lon=0
         Lon=90
    Although simple in theory, this task requires both mathematics to change the numerical value of the longitude coordinate, data hyperslabbing to split the input on-disk arrays at Greenwich, and data re-ordering within to stitch the western hemisphere onto the eastern hemisphere at the date-line. The ‘-msa’ user-order switch overrides the default that data are output in the same order in which they are stored on-disk in the input file, and instead stores them in the same order as the multi-slabs are given to the command line. This default is intuitive and is not important in most uses. However, the MSA user-order switch allows users to meet their output order needs by specifying multi-slabs in a certain order. Compare the results of default ordering to user-ordering for longitude:

         % ncks -O -H       -v Lon -d Lon,0.,180. -d Lon,-180.,-1.0 in.nc
         Lon=-180
         Lon=-90
         Lon=0
         Lon=90
         % ncks -O -H -msa -v Lon -d Lon,0.,180. -d Lon,-180.,-1.0 in.nc
         Lon=0
         Lon=90
         Lon=-180
         Lon=-90
    The two multi-slabs are the same but they can be presented to screen, or to an output file, in either order. The second example shows how to place the western hemisphere after the eastern hemisphere, although they are stored in the opposite order in the input file.

    With this background, one sees that the following commands suffice to rotate the input file by 180 degrees longitude:

         % ncks -O -v LatLon -msa -d Lon,0.,180. -d Lon,-180.,-1.0 in.nc out.nc
         % ncap2 -O -s 'where(Lon < 0) Lon=Lon+360' out.nc out.nc
         % ncks -C -H -v LatLon ~/nco/data/in.nc
         Lat=-45 Lon=-180 LatLon=0
         Lat=-45 Lon=-90 LatLon=1
         Lat=-45 Lon=0 LatLon=2
         Lat=-45 Lon=90 LatLon=3
         Lat=45 Lon=-180 LatLon=4
         Lat=45 Lon=-90 LatLon=5
         Lat=45 Lon=0 LatLon=6
         Lat=45 Lon=90 LatLon=7
         % ncks -C -H -v LatLon ~/out.nc
         Lat=-45 Lon=0 LatLon=2
         Lat=-45 Lon=90 LatLon=3
         Lat=-45 Lon=180 LatLon=0
         Lat=-45 Lon=270 LatLon=1
         Lat=45 Lon=0 LatLon=6
         Lat=45 Lon=90 LatLon=7
         Lat=45 Lon=180 LatLon=4
         Lat=45 Lon=270 LatLon=5
    There are other workable, valid methods to accomplish this rotation, yet none are simpler nor more efficient than utilizing MSA user-ordering. Some final comments on applying this algorithm: Be careful to specify hemispheres that do not overlap, e.g., by inadvertently specifying coordinate ranges that both include Greenwich. Some users will find using index-based rather than coordinate-based hyperslabs makes this clearer.

     
  • Phil Bentley

    Phil Bentley - 2013-07-27

    Hi guys,

    Thanks for the superb answers. That works a treat. Great to learn about that previously hidden -msa option - I for one am going to find that extremely useful in future.

    NCO is one awesome box of tricks :-)

    Thanks again.

    Phil

     
  • Phil Bentley

    Phil Bentley - 2013-07-27

    Just discovered a minor caveat to this solution. If, like me, your netcdf files also contain a longitude bounds variable then that too can also be adjusted using the ncap2 command. My bounds variable has the customary name of 'lon_bnds', so my ncap2 command now appears as follows:

    ncap2 -O -s 'where(lon<0) lon_bnds=lon_bnds+360; where(lon<0) lon=lon+360' $outfile $outfile
    

    The bounds variable must be adjusted before the longitude coordinate variable, otherwise the (lon<0) test returns zero matches.

    Phil

     
  • Emma John

    Emma John - 2016-11-23

    Hi guys, this is a really helpful thread (and I realise probably considered closed) but I'm having issues trying to do the reverse of this.

    I have data on a lonlat grid from lon 0/360 and I want to convert it to lon -180/180. I've tried to use the above code (see below) and it is giving me the right lon values but it isn't changing the position of the East and West hemispheres, so when I look at the grid values (using 'cdo sinfon ifile' and 'cdo griddes ifile') I have 0-180 immediately followed by -179.75 - -0.25. I've read the userguide and help pages and I can't see why it isn't working.

    ncks -O -msa -d lon,180.25,359.75 -d lon,0.0,180.0 2000-10_MON_WIND.nc 2000-10_MON_WIND_test.nc
    
    ncap2 -O -s "where(lon>180)lon=lon-360" 2000-10_MON_WIND_test.nc 2000-10_MON_WIND_test.nc
    

    I also find that my grid goes from being lonlat to generic, and all the metadata/attributes for my variables and dimensions aren't included in the new output file. I can code these back in, reset the grid and loop it through all my files but if there's a way to not lose them in the first place that would be great.

    I'm also going to have to do bilinear interpolation on this, I have 4 satellites that I need to get on the same spatial resolution, would I be better doing that first or does the order make no difference?

    Any help greatly appreciated

    Cheers,
    Emma

     
  • henry Butowsky

    henry Butowsky - 2016-11-23

    Hi Emma,
    you are nearly there. try the following. (the arg -msa is wrong)

    ncks -O --msa_usr_rdr -d lon,180.25,359.75 -d lon,0.0,180.0 2000-10_MON_WIND.nc 2000-10_MON_WIND_test.nc
    
    ncap2 -O -s "where(lon>180)lon=lon-360" 2000-10_MON_WIND_test.nc 2000-10_MON_WIND_test.nc
    

    The bilinear interpolation function should work on either (-180,180) or (0,260)
    ...Henry

     
  • Emma John

    Emma John - 2016-11-23

    Hi Henry,

    You are a wonderful human being, that worked a treat! I'd thought the abbreviated -msa switch did the same as the extended one.

    Thanks so much!
    Emma

     
  • zizzino

    zizzino - 2016-12-15

    Hi everyone,

    I am in a huge pickle because I have a similar problem. I am tooootally new to cdo, nco, and netcdf. I have been given a netcdf file from the ERA-Interim climate model with resolution 1.5 km. The data grid uses a rotated pole (lat=39.25, long=198.0), and has a lat range of [-0.139, 5.61] and a long range of [344, 352]. The dataset has 365 latitudes, 466 longitudes, and is centered over England.

    For ArcGIS post-processing purposes, I need to convert the lat/long of my dataset into regular lat/long with Datum OSGB36. This task has to be done only once. I have had a few difficulties trying to find a tutorial or something similar, so maybe I am in the right place talking to the right people.

    The attached file shows how the dataset is plotted in Panoply. Btw, how does that software plot it in the correct lat/long? Does it perform any transformation in the background?

    I would appreciate any kind of help.

    Thank you so much!!!

     
    Last edit: zizzino 2016-12-15
  • zizzino

    zizzino - 2016-12-15

    Hi everyone,

    I am in a huge pickle because I have a similar problem. I am tooootally new to cdo, nco, and netcdf. I have been given a netcdf file from the ERA-Interim climate model with resolution 1.5 km. The data grid uses a rotated pole (lat=39.25, long=198.0), and has a lat range of [-0.139, 5.61] and a long range of [344, 352]. The dataset has 365 latitudes, 466 longitudes, and is centered over England.

    For post-processing purposes, I need to convert the lat/long of my dataset into regular lat/long with Datum OSGB36. This task has to be done only once. I have had a few difficulties trying to find a tutorial or something similar, so maybe I am in the right place talking to the right people.

    The attached file shows how the dataset is plotted in Panoply. Btw, how does that software plot it in the correct lat/long? Does it perform any transformation in the background?

    I would appreciate any kind of help.

    Thank you so much!!!

     
    • Charlie Zender

      Charlie Zender - 2016-12-16

      You might use ncremap to define a suitably sized destination grid and then remap these data to that grid. Reading the manual http://nco.sf.net/nco.html#ncremap is the next step.

       
  • Cian Woods

    Cian Woods - 2017-01-24

    Hi all,

    I was struggling with this issue today and found a solution. Code pasted below. Hopefully it will be of help.

    file.nc has longitude dimension between -180 to 180 and I want to shift the file so that the data is on 0 to 360.

    1) Extract eastern (fileE.nc) and western (fileW.nc) longitudes like OP
    ncks -d lon,0.0,180.0 file.nc fileE.nc
    ncks -d lon,-180.0,-0.0 file.nc fileW.nc

    2) Change western (negative) longitudes to be from 180 to 360 like OP
    ncap2 -s "lon=lon+360.0" fileW.nc fileW.360.nc

    3) Chnage record dimension from time to lon in the new files
    ncpdq -O -a lon,time fileW.360.nc fileW.360.nc
    ncpdq -O -a lon,time fileE.nc fileE.nc

    4) Concatenate files along new record dimension (-O to overwrite original file.nc)
    ncrcat fileE.nc fileW.360.nc -O file.nc

    5) Change record dimension back to time
    ncpdq -O -a time,lon file.nc file.nc

    6) remove residual files
    rm fileE.nc fileW.nc fileW.360.nc

    file.nc should now be shifted to a 0 to 360 lon grid

     
  • Charlie Zender

    Charlie Zender - 2017-01-24

    Thanks for posting your solution. However, a much more elegant method is described
    here

     

Log in to post a comment.

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:





No, thanks