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

Renamed NDAY to MAX_DAYS.

parent bb837c81
Hide whitespace changes
Inline Side-by-side
src/Ydaypctl.cc
View file @ 192dd9f0
......@@ -23,17 +23,16 @@
#include <cdi.h>
#include "cdo_options.h"
#include "functs.h"
#include "process_int.h"
#include "cdo_vlist.h"
#include "param_conversion.h"
#include "percentiles_hist.h"
#include "percentiles.h"
#define NDAY 373
constexpr int MAX_DAYS = 373;
int getmonthday(int64_t date);
......@@ -47,20 +46,20 @@ Ydaypctl(void *process)
int levelID;
size_t nmiss;
int nlevels;
int64_t vdates1[NDAY], vdates2[NDAY];
int vtimes1[NDAY];
long nsets[NDAY];
FieldVector2D vars1[NDAY];
HISTOGRAM_SET *hsets[NDAY];
int64_t vdates1[MAX_DAYS], vdates2[MAX_DAYS];
int vtimes1[MAX_DAYS];
long nsets[MAX_DAYS];
FieldVector2D vars1[MAX_DAYS];
HISTOGRAM_SET *hsets[MAX_DAYS];
cdoInitialize(process);
cdoOperatorAdd("ydaypctl", func_pctl, 0, nullptr);
operatorInputArg("percentile number");
const double pn = parameter2double(operatorArgv()[0]);
const auto pn = parameter2double(operatorArgv()[0]);
percentile_check_number(pn);
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
{
hsets[dayoy] = nullptr;
nsets[dayoy] = 0;
......@@ -70,17 +69,17 @@ Ydaypctl(void *process)
const auto streamID2 = cdoOpenRead(1);
const auto streamID3 = cdoOpenRead(2);
const int vlistID1 = cdoStreamInqVlist(streamID1);
const int vlistID2 = cdoStreamInqVlist(streamID2);
const int vlistID3 = cdoStreamInqVlist(streamID3);
const int vlistID4 = vlistDuplicate(vlistID1);
const auto vlistID1 = cdoStreamInqVlist(streamID1);
const auto vlistID2 = cdoStreamInqVlist(streamID2);
const auto vlistID3 = cdoStreamInqVlist(streamID3);
const auto vlistID4 = vlistDuplicate(vlistID1);
vlistCompare(vlistID1, vlistID2, CMP_ALL);
vlistCompare(vlistID1, vlistID3, CMP_ALL);
const int taxisID1 = vlistInqTaxis(vlistID1);
const int taxisID2 = vlistInqTaxis(vlistID2);
const int taxisID3 = vlistInqTaxis(vlistID3);
const auto taxisID1 = vlistInqTaxis(vlistID1);
const auto taxisID2 = vlistInqTaxis(vlistID2);
const auto taxisID3 = vlistInqTaxis(vlistID3);
/* TODO - check that time axes 2 and 3 are equal */
const int taxisID4 = taxisDuplicate(taxisID1);
......@@ -90,9 +89,9 @@ Ydaypctl(void *process)
const auto streamID4 = cdoOpenWrite(3);
cdoDefVlist(streamID4, vlistID4);
const int nvars = vlistNvars(vlistID1);
const auto nvars = vlistNvars(vlistID1);
const int maxrecs = vlistNrecs(vlistID1);
const auto maxrecs = vlistNrecs(vlistID1);
std::vector<RecordInfo> recList(maxrecs);
const auto gridsizemax = vlistGridsizeMax(vlistID1);
......@@ -107,7 +106,7 @@ Ydaypctl(void *process)
cdoAbort("Number of records at time step %d of %s and %s differ!", tsID + 1, cdoGetStreamName(1), cdoGetStreamName(2));
const auto vdate = taxisInqVdate(taxisID2);
const int vtime = taxisInqVtime(taxisID2);
const auto vtime = taxisInqVtime(taxisID2);
if (vdate != taxisInqVdate(taxisID3))
cdoAbort("Verification dates at time step %d of %s and %s differ!", tsID + 1, cdoGetStreamName(1), cdoGetStreamName(2));
......@@ -117,7 +116,7 @@ Ydaypctl(void *process)
cdiDecodeDate(vdate, &year, &month, &day);
const int dayoy = (month >= 1 && month <= 12) ? (month - 1) * 31 + day : 0;
if (dayoy < 0 || dayoy >= NDAY) cdoAbort("Day %d out of range!", dayoy);
if (dayoy < 0 || dayoy >= MAX_DAYS) cdoAbort("Day %d out of range!", dayoy);
vdates2[dayoy] = vdate;
......@@ -159,14 +158,14 @@ Ydaypctl(void *process)
while ((nrecs = cdoStreamInqTimestep(streamID1, tsID)))
{
const auto vdate = taxisInqVdate(taxisID1);
const int vtime = taxisInqVtime(taxisID1);
const auto vtime = taxisInqVtime(taxisID1);
if (Options::cdoVerbose) cdoPrint("process timestep: %d %d %d", tsID + 1, vdate, vtime);
cdiDecodeDate(vdate, &year, &month, &day);
const int dayoy = (month >= 1 && month <= 12) ? (month - 1) * 31 + day : 0;
if (dayoy < 0 || dayoy >= NDAY) cdoAbort("Day %d out of range!", dayoy);
if (dayoy < 0 || dayoy >= MAX_DAYS) cdoAbort("Day %d out of range!", dayoy);
vdates1[dayoy] = vdate;
vtimes1[dayoy] = vtime;
......@@ -195,7 +194,7 @@ Ydaypctl(void *process)
}
int otsID = 0;
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
if (nsets[dayoy])
{
if (getmonthday(vdates1[dayoy]) != getmonthday(vdates2[dayoy]))
......@@ -219,8 +218,8 @@ Ydaypctl(void *process)
{
if (otsID && recList[recID].lconst) continue;
const int varID = recList[recID].varID;
const int levelID = recList[recID].levelID;
const auto varID = recList[recID].varID;
const auto levelID = recList[recID].levelID;
cdoDefRecord(streamID4, varID, levelID);
cdoWriteRecord(streamID4, vars1[dayoy][varID][levelID].vec.data(), vars1[dayoy][varID][levelID].nmiss);
}
......@@ -228,12 +227,9 @@ Ydaypctl(void *process)
otsID++;
}
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
{
if (vars1[dayoy].size())
{
hsetDestroy(hsets[dayoy]);
}
if (vars1[dayoy].size()) hsetDestroy(hsets[dayoy]);
}
cdoStreamClose(streamID4);
......
src/Ydrunstat.cc
View file @ 192dd9f0
......@@ -31,29 +31,28 @@
#include <cdi.h>
#include "dmemory.h"
#include "functs.h"
#include "process_int.h"
#include "param_conversion.h"
#include "calendar.h"
#include "datetime.h"
#define NDAY 373
constexpr int MAX_DAYS = 373;
class YDAY_STATS
{
public:
int64_t vdate[NDAY];
int vtime[NDAY];
FieldVector2D vars1[NDAY];
FieldVector2D vars2[NDAY];
int nsets[NDAY];
int64_t vdate[MAX_DAYS];
int vtime[MAX_DAYS];
FieldVector2D vars1[MAX_DAYS];
FieldVector2D vars2[MAX_DAYS];
int nsets[MAX_DAYS];
int vlist;
YDAY_STATS(int vlistID) : vlist(vlistID)
{
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
{
vdate[dayoy] = 0;
vtime[dayoy] = 0;
......@@ -68,13 +67,13 @@ ydstatUpdate(YDAY_STATS &stats, int64_t vdate, int vtime, const FieldVector2D &v
{
bool lvarstd = vars2.size();
auto nvars = vlistNvars(stats.vlist);
const auto nvars = vlistNvars(stats.vlist);
int year, month, day;
cdiDecodeDate(vdate, &year, &month, &day);
const int dayoy = (month >= 1 && month <= 12) ? (month - 1) * 31 + day : 0;
if (dayoy < 0 || dayoy >= NDAY) cdoAbort("day %d out of range!", dayoy);
if (dayoy < 0 || dayoy >= MAX_DAYS) cdoAbort("day %d out of range!", dayoy);
stats.vdate[dayoy] = vdate;
stats.vtime[dayoy] = vtime;
......@@ -129,7 +128,7 @@ ydstatFinalize(YDAY_STATS &stats, int operfunc)
const auto nvars = vlistNvars(stats.vlist);
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
if (stats.nsets[dayoy])
{
switch (operfunc)
......@@ -170,17 +169,9 @@ ydstatFinalize(YDAY_STATS &stats, int operfunc)
}
}
void *
Ydrunstat(void *process)
static void
addOperators(void)
{
int varID;
int nrecs;
int levelID;
int tsID;
int inp, its;
cdoInitialize(process);
// clang-format off
cdoOperatorAdd("ydrunmin", func_min, 0, nullptr);
cdoOperatorAdd("ydrunmax", func_max, 0, nullptr);
......@@ -192,6 +183,16 @@ Ydrunstat(void *process)
cdoOperatorAdd("ydrunstd", func_std, 0, nullptr);
cdoOperatorAdd("ydrunstd1", func_std1, 0, nullptr);
// clang-format on
}
void *
Ydrunstat(void *process)
{
int varID, levelID;
cdoInitialize(process);
addOperators();
const auto operatorID = cdoOperatorID();
const auto operfunc = cdoOperatorF1(operatorID);
......@@ -211,8 +212,7 @@ Ydrunstat(void *process)
if (taxisHasBounds(taxisID2)) taxisDeleteBounds(taxisID2);
vlistDefTaxis(vlistID2, taxisID2);
const auto calendar = taxisInqCalendar(taxisID1);
const auto dpy = calendar_dpy(calendar);
const auto dpy = calendar_dpy(taxisInqCalendar(taxisID1));
const auto streamID2 = cdoOpenWrite(1);
cdoDefVlist(streamID2, vlistID2);
......@@ -225,15 +225,16 @@ Ydrunstat(void *process)
YDAY_STATS stats(vlistID1);
FieldVector3D vars1(ndates + 1), vars2(ndates + 1);
for (its = 0; its < ndates; its++)
for (int its = 0; its < ndates; its++)
{
fieldsFromVlist(vlistID1, vars1[its], FIELD_VEC);
if (lvarstd) fieldsFromVlist(vlistID1, vars2[its], FIELD_VEC);
}
int tsID = 0;
for (tsID = 0; tsID < ndates; tsID++)
{
nrecs = cdoStreamInqTimestep(streamID1, tsID);
const auto nrecs = cdoStreamInqTimestep(streamID1, tsID);
if (nrecs == 0) cdoAbort("File has less then %d timesteps!", ndates);
datetime[tsID].date = taxisInqVdate(taxisID1);
......@@ -256,7 +257,7 @@ Ydrunstat(void *process)
if (lvarstd)
{
vfarmoq(vars2[tsID][varID][levelID], vars1[tsID][varID][levelID]);
vfarmoq(vars2[tsID][varID][levelID], rvars1);
for (int inp = 0; inp < tsID; inp++)
{
vfarsumq(vars2[inp][varID][levelID], rvars1);
......@@ -286,14 +287,14 @@ Ydrunstat(void *process)
vars1[ndates] = vars1[0];
if (lvarstd) vars2[ndates] = vars2[0];
for (inp = 0; inp < ndates; inp++)
for (int inp = 0; inp < ndates; inp++)
{
datetime[inp] = datetime[inp + 1];
vars1[inp] = vars1[inp + 1];
if (lvarstd) vars2[inp] = vars2[inp + 1];
}
nrecs = cdoStreamInqTimestep(streamID1, tsID);
const auto nrecs = cdoStreamInqTimestep(streamID1, tsID);
if (nrecs == 0) break;
datetime[ndates - 1].date = taxisInqVdate(taxisID1);
......@@ -309,16 +310,16 @@ Ydrunstat(void *process)
if (lvarstd)
{
for (inp = 0; inp < ndates - 1; inp++)
vfarmoq(vars2[ndates - 1][varID][levelID], rvars1);
for (int inp = 0; inp < ndates - 1; inp++)
{
vfarsumq(vars2[inp][varID][levelID], rvars1);
vfarsum(vars1[inp][varID][levelID], rvars1);
}
vfarmoq(vars2[ndates - 1][varID][levelID], rvars1);
}
else
{
for (inp = 0; inp < ndates - 1; inp++)
for (int inp = 0; inp < ndates - 1; inp++)
{
vfarfun(vars1[inp][varID][levelID], rvars1, operfunc);
}
......@@ -328,31 +329,11 @@ Ydrunstat(void *process)
tsID++;
}
/*
// set the year to the minimum of years found on output timestep
int outyear = 1e9;
int year, month, day;
for ( int dayoy = 0; dayoy < NDAY; dayoy++ )
if ( stats.nsets[dayoy] )
{
cdiDecodeDate(stats.vdate[dayoy], &year, &month, &day);
if ( year < outyear ) outyear = year;
}
for ( int dayoy = 0; dayoy < NDAY; dayoy++ )
if ( stats.nsets[dayoy] )
{
cdiDecodeDate(stats.vdate[dayoy], &year, &month, &day);
// printf("vdates[%d] = %d nsets = %d\n", dayoy, stats.vdate[dayoy],
stats.nsets[dayoy]); if ( year > outyear ) stats.vdate[dayoy] =
cdiEncodeDate(outyear, month, day);
}
*/
ydstatFinalize(stats, operfunc);
int otsID = 0;
for (int dayoy = 0; dayoy < NDAY; dayoy++)
for (int dayoy = 0; dayoy < MAX_DAYS; dayoy++)
if (stats.nsets[dayoy])
{
taxisDefVdate(taxisID2, stats.vdate[dayoy]);
......@@ -363,8 +344,8 @@ Ydrunstat(void *process)
{
if (otsID && recList[recID].lconst) continue;
auto varID = recList[recID].varID;
auto levelID = recList[recID].levelID;
const auto varID = recList[recID].varID;
const auto levelID = recList[recID].levelID;
Field &rvars1 = stats.vars1[dayoy][varID][levelID];
cdoDefRecord(streamID2, varID, levelID);
......
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