Hello there,
I've got a netCDF file which I want to read, modify, and save, which works
well with 'ncdf4'-package.
The problem is, that the grid representing the data (temperature values) is
irregular, and I am not able to output an irregular grid. In the original
file I have two variables ('lat' and 'lon'), which are both
two-dimensional
(356 x 236), with each entry representing the latitude or longitude,
respectively, of each cell in the original grid (356 x 236 x 365 [days]).
In the input file, the data grid (TMAX_2M, please see blow) references the
'lat' and 'lon' variables by the argument 'coordinates =
"lon lat"', which,
as far as I understood, make my viewer application (Panoply) make the data
appear in the correct locations.
However, I am not able to hand this argument to the 'ncvar_def' method
(the
error is: 'unused argument'), so at the moment I use only one dimension
of
the 'lat' and 'lon' variables and hand them to the
'ncvar_def'-method, which
produces an distorted (regular) output.
I would love to know how to make the netCDF file aware that the
latitude/longitude coordinates of the data grid cells are stored in the
"external" 'lat' and 'lon' variables, and not in the
grid object itself.
Below I'm posting the structure of my input file, which I hope to reproduce,
as well as my source code so far.
Any help is highly appreciated, thank you very much in advance.
Best regards, Michael
this is the structure of the source netCDF file:
-----------------------------------------------------
dimensions:
x = 356;
y = 236;
height_2m = 1;
time = UNLIMITED; // (365 currently
nb2 = 2;
variables:
float lon(y=236, x=356);
:standard_name = "longitude";
:long_name = "longitude";
:units = "degrees_east";
:_CoordinateAxisType = "Lon";
float lat(y=236, x=356);
:standard_name = "latitude";
:long_name = "latitude";
:units = "degrees_north";
:_CoordinateAxisType = "Lat";
float height_2m(height_2m=1);
:standard_name = "height";
:long_name = "height above the surface";
:units = "m";
:positive = "up";
:axis = "Z";
double time(time=365);
:standard_name = "time";
:long_name = "time";
:bounds = "time_bnds";
:units = "seconds since 1979-01-01 00:00:00";
:calendar = "proleptic_gregorian";
double time_bnds(time=365, nb2=2);
:units = "seconds since 1979-01-01 00:00:00";
:calendar = "proleptic_gregorian";
float TMAX_2M(time=365, height_2m=1, y=236, x=356);
:standard_name = "air_temperature";
:long_name = "2m maximum temperature";
:units = "K";
:coordinates = "lon lat";
:cell_methods = "time: maximum";
:_FillValue = -9.0E33f; // float
this is my source code so far:
---------------------------------
library(ncdf4)
library(maptools)
library(sp)
ncdfname <- "TMAX_2M.nc"
ncdffile <- nc_open(ncdfname)
latorig <- ncvar_get(ncdffile,"lat",verbose=F)
lonorig <- ncvar_get(ncdffile,"lon")
timeorig <- ncvar_get(ncdffile,"time",verbose=F)
tunits <- ncatt_get(ncdffile,"time","units")
# get global attributes
titleorig <- ncatt_get(ncdffile,0,"title")
institutionorig <- ncatt_get(ncdffile,0,"institution")
datasourceorig <- ncatt_get(ncdffile,0,"source")
referencesorig <- ncatt_get(ncdffile,0,"references")
historyorig <- ncatt_get(ncdffile,0,"history")
Conventionsorig <- ncatt_get(ncdffile,0,"Conventions")
temp.array <- ncvar_get(ncdffile,"TMAX_2M")
dlname <- ncatt_get(ncdffile,"TMAX_2M","long_name")
dunits <- ncatt_get(ncdffile,"TMAX_2M","units")
fillvalue <- ncatt_get(ncdffile,"TMAX_2M","_FillValue")
temp.array[temp.array==fillvalue$value] <- NA
temp.array[temp.array==NULL] <- NA
nc_close(ncdffile)
temp.slice <- temp.array[,,1]
lonorig1d = lonorig[,1]
latorig1d = latorig[1,]
lonorig1d_sort = sort(lonorig1d, decreasing=FALSE)
latorig1d_sort = sort(latorig1d, decreasing=FALSE)
dim_lon <- ncdim_def('x', 'degrees_east',
as.double(lonorig1d_sort),
create_dimvar=TRUE)
dim_lat <- ncdim_def('y', 'degrees_north',
as.double(latorig1d_sort),
create_dimvar=TRUE)
dim_height <- ncdim_def('height_2m', 'height', c(0:0))
dim_time <- ncdim_def('time', "seconds since 1979-01-01
00:00:00",
calendar="proleptic_gregorian", unlim=T, c(364:0))
var_out_lon <- ncvar_def('lon', 'degrees_east',
list(dim_lon),
longname="longitude")
var_out_lat <- ncvar_def('lat', 'degrees_north',
list(dim_lat),
longname="latitude")
var_out_timebnds <- ncvar_def('time_bnds', "seconds since
1979-01-01
00:00:00", list(dim_time))
var_out <- ncvar_def('TMAX_2M', 'K',
list(dim_lon,dim_lat,dim_time), 9.e20,
longname="2m maximum temperature")
ncid_out <- nc_create('updated_netcdf.nc', list(var_out, var_out_lon,
var_out_lat))
ncvar_put(ncid_out, var_out, temp.array[,,])
nc_close(ncid_out)
--
View this message in context:
http://r.789695.n4.nabble.com/writing-netCDF-with-irregular-grid-tp4685584.html
Sent from the R help mailing list archive at Nabble.com.