Commit 039728b3 authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

New operator: mrotuvb

parent 9f56f840
2006-11-?? Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
* using CDI library version 1.0.3
* New operator: intntime [request: Michael Boettinger]
* New operator: intntime [request: Michael Boettinger]
* New operator: mrotuvb [request: Helmuth Haak]
* New operator: outputvector [request: Helmuth Haak]
* inttime: extention for months and years [request: Holger Goettel]
* outputcenter: support of zonal and meridional fields
* zonavg: bug fix for non gaussian or lonlat grids
......
......@@ -25,14 +25,17 @@
Plotting example:
- outputcenter
- outputbounds:
- outputboundscpt:
- outputbounds
- outputboundscpt
- outputvector
*/
#if defined (HAVE_CONFIG_H)
# include "config.h" /* VERSION */
#endif
#include <string.h>
#include "cdi.h"
#include "cdo.h"
#include "cdo_int.h"
......@@ -43,9 +46,9 @@
void *Outputgmt(void *argument)
{
static char func[] = "Outputgmt";
int OUTPUTCENTER, OUTPUTBOUNDS, OUTPUTBOUNDSCPT;
int OUTPUTCENTER, OUTPUTBOUNDS, OUTPUTBOUNDSCPT, OUTPUTVECTOR;
int operatorID;
int i;
int i, j;
int varID0, varID, recID;
int gridsize = 0;
int gridID, code;
......@@ -61,12 +64,14 @@ void *Outputgmt(void *argument)
int status;
int lgrid_gen_bounds = FALSE, luse_grid_corner = FALSE;
int zaxisID, taxisID;
int ninc = 1;
int vdate, vtime;
int year, month, day, hour, minute;
char varname[256];
double level;
double missval;
double *array = NULL;
double *uf = NULL, *vf = NULL, *alpha = NULL, *auv = NULL;
double *grid_center_lat = NULL, *grid_center_lon = NULL;
double *grid_corner_lat = NULL, *grid_corner_lon = NULL;
FILE *cpt_fp;
......@@ -77,9 +82,17 @@ void *Outputgmt(void *argument)
OUTPUTCENTER = cdoOperatorAdd("outputcenter", 0, 0, NULL);
OUTPUTBOUNDS = cdoOperatorAdd("outputbounds", 0, 0, NULL);
OUTPUTBOUNDSCPT = cdoOperatorAdd("outputboundscpt", 0, 0, NULL);
OUTPUTVECTOR = cdoOperatorAdd("outputvector", 0, 0, NULL);
operatorID = cdoOperatorID();
if ( operatorID == OUTPUTVECTOR )
{
operatorInputArg("increment");
operatorCheckArgc(1);
ninc = atoi(operatorArgv()[0]);
}
if ( operatorID == OUTPUTBOUNDS || operatorID == OUTPUTBOUNDSCPT )
luse_grid_corner = TRUE;
......@@ -142,6 +155,10 @@ void *Outputgmt(void *argument)
if ( nlon == 1 && nlat > 1 && nlev > 1 ) lzon = TRUE;
if ( nlon > 1 && nlat == 1 && nlev > 1 ) lmer = TRUE;
}
else
{
nlat = 1;
}
if ( cdoVerbose && lzon ) cdoPrint("Process zonal data");
if ( cdoVerbose && lmer ) cdoPrint("Process meridional data");
......@@ -188,6 +205,13 @@ void *Outputgmt(void *argument)
}
array = (double *) malloc(gridsize*sizeof(double));
if ( operatorID == OUTPUTVECTOR )
{
uf = (double *) malloc(gridsize*sizeof(double));
vf = (double *) malloc(gridsize*sizeof(double));
alpha = (double *) malloc(gridsize*sizeof(double));
auv = (double *) malloc(gridsize*sizeof(double));
}
tsID = 0;
nrecs = streamInqTimestep(streamID, tsID);
......@@ -201,6 +225,8 @@ void *Outputgmt(void *argument)
fprintf(stdout, "# Generated by CDO version %s\n", VERSION);
fprintf(stdout, "#\n");
fprintf(stdout, "# Operator = %s\n", cdoOperatorName(operatorID));
fprintf(stdout, "#\n");
fprintf(stdout, "# File = %s\n", cdoStreamName(0));
fprintf(stdout, "# Date = %4.4d-%2.2d-%2.2d\n", year, month, day);
fprintf(stdout, "# Time = %2.2d:%2.2d\n", hour, minute);
......@@ -216,13 +242,13 @@ void *Outputgmt(void *argument)
if ( varID != varID0 ) continue;
if ( recID > 0 && !lzon && !lmer ) continue;
streamReadRecord(streamID, array, &nmiss);
level = zaxisInqLevel(zaxisID, levelID);
fprintf(stdout, "# Level = %g\n", level);
fprintf(stdout, "#\n");
streamReadRecord(streamID, array, &nmiss);
if ( operatorID == OUTPUTCENTER )
{
for ( i = 0; i < gridsize; i++ )
......@@ -239,6 +265,36 @@ void *Outputgmt(void *argument)
}
fprintf(stdout, "#\n");
}
else if ( operatorID == OUTPUTVECTOR )
{
if ( nrecs < 2 ) cdoAbort("No enough fields!");
memcpy(uf, array, gridsize*sizeof(double));
streamInqRecord(streamID, &varID, &levelID);
streamReadRecord(streamID, vf, &nmiss);
for ( j = 0; j < nlat; j += ninc )
for ( i = 0; i < nlon; i += ninc )
{
/* compute length of velocity vector */
auv[IX2D(j,i,nlon)] = sqrt(uf[IX2D(j,i,nlon)]*uf[IX2D(j,i,nlon)] +
vf[IX2D(j,i,nlon)]*vf[IX2D(j,i,nlon)]);
alpha[IX2D(j,i,nlon)] = atan2(vf[IX2D(j,i,nlon)],uf[IX2D(j,i,nlon)]);
alpha[IX2D(j,i,nlon)] = 90. - alpha[IX2D(j,i,nlon)]*RAD2DEG;
if ( alpha[IX2D(j,i,nlon)] < 0 ) alpha[IX2D(j,i,nlon)] += 360;
if ( alpha[IX2D(j,i,nlon)] > 360 ) alpha[IX2D(j,i,nlon)] -= 360;
if ( fabs(auv[IX2D(j,i,nlon)]) > 0 )
fprintf(stdout, " %g %g %g %g\n",
grid_center_lon[IX2D(j,i,nlon)], grid_center_lat[IX2D(j,i,nlon)],
alpha[IX2D(j,i,nlon)], auv[IX2D(j,i,nlon)]);
}
fprintf(stdout, "#\n");
break;
}
else if ( operatorID == OUTPUTBOUNDS || operatorID == OUTPUTBOUNDSCPT )
{
for ( i = 0; i < gridsize; i++ )
......
......@@ -43,21 +43,33 @@ char *strdup(const char *s);
*/
#endif
#ifndef DBL_IS_EQUAL
#ifndef DBL_IS_EQUAL
/*
#define DBL_IS_EQUAL(x,y) (fabs(x - y) <= 2.0*(y*DBL_EPSILON + DBL_MIN))
#define DBL_IS_EQUAL(x,y) (fabs(x - y) <= 2.0*(y*DBL_EPSILON + DBL_MIN))
*/
#define DBL_IS_EQUAL(x,y) (!(fabs(x - y) > 0))
#define DBL_IS_EQUAL(x,y) (!(fabs(x - y) > 0))
#endif
#ifndef M_PI
#define M_PI 3.14159265358979323846 /* pi */
#ifndef M_PI
#define M_PI 3.14159265358979323846 /* pi */
#endif
#ifndef M_LN10
#define M_LN10 2.30258509299404568402 /* log_e 10 */
#ifndef M_LN10
#define M_LN10 2.30258509299404568402 /* log_e 10 */
#endif
#ifndef RAD2DEG
#define RAD2DEG (180./M_PI) /* conversion for rad to deg */
#endif
#ifndef DEG2RAD
#define DEG2RAD (M_PI/180.) /* conversion for deg to rad */
#endif
#define IX2D(y,x,nx) ((y)*(nx)+(x))
#define CDO_EXP_LOCAL 1
#define CDO_EXP_REMOTE 2
......
......@@ -178,7 +178,7 @@ void *Zonstat(void *argument);
#define NinfoOperators {"nyear", "nmon", "ndate", "ntime", "ncode", "nvar", "nlevel"}
#define NmltestOperators {"nmltest"}
#define OutputOperators {"output", "outputint", "outputsrv", "outputext", "outputf", "outputts", "outputfld"}
#define OutputgmtOperators {"outputcenter", "outputbounds", "outputboundscpt"}
#define OutputgmtOperators {"outputcenter", "outputbounds", "outputboundscpt", "outputvector"}
#define PinfoOperators {"pinfo", "pinfov"}
#define RemapOperators {"remap"}
#define RemapgridOperators {"remapcon", "remapbil", "remapbic", "remapdis", "remapdis1", "remapcon1"}
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment