Skip to content
GitLab
Commit b525e00d authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

added operator runvar1 and runstd1

parent 0bedf487
Hide whitespace changes
Inline Side-by-side
ChangeLog
View file @ b525e00d
......@@ -3,6 +3,11 @@
* using CDI library version 1.6.0
* Version 1.6.0 released
2013-02-06 Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
* New operator: runvar1 - Running variance [Divisor is (n-1)]
* New operator: runstd1 - Running standard deviation [Divisor is (n-1)]
2013-02-05 Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
* ensrkhisttime: fixed memory fault [https://code.zmaw.de/boards/1/topics/1657]
......
src/Runstat.c
View file @ b525e00d
......@@ -24,7 +24,9 @@
Runstat runmean Running mean
Runstat runavg Running average
Runstat runvar Running variance
Runstat runvar1 Running variance [Divisor is (n-1)]
Runstat runstd Running standard deviation
Runstat runstd1 Running standard deviation [Divisor is (n-1)]
*/
#include <cdi.h>
......@@ -123,7 +125,6 @@ void *Runstat(void *argument)
int i;
int varID;
int recID;
int gridID;
int nrecs, nrecords;
int levelID;
int tsID;
......@@ -134,8 +135,10 @@ void *Runstat(void *argument)
int nmiss;
int nvars, nlevel;
int *recVarID, *recLevelID;
int lmean = FALSE, lvarstd = FALSE, lstd = FALSE;
int *imask;
double missval;
double divisor;
field_t ***vars1 = NULL, ***vars2 = NULL, ***samp1 = NULL;
dtinfo_t *dtinfo;
int taxisID1, taxisID2;
......@@ -182,7 +185,9 @@ void *Runstat(void *argument)
cdoOperatorAdd("runmean", func_mean, 0, NULL);
cdoOperatorAdd("runavg", func_avg, 0, NULL);
cdoOperatorAdd("runvar", func_var, 0, NULL);
cdoOperatorAdd("runvar1", func_var1, 0, NULL);
cdoOperatorAdd("runstd", func_std, 0, NULL);
cdoOperatorAdd("runstd1", func_std1, 0, NULL);
operatorID = cdoOperatorID();
operfunc = cdoOperatorF1(operatorID);
......@@ -190,6 +195,11 @@ void *Runstat(void *argument)
operatorInputArg("number of timesteps");
ndates = atoi(operatorArgv()[0]);
lmean = operfunc == func_mean || operfunc == func_avg;
lstd = operfunc == func_std || operfunc == func_std1;
lvarstd = operfunc == func_std || operfunc == func_var || operfunc == func_std1 || operfunc == func_var1;
divisor = operfunc == func_std1 || operfunc == func_var1;
streamID1 = streamOpenRead(cdoStreamName(0));
vlistID1 = streamInqVlist(streamID1);
......@@ -215,53 +225,16 @@ void *Runstat(void *argument)
vars1 = (field_t ***) malloc((ndates+1)*sizeof(field_t **));
if ( !runstat_nomiss )
samp1 = (field_t ***) malloc((ndates+1)*sizeof(field_t **));
if ( operfunc == func_std || operfunc == func_var )
if ( lvarstd )
vars2 = (field_t ***) malloc((ndates+1)*sizeof(field_t **));
for ( its = 0; its < ndates; its++ )
{
vars1[its] = (field_t **) malloc(nvars*sizeof(field_t *));
vars1[its] = field_malloc(vlistID1, FIELD_PTR);
if ( !runstat_nomiss )
samp1[its] = (field_t **) malloc(nvars*sizeof(field_t *));
if ( operfunc == func_std || operfunc == func_var )
vars2[its] = (field_t **) malloc(nvars*sizeof(field_t *));
for ( varID = 0; varID < nvars; varID++ )
{
gridID = vlistInqVarGrid(vlistID1, varID);
gridsize = gridInqSize(gridID);
nlevel = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
missval = vlistInqVarMissval(vlistID1, varID);
vars1[its][varID] = (field_t *) malloc(nlevel*sizeof(field_t));
if ( !runstat_nomiss )
samp1[its][varID] = (field_t *) malloc(nlevel*sizeof(field_t));
if ( operfunc == func_std || operfunc == func_var )
vars2[its][varID] = (field_t *) malloc(nlevel*sizeof(field_t));
for ( levelID = 0; levelID < nlevel; levelID++ )
{
vars1[its][varID][levelID].grid = gridID;
vars1[its][varID][levelID].nmiss = 0;
vars1[its][varID][levelID].missval = missval;
vars1[its][varID][levelID].ptr = (double *) malloc(gridsize*sizeof(double));
if ( !runstat_nomiss )
{
samp1[its][varID][levelID].grid = gridID;
samp1[its][varID][levelID].nmiss = 0;
samp1[its][varID][levelID].missval = missval;
samp1[its][varID][levelID].ptr = (double *) malloc(gridsize*sizeof(double));
}
if ( operfunc == func_std || operfunc == func_var )
{
vars2[its][varID][levelID].grid = gridID;
vars2[its][varID][levelID].nmiss = 0;
vars2[its][varID][levelID].missval = missval;
vars2[its][varID][levelID].ptr = (double *) malloc(gridsize*sizeof(double));
}
}
}
samp1[its] = field_malloc(vlistID1, FIELD_PTR);
if ( lvarstd )
vars2[its] = field_malloc(vlistID1, FIELD_PTR);
}
gridsizemax = vlistGridsizeMax(vlistID1);
......@@ -314,7 +287,7 @@ void *Runstat(void *argument)
}
}
if ( operfunc == func_std || operfunc == func_var )
if ( lvarstd )
{
farmoq(&vars2[tsID][varID][levelID], vars1[tsID][varID][levelID]);
#if defined (_OPENMP)
......@@ -342,7 +315,7 @@ void *Runstat(void *argument)
otsID = 0;
while ( TRUE )
{
if ( operfunc == func_mean || operfunc == func_avg )
if ( lmean )
for ( varID = 0; varID < nvars; varID++ )
{
if ( vlistInqVarTsteptype(vlistID1, varID) == TSTEP_CONSTANT ) continue;
......@@ -355,7 +328,7 @@ void *Runstat(void *argument)
fardiv(&vars1[0][varID][levelID], samp1[0][varID][levelID]);
}
}
else if ( operfunc == func_std || operfunc == func_var )
else if ( lvarstd )
for ( varID = 0; varID < nvars; varID++ )
{
if ( vlistInqVarTsteptype(vlistID1, varID) == TSTEP_CONSTANT ) continue;
......@@ -364,18 +337,17 @@ void *Runstat(void *argument)
{
if ( runstat_nomiss )
{
if ( operfunc == func_std )
farcstd(&vars1[0][varID][levelID], vars2[0][varID][levelID], 1.0/ndates);
if ( lstd )
farcstdx(&vars1[0][varID][levelID], vars2[0][varID][levelID], ndates, divisor);
else
farcvar(&vars1[0][varID][levelID], vars2[0][varID][levelID], 1.0/ndates);
farcvarx(&vars1[0][varID][levelID], vars2[0][varID][levelID], ndates, divisor);
}
else
{
farinv(&samp1[0][varID][levelID]);
if ( operfunc == func_std )
farstd(&vars1[0][varID][levelID], vars2[0][varID][levelID], samp1[0][varID][levelID]);
if ( lstd )
farstdx(&vars1[0][varID][levelID], vars2[0][varID][levelID], samp1[0][varID][levelID], divisor);
else
farvar(&vars1[0][varID][levelID], vars2[0][varID][levelID], samp1[0][varID][levelID]);
farvarx(&vars1[0][varID][levelID], vars2[0][varID][levelID], samp1[0][varID][levelID], divisor);
}
}
}
......@@ -410,7 +382,7 @@ void *Runstat(void *argument)
vars1[ndates] = vars1[0];
if ( !runstat_nomiss )
samp1[ndates] = samp1[0];
if ( operfunc == func_std || operfunc == func_var )
if ( lvarstd )
vars2[ndates] = vars2[0];
for ( inp = 0; inp < ndates; inp++ )
......@@ -419,7 +391,7 @@ void *Runstat(void *argument)
vars1[inp] = vars1[inp+1];
if ( !runstat_nomiss )
samp1[inp] = samp1[inp+1];
if ( operfunc == func_std || operfunc == func_var )
if ( lvarstd )
vars2[inp] = vars2[inp+1];
}
......@@ -462,7 +434,7 @@ void *Runstat(void *argument)
}
}
if ( operfunc == func_std || operfunc == func_var )
if ( lvarstd )
{
farmoq(&vars2[ndates-1][varID][levelID], vars1[ndates-1][varID][levelID]);
#if defined (_OPENMP)
......@@ -491,30 +463,15 @@ void *Runstat(void *argument)
for ( its = 0; its < ndates; its++ )
{
for ( varID = 0; varID < nvars; varID++ )
{
nlevel = zaxisInqSize(vlistInqVarZaxis(vlistID1, varID));
for ( levelID = 0; levelID < nlevel; levelID++ )
{
free(vars1[its][varID][levelID].ptr);
if ( !runstat_nomiss ) free(samp1[its][varID][levelID].ptr);
if ( operfunc == func_std || operfunc == func_var ) free(vars2[its][varID][levelID].ptr);
}
free(vars1[its][varID]);
if ( !runstat_nomiss ) free(samp1[its][varID]);
if ( operfunc == func_std || operfunc == func_var ) free(vars2[its][varID]);
}
free(vars1[its]);
if ( !runstat_nomiss ) free(samp1[its]);
if ( operfunc == func_std || operfunc == func_var ) free(vars2[its]);
field_free(vars1[its], vlistID1);
if ( !runstat_nomiss ) field_free(samp1[its], vlistID1);
if ( lvarstd ) field_free(vars2[its], vlistID1);
}
free(dtinfo);
free(vars1);
if ( !runstat_nomiss ) free(samp1);
if ( operfunc == func_std || operfunc == func_var ) free(vars2);
if ( lvarstd ) free(vars2);
if ( recVarID ) free(recVarID);
if ( recLevelID ) free(recLevelID);
......
src/field2.c
View file @ b525e00d
......@@ -452,9 +452,6 @@ void farvar(field_t *field1, field_t field2, field_t field3)
{
for ( i = 0; i < len; i++ )
{
// array3[i] = 0.111111;
if ( i < 10 )
printf("array123 %g %g %g", array1[i], array2[i], array3[i]);
if ( !DBL_IS_EQUAL(array1[i], missval1) && !DBL_IS_EQUAL(array2[i], missval2) && !DBL_IS_EQUAL(array3[i], missval3) )
{
array1[i] = array2[i]*array3[i] - (array1[i]*array3[i])*(array1[i]*array3[i]);
......@@ -462,8 +459,6 @@ void farvar(field_t *field1, field_t field2, field_t field3)
}
else
array1[i] = missval1;
if ( i < 10 )
printf(" %g\n", array1[i]);
}
}
else
......
src/modules.c
View file @ b525e00d
......@@ -375,7 +375,7 @@ void *Maggraph(void *argument);
#define RhopotOperators {"rhopot"}
#define RotuvOperators {"rotuvb"}
#define RunpctlOperators {"runpctl"}
#define RunstatOperators {"runmin", "runmax", "runsum", "runmean", "runavg", "runvar", "runstd"}
#define RunstatOperators {"runmin", "runmax", "runsum", "runmean", "runavg", "runstd", "runstd1", "runvar", "runvar1"}
#define SeascountOperators {"seascount"}
#define SeaspctlOperators {"seaspctl"}
#define SeasstatOperators {"seasmin", "seasmax", "seassum", "seasmean", "seasavg", "seasvar", "seasstd"}
......
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

Imprint | Privacy Policy