Regridding netCDF

Help
2013-12-04
2013-12-05
  • Maria Catala
    Maria Catala
    2013-12-04

    Hello all,

    I have a netcdf file with one variable called Z500 (geopotential height at 500 hPa) and three dimensions (latitude, longitude and time). The longitude size is 26 (from 28 to 60 by 2), the latitude size is 17 (from -20 to 30 by 2) and the time size is 205996 (150 years with 6 hourly measurements).

    My grid has a spatial resolution of 2 degrees and I want to regrid it in order to have a spatial resolution of 2.5 degrees. The new latitude should have a size of 13 (from 30 to 60 by 2.5) and the longitude should have a size of 19 (from -17.5 to 27.5 by 2.5).

    In order to do that I am using ncap2 with the folowing script (called Interpol.nco):

    // Bilinear interpolation
    // Included by rgr.sh

    defdim("latn",13);
    defdim("lonn",19);
    latn[$latn] = {30, 32.5, 35, 37.5, 40, 42.5, 45, 47.5, 50, 52.5, 55, 57.5, 60};
    lonn[$lonn] = {-17.5, -15, -12.5, -10, -7.5, -5, -2.5, 0, 2.5, 5, 7.5, 10, 12.5, 15, 17.5, 20 ,22.5, 25, 27.5};

    out[$time,$latn,$lonn]=0.0;

    // Bi-linear interpolation
    Z500f = bilinear_interp_wrap(Z500,out,latn,lonn,lat,lon);

    After that, in the console:

    ncap2 -0 -S Interpol.nco Hgt_500_1871-2011_3d.nc test.nc

    This is the error message that I get:

    ncap2: ERROR bilinear_interp_wrap(): Dimension size mismatch with input variables

    I thought it was because of the time dimension so I tried to write the following code into the script:
    i=1
    until ($i == dim(time))
    do Z500f[$i,,] = bilinear_interp_wrap(Z500[$i,,],out[$i,,],latn,lonn,lat,lon);
    $i = $i + 1
    done

    But I don't know if it is a good idea or the good syntax.

    Can someone help me?

    Thanks!

     
    Last edit: Maria Catala 2013-12-04
  • henry Butowsky
    henry Butowsky
    2013-12-04

    Hi Maria,
    you are nearly there: try the following script.

    /******* interpol.nco **********/
    defdim("latn",13);
    defdim("lonn",19);
    latn[$latn] = {30, 32.5, 35, 37.5, 40, 42.5, 45, 47.5, 50, 52.5, 55, 57.5, 60};
    lonn[$lonn] = {-17.5, -15, -12.5, -10, -7.5, -5, -2.5, 0, 2.5, 5, 7.5, 10, 12.5, 15, 17.5, 20 ,22.5, 25, 27.5};

    out[$latn,$lonn]=0.0;
    Z500_new[$time,$latn,$lonn]=0.0;

    for( idx=0 ; idx<$time.size ; idx++){
    Z500_new(idx,:,:)=bilinear_interp(Z500(idx,:,:),out,latn,lonn,lat,lon);
    }

    /*******************/

    Then run the script as follows
    ncap2 -v -S interpol.nco in.nc out.nc

    ... Henry

     
  • Maria Catala
    Maria Catala
    2013-12-05

    Thanks for your answer! This is exactly what I tried to do (but I am only used to R language)

    I got several problems, this is what I have with ncdump -h in.nc:

    dimensions:

    lon = 26 ;
    lat = 17 ;
    time = UNLIMITED ; // (5844 currently)
    

    variables:

    float lon(lon) ;
        lon:standard_name = "longitude" ;
        lon:long_name = "Longitude" ;
        lon:units = "degrees_east" ;
        lon:axis = "X" ;
    
    float lat(lat) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "Latitude" ;
        lat:units = "degrees_north" ;
        lat:axis = "Y" ;
    
    double time(time) ;
        time:standard_name = "time" ;
        time:units = "hours since 1800-01-01 00:00:00" ;
    
        time:calendar = "standard" ;
    float Z500(time, lat, lon) ;
        Z500:standard_name = "geopotential_height" ;
        Z500:long_name = "6-hourly Geopotential Heights on Pressure Levels" ;
        Z500:units = "m" ;
        Z500:_FillValue = -32767.f ;
        Z500:unpacked_valid_range = -1500.f, 60000.f ;
        Z500:actual_range = 4866.f, 5998.f ;
        Z500:precision = 0s ;
        Z500:GRIB_id = 7s ;
        Z500:GRIB_name = "HGT" ;
        Z500:var_desc = "Geopotential height" ;
        Z500:dataset = "NOAA-CIRES 20th Century Reanalysis version 2" ;
        Z500:level_desc = "Pressure Levels" ;
        Z500:statistic = "Ensemble Mean" ;
        Z500:parent_stat = "Other" ;
    

    This is what I have with the outpout file:

    dimensions:

    latn = 13 ;
    lonn = 19 ;
    time = UNLIMITED ; // (5844 currently)
    lat = 17 ;
    lon = 26 ;
    

    variables:

    double Z500_new(time, latn, lonn) ;
    
    double latn(latn) ;
    
    double lonn(lonn) ;
    
    double out(latn, lonn) ;
    
    int idx ;
    
    float lon(lon) ;
    
        lon:standard_name = "longitude" ;
        lon:long_name = "Longitude" ;
        lon:units = "degrees_east" ;
        lon:axis = "X" ;
    
    float lat(lat) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "Latitude" ;
        lat:units = "degrees_north" ;
        lat:axis = "Y" ;
    
    double time(time) ;
        time:standard_name = "time" ;
        time:units = "hours since 1800-01-01 00:00:00" ;
        time:calendar = "standard" ;
    

    In fact I want to find the same format as the input file. lonn, latn and Z500 are considered as double and not as float, can I change that ? Moreover Z500_new[,,1] is a grid with 0 (after this time step I have normal values).

     
    Last edit: Maria Catala 2013-12-05
  • henry Butowsky
    henry Butowsky
    2013-12-05

    Hi Maria,
    your editing your previous posts - its getting confusing. Can just reply/post instead of editing your posts.

    The existence of the tmp file indicates that the script has run incorrectly.
    Can you try doing a bilineariterp on the first time step in the series. Just to check things ate working :

    /***********/
    defdim("latn",13);
    defdim("lonn",19);
    // change to your new range here
    latn[$latn] = {30, 32.5, 35, 37.5, 40, 42.5, 45, 47.5, 50, 52.5, 55, 57.5, 60};
    // change to your new range here
    lonn[$lonn] = {-17.5, -15, -12.5, -10, -7.5, -5, -2.5, 0, 2.5, 5, 7.5, 10, 12.5, 15, 17.5, 20 ,22.5, 25, 27.5};
    out[$latn,$lonn]=0.0;

    // no time dimension
    Z500_new[$latn,$lonn]=0.0;

    Z500_new=bilinear_interp(Z500(0,:,:),out,latn,lonn,lat,lon);

    /***/

    ... Henry

     
    Last edit: henry Butowsky 2013-12-05
  • Maria Catala
    Maria Catala
    2013-12-05

    I am sorry Henry, I didn't see your post!

    I had 150 years data so it really took time. I used another file to test (with only three years data) and it works well. The weird thing is when I add the line "print(Z500_new)" the "for" never ends ! But when the line is not here, it works well...

    The output I get is the following:

    dimensions:

    latn = 13 ;
    
    lonn = 19 ;
    
    time = UNLIMITED ; // (5844 currently)
    
    lat = 17 ;
    
    lon = 26 ;
    

    variables:

    double Z500_new(time, latn, lonn) ;
    
    double latn(latn) ;
    
    double lonn(lonn) ;
    
    double out(latn, lonn) ;
    
    int idx ;
    
    float lon(lon) ;
        lon:standard_name = "longitude" ;
        lon:long_name = "Longitude" ;
        lon:units = "degrees_east" ;
        lon:axis = "X" ;
    
    float lat(lat) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "Latitude" ;
        lat:units = "degrees_north" ;
        lat:axis = "Y" ;
    
    double time(time) ;
        time:standard_name = "time" ;
        time:units = "hours since 1800-01-01 00:00:00" ;
        time:calendar = "standard" ;
    

    I used ncrcat -v var -x out.nc out2.nc to delete lon, lat, idx and out variables
    I also used ncwa -a dim out2.nc out3.nc to delete lon and lat dimensions. After that I have:

    dimensions:

    latn = 13 ;
    
    lonn = 19 ;
    
    time = UNLIMITED ; // (5844 currently)
    

    variables:

    double Z500_new(time, latn, lonn) ;
    
    double latn(latn) ;
    
    double lonn(lonn) ;
    
    double time(time) ;
        time:standard_name = "time" ;
        time:units = "hours since 1800-01-01 00:00:00" ;
        time:calendar = "standard" ;
    

    I will use this file in a program and usually it is float or short variables. But I think it will be ok with a double type.

    Now I will try with the "big" file, I'll tell in one hour...

     
    Last edit: Maria Catala 2013-12-05
    • Maria Catala
      Maria Catala
      2013-12-05

      Wow in fact, when I do it with 3 years it takes 30 seconds. But when I do it with 158 years, I've calculated that it will take 17 hours!! I don't know why there is such a difference... I'll let the function run this night...

      Edit: It works well thanks but is there a way to keep the type of the variables? (Short or float for coordinates)

       
      Last edit: Maria Catala 2013-12-06