/* ncoutput.c
* *
* Emilio 4/11/2001
* Added some comments, prettified some of the layout.
* *
* Functions used to define and write netcdf files, and to define generalized axes.
* This is a part of the Environmental Modeling Utilities (EMU) toolbox.
******************************************************************************/
#include "emu.h"
#ifdef NONETCDF
#include "nonetcdf.c"
#endif
/* external variable defined in function.c */
extern char * pszOutputFilePath;
int errnc;
/* ================================================================= */
void SetupDefaultAxes(BOOL UseTimeAxis, BOOL monmean)
{
int i, mo;
double *lon, *lat, tVal = 1.0, tVals[12];
double lonOffset = 360.0, lonLeftOffset, lonRightOffset;
/* longitude axis */
lon = alloc1d_d(0,nCols-1);
for(i=0; i<nCols; i++)
{
lon[i] = fLonLeft + (i+0.5)*fHScale;
if (lon[i] < 0)
lon[i] += lonOffset;
}
lonLeftOffset = (fLonLeft < 0) ? (lonOffset + fLonLeft) : fLonLeft;
lonRightOffset = (fLonRight < 0) ? (lonOffset + fLonRight) : fLonRight;
FillAxisProp("lon", "longitude", "degrees_east", TRUE, (float) lonLeftOffset, (float) lonRightOffset, nCols, lon, FALSE);
/* latitude axis */
lat = alloc1d_d(0,nRows-1);
for(i=0; i<nRows; i++)
lat[i] = fLatBtm + (i+0.5)*fVScale;
FillAxisProp("lat", "latitude", "degrees_north", TRUE, fLatBtm, fLatTop, nRows, lat, FALSE);
/* time axis */
if (UseTimeAxis)
if (!monmean)
FillTAxisProp("time", "time", "default", TRUE,0.0F,1.0F, NC_UNLIMITED, &tVal, "default", FALSE);
else
{
for (mo=0; mo < 12; mo++)
tVals[mo] = (double) MONTH2HOUR*(mo + 0.501);
FillTAxisProp("time", "time", "default", TRUE, (float) tVals[0], (float) tVals[11], 12, tVals, "default", TRUE);
}
/* WATCH OUT, THE ASSIGNMENT OF THE TVALUE tVal MAY CAUSE PROBLEMS, B/C I USE A MEMORY POINTER ! */
/* free allocated memory */
free1d_d(lon,0,nCols-1);
free1d_d(lat,0,nRows-1);
}
/* ================================================================= */
/* Fill the axis dimension structure with all necessary information; *
* used for creating netcdf files */
void FillAxisProp(char *name, char *lname, char *units, short point_spacing_even,
float min, float max, int length, double *values, BOOL positiveDown)
{
short i, axis;
static char positDown[] = "down", positUp[] = "up";
if (istype(name, "x") || istype(name, "lon"))
axis = 0;
else if (istype(name, "y") || istype(name, "lat"))
axis = 1;
else if (istype(name, "z") || istype(name, "depth") || istype(name, "elev"))
{
axis = 2;
if (positiveDown)
strcpy(stAxes[axis].positiveOrient, positDown);
else
strcpy(stAxes[axis].positiveOrient, positUp);
}
else if (istype(name, "t") || istype(name, "time"))
axis = 3;
stAxes[axis].active = TRUE;
strcpy(stAxes[axis].name, name);
strcpy(stAxes[axis].lname, lname);
strcpy(stAxes[axis].units, units);
stAxes[axis].point_spacing_even = point_spacing_even;
stAxes[axis].valid_range[0] = min;
stAxes[axis].valid_range[1] = max;
stAxes[axis].length = length;
/* this array is freed during clean up of the application, in VariableCleanUp */
stAxes[axis].values = alloc1d_d(0,length-1);
for (i=0; i < length; i++)
stAxes[axis].values[i] = values[i];
/* update axes */
for (axis=0; axis < 4; axis++)
if (stAxes[axis].active != TRUE) stAxes[axis].active = FALSE;
AxesType = 1000*stAxes[0].active + 100*stAxes[1].active
+ 10*stAxes[2].active + stAxes[3].active;
}
/* ================================================================= */
/* Fill the time axis dimension structure with all necessary information; *
* used for creating netcdf files *
* If an UNLIMITED time axis is desired, use NC_UNLIMITED for length */
void FillTAxisProp(char *name, char *lname, char *units, short point_spacing_even,
float min, float max, int length, double *values, char *time_origin, BOOL Tmonmean)
{
static char Tunits1[] = "hour since ", Tunits2[] = "-12-15 00:00:00";
char Tunits[40];
static char Torig1[] = "15-DEC-", Torig2[] = " 00:00:00";
char Torigin[40];
static char TunitsMonmean[] = "hour since 0000-01-01 00:00:00";
static char ToriginMonmean[] = "01-JAN-0000 00:00:00";
static char modulo[] = " "; /* modulo attribute string */
/* char s[20]; */
sprintf(Tunits, "%s%d%s", Tunits1,YEARORIGIN,Tunits2);
sprintf(Torigin, "%s%d%s", Torig1,YEARORIGIN,Torig2);
stAxes[3].monmean = Tmonmean;
if (istype(units, "default"))
if (Tmonmean)
units = TunitsMonmean;
else
units = Tunits;
if (istype(time_origin, "default"))
{
if (Tmonmean)
{
strcpy(stAxes[3].time_origin, ToriginMonmean);
strcpy(stAxes[3].time_modulo, modulo);
}
else
strcpy(stAxes[3].time_origin, Torigin);
}
else
strcpy(stAxes[3].time_origin, time_origin);
FillAxisProp(name, lname, units, point_spacing_even, min, max, length, values, FALSE);
}
/* ================================================================= */
int newNCsetup(int *tID, char *pszOutputFile, char *title, char *history)
{
int ncID; /* netCDF file ID */
char ncfile[FILENAMELEN]; /* NC file to be written */
int xdimID, ydimID, zdimID, tdimID, xID, yID, zID;
size_t timeindx = 0;
/* construct file path and file name for nc file */
strcpy(ncfile, pszOutputFilePath);
strcat(ncfile, "/");
strcat(ncfile, pszOutputFile); /* file name maybe should be handled using Bill's */
strcat(ncfile, ".nc"); /* OUTPUT VARS feature in vars.cfg ! */
errnc = nc_create(ncfile, NC_NOCLOBBER, &ncID);
if (errnc == NC_EEXIST)
ExitError("Couldn't create NC file %s, NC file with that name already exists\n", ncfile);
/* write the Global Attributes */
errnc = nc_put_att_text(ncID, NC_GLOBAL, "title", (short) strlen(title), title);
errnc = nc_put_att_text(ncID, NC_GLOBAL, "history", (short) strlen(history), history);
/* define dimensions and coordinate variables (following conventions) */
if (stAxes[0].active) /* X AXIS */
{
errnc = nc_def_dim(ncID, stAxes[0].name, (size_t)stAxes[0].length, &xdimID);
errnc = nc_def_var(ncID, stAxes[0].name, NC_FLOAT, 1, &xdimID, &xID);
errnc = nc_put_att_text(ncID, xID, "long_name", (size_t)strlen(stAxes[0].lname), stAxes[0].lname);
errnc = nc_put_att_text(ncID, xID, "units", (size_t)strlen(stAxes[0].units), stAxes[0].units);
errnc = nc_put_att_float(ncID, xID, "valid_range", NC_FLOAT, 2, stAxes[0].valid_range);
if (stAxes[0].point_spacing_even)
errnc = nc_put_att_text(ncID, xID, "point_spacing", 4, "even");
}
if (stAxes[1].active) /* Y AXIS */
{
errnc = nc_def_dim(ncID, stAxes[1].name, (size_t)stAxes[1].length, &ydimID);
errnc = nc_def_var(ncID, stAxes[1].name, NC_FLOAT, 1, &ydimID, &yID);
errnc = nc_put_att_text(ncID, yID, "long_name", (size_t)strlen(stAxes[1].lname), stAxes[1].lname);
errnc = nc_put_att_text(ncID, yID, "units", (size_t)strlen(stAxes[1].units), stAxes[1].units);
errnc = nc_put_att_float(ncID, yID, "valid_range", NC_FLOAT, 2, stAxes[1].valid_range);
if (stAxes[1].point_spacing_even)
errnc = nc_put_att_text(ncID, yID, "point_spacing", 4, "even");
}
if (stAxes[2].active) /* Z AXIS */
{
errnc = nc_def_dim(ncID, stAxes[2].name, (size_t)stAxes[2].length, &zdimID);
errnc = nc_def_var(ncID, stAxes[2].name, NC_FLOAT, 1, &zdimID, &zID);
errnc = nc_put_att_text(ncID, zID, "long_name", (size_t)strlen(stAxes[2].lname), stAxes[2].lname);
errnc = nc_put_att_text(ncID, zID, "units", (size_t)strlen(stAxes[2].units), stAxes[2].units);
errnc = nc_put_att_float(ncID, zID, "valid_range", NC_FLOAT, 2, stAxes[2].valid_range);
if (stAxes[2].point_spacing_even)
errnc = nc_put_att_text(ncID, zID, "point_spacing", 4, "even");
errnc = nc_put_att_text(ncID, zID, "positive", (size_t)strlen(stAxes[2].positiveOrient),
stAxes[2].positiveOrient);
}
if (stAxes[3].active) /* T AXIS */
{
errnc = nc_def_dim(ncID, stAxes[3].name, (size_t)stAxes[3].length, &tdimID);
errnc = nc_def_var(ncID, stAxes[3].name, NC_DOUBLE, 1, &tdimID, tID);
errnc = nc_put_att_text(ncID, *tID, "long_name", (size_t)strlen(stAxes[3].lname), stAxes[3].lname);
errnc = nc_put_att_text(ncID, *tID, "units", (size_t)strlen(stAxes[3].units), stAxes[3].units);
errnc = nc_put_att_text(ncID, *tID, "time_origin", (size_t)strlen(stAxes[3].time_origin),
stAxes[3].time_origin);
if (stAxes[3].monmean)
errnc = nc_put_att_text(ncID, *tID, "modulo", (size_t)strlen(stAxes[3].time_modulo),
stAxes[3].time_modulo);
}
/* leave define mode for the new NC file */
errnc = nc_enddef(ncID);
/* write values for coordinate variables; a single value for time */
if (stAxes[0].active)
errnc = nc_put_var_double(ncID, xID, stAxes[0].values);
if (stAxes[1].active)
errnc = nc_put_var_double(ncID, yID, stAxes[1].values);
if (stAxes[2].active)
errnc = nc_put_var_double(ncID, zID, stAxes[2].values);
if (stAxes[3].active)
{
if (stAxes[3].length == NC_UNLIMITED)
errnc = nc_put_var1_double(ncID, *tID, &timeindx, stAxes[3].values);
else
errnc = nc_put_var_double(ncID, *tID, stAxes[3].values);
}
return(ncID);
}
/* ================================================================= */
/* Setup new netcdf file from nc headers */
void CreateCDFFromCDL(char * pszCDL, char * pszCDF, char * pszTitle, char * pszHistory,
char * timeunits, char * timeorigin, char * pszType)
{
int ncid; /* netCDF file ID */
static char szModulo[] = " "; /* modulo attribute string */
char s[FILENAMELEN];
int nTimeID;
sprintf(s, "%s -o %s %s", pszNCGenFile, pszCDF, pszCDL);
system(s);
ncid = ncopen(pszCDF, NC_WRITE);
ncredef(ncid);
/* set up time axis */
nTimeID = ncvarid(ncid, "time");
ncattput(ncid, nTimeID, "units", NC_CHAR, (short) strlen(timeunits), timeunits);
ncattput(ncid, nTimeID, "time_origin", NC_CHAR, (short) strlen(timeorigin), timeorigin);
if (istype(pszType, "monmean"))
ncattput(ncid, nTimeID, "modulo", NC_CHAR, (short) strlen(szModulo), szModulo);
/* write the Global Attributes 'title' and 'history' (char arrays) */
if (pszTitle != NULL && !strcmp(pszTitle, ""))
ncattput(ncid, NC_GLOBAL, "title", NC_CHAR, (short) strlen(pszTitle)+1, pszTitle);
if (pszHistory != NULL && !strcmp(pszHistory, ""))
ncattput(ncid, NC_GLOBAL, "history", NC_CHAR, (short) strlen(pszHistory)+1, pszHistory);
ncendef(ncid); /* leave define mode for ncid file */
ncclose(ncid); /* leave define mode for ncid file */
}
/* ================================================================= */
/* DEFINE VARIABLE TO NC FILE */
int DefNCvar(int ncID, char *name, emu_type EmuType, char *dims, char *longname, char *units)
{
/* dims: xyt or xyz or xy or xzt or zt or t */
int varID, old_fill_mode; int bug=0;
int dimids210[] = {2,1,0}, dimids21[] = {2,1}, dimids10[] = {1,0}, dimid2=2, dimid1=1, dimid0=0;
nc_type ncType;
float fFillVal = -999.0;
short sFillVal = -999;
BYTE bFillVal = 255;
errnc = nc_redef(ncID);
errnc = nc_set_fill(ncID, NC_NOFILL, &old_fill_mode);
ncType = EMUtoNCType(EmuType);
if (istype(dims, "xyt") || istype(dims, "xyz") || istype(dims, "xzt"))
errnc = nc_def_var(ncID, name, ncType, 3, dimids210, &varID);
else if (istype(dims, "xy"))
{
if (AxesType == XYT">XYT || AxesType == XYZ">XYZ || AxesType == XY)
errnc = nc_def_var(ncID, name, ncType, 2, dimids10, &varID);
}
else if (istype(dims, "zt"))
{
/* is the use for XZT axes wrong?? Should it be dimids21 ?? */
if (AxesType == ZT">XZT || AxesType == ZT)
errnc = nc_def_var(ncID, name, ncType, 2, dimids10, &varID);
}
else if (istype(dims, "z"))
{
if (AxesType == XYZ)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid2, &varID);
else if (AxesType == XZT)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid1, &varID);
else if (AxesType == ZT)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid0, &varID);
}
else if (istype(dims, "t"))
{
if (AxesType == XYT || AxesType == XZT)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid2, &varID);
else if (AxesType == ZT)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid1, &varID);
T COLOR="#FF0000">elseT> T COLOR="#FF0000">ifT> (Type">AxesType == T)
errnc = nc_def_var(ncID, name, ncType, 1, &dimid0, &varID);
}
/* add variable attributes */
errnc = nc_put_att_text(ncID, varID, "long_name", (size_t) strlen(longname), longname);
errnc = nc_put_att_text(ncID, varID, "units", (size_t) strlen(units), units);
/* maybe _FillValue should be optional, to be inserted only if it is known
that there are no data values in the variable */
if (EmuType == EMU_SHORT || EmuType == EMU_USHORT)
errnc = nc_put_att_short(ncID, varID, "_FillValue", NC_SHORT, 1, &sFillVal);
else if (EmuType == EMU_BYTE)
errnc = nc_put_att_uchar(ncID, varID, "_FillValue", NC_BYTE, 1, &bFillVal);
else
errnc = nc_put_att_float(ncID, varID, "_FillValue", NC_FLOAT, 1, &fFillVal);
errnc = nc_enddef(ncID);
return(varID);
}
/* ================================================================= */
/* Add variables to nc file setup from a template (dimvars.cdl) */
void DefNCvarsFromNCs(int ncIDin, char * pszNC, emu_type EmuType)
{
int v, i; int bug=0;
int fromncID; /* netCDF file IDs */
int varID, fromvarID; /* variable ID */
int dimids[] = {2,1,0}; /* dimension ids: time, lat, lon */
int old_fill_mode;
nc_type ncType;
int nAtts;
char attname[MAX_NC_NAME];
char fromncfile[FILENAMELEN]; /* netCDF header file */
int ncID;
if (istype(pszNC, "none") || EmuType == 0)
ncID = ncIDin;
else
errnc = nc_open(pszNC, NC_WRITE, &ncID);
errnc = nc_redef(ncID);
errnc = nc_set_fill(ncID, NC_NOFILL, &old_fill_mode);
/* in netCDF file, create all variables and copy attributes to them */
for (v = 0; v < nVars; v++)
{
if (bDebug) printf("HERE: %d var-parameter: %s\n",bug++, pstVars[v].name);
/* open variable header netCDF file, "fromncfile" */
if (istype(pstVars[v].filetype, "nc"))
{
fromncID = pstVars[v].ncid;
fromvarID = pstVars[v].var_id;
}
else
{ /* use the nc header file found together with the data files */
strcpy(fromncfile, pstVars[v].path);
strcat(fromncfile, pstVars[v].name);
strcat(fromncfile, ".nc");
if (bDebug) printf("HERE: %d fromncfile path: %s\n",bug++, fromncfile);
errnc = nc_open(fromncfile, NC_NOWRITE, &fromncID);
errnc = nc_inq_varid(fromncID, pstVars[v].name, &fromvarID);
}
/* define the variable in ncID */
if (EmuType == 0)
errnc = nc_inq_vartype(fromncID, fromvarID, &ncType);
else
ncType = EMUtoNCType(EmuType);
errnc = nc_inq_varnatts(fromncID, fromvarID, &nAtts);
errnc = nc_def_var(ncID, pstVars[v].name, ncType, 3, dimids, &varID);
/* copy variable attributes from "fromvarID" to "varID" */
for (i=0; i<nAtts; i++)
{
errnc = nc_inq_attname(fromncID, fromvarID, i, attname);
errnc = nc_copy_att(fromncID, fromvarID, attname, ncID, varID);
}
if (!istype(pstVars[v].filetype, "nc"))
errnc = nc_close(fromncID);
} /* END v FOR-LOOP */
if (istype(pszNC, "none") || EmuType == 0)
errnc = nc_enddef(ncID);
else
errnc = nc_close(ncID);
}
/* ================================================================= */
void AddNCCompressionAttributes(int ncID, int varID, float scale_factor, float add_offset)
{
float sFillVal = -999;
/* Open Define Mode */
errnc = nc_redef(ncID);
/* NC data type of _FillValue & missing_value and of the variable should be the same */
/* need to use "missing_value" too, for Ferret */
errnc = nc_put_att_float(ncID, varID, "missing_value", NC_SHORT, 1, &sFillVal);
/* unpacked_data = scale_factor * packed_data + add_offset */
errnc = nc_put_att_float(ncID, varID, "scale_factor", NC_FLOAT, 1, &scale_factor);
errnc = nc_put_att_float(ncID, varID, "add_offset", NC_FLOAT, 1, &add_offset);
errnc = nc_enddef(ncID);
}
/* ================================================================= */
void AddAttribute(int ncIDin, char *pszNC, int nVarIDin, char *pszVar,
char *pszAtt, emu_type EmuType, int nLen, void *pvValues)
{
/* use "none" for pszNC and/or pszVar if the IDs are supplied instead of the names */
int nVarID, ncID;
if (istype(pszNC, "none"))
ncID = ncIDin;
else
errnc = nc_open(pszNC, NC_WRITE, &ncID);
/* use NC_GLOBAL in nVarIDin to add a global attribute */
if (istype(pszVar, "none"))
nVarID = nVarIDin;
else
errnc = nc_inq_varid(ncID, pszVar, &nVarID);
errnc = nc_redef(ncID);
ncattput(ncID, nVarID, pszAtt, EMUtoNCType(EmuType), nLen, pvValues);
/* errnc = nc_put_att_text(ncID, nVarID, pszAtt, (size_t) strlen(units), units); */
if (istype(pszNC, "none") || EmuType == 0)
errnc = nc_enddef(ncID);
else
errnc = nc_close(ncID);
}