diff --git a/ChangeLog b/ChangeLog index 59a1af4ebe27bf5f8700fa64557d988a74ccc8b6..449802dd900622c8fb20a4a61c420b36e6827cc6 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,7 @@ +2025-03-15 Uwe Schulzweida + + * timpctl: Missing error handling for CDO_PCTL_NBINS>32768 (bug fix) + 2025-03-14 Uwe Schulzweida * cat, mergetime: unpack input data diff --git a/src/Timpctl.cc b/src/Timpctl.cc index 4cf13f1d4f2f54b9b59e6ccb3538cb0fcf5bfa00..8dd8fcc61d92ca482cfa46f7b54e6045cd490f02 100644 --- a/src/Timpctl.cc +++ b/src/Timpctl.cc @@ -185,7 +185,6 @@ public: { field1.init(var1); cdo_read_field(streamID1, field1); - hset.addVarLevelValues(varID, levelID, field1); } } diff --git a/src/percentiles_hist.cc b/src/percentiles_hist.cc index 8705ee1debc7343366bb07a7a7a94139915f327e..4a5f01151ef1ebb84e6d3767f18f99a0e35464db 100644 --- a/src/percentiles_hist.cc +++ b/src/percentiles_hist.cc @@ -22,7 +22,7 @@ #define SHR_PTR(p) ((unsigned short *) (p)) static int -histGetEnvNBins() +histGetEnvNumBins() { constexpr int NBINS_DEFAULT = 101; constexpr int NBINS_MINIMUM = 11; @@ -34,47 +34,47 @@ histGetEnvNBins() template <typename T> static void -hist_init_bins(int nbins, T *bins) +hist_init_bins(int numBins, T *bins) { - for (int i = 0; i < nbins; ++i) bins[i] = 0; + for (int i = 0; i < numBins; ++i) bins[i] = 0; } static void histDefBounds(HistogramEntry &hist, float a, float b) { - assert(hist.nbins > 0); + assert(hist.numBins > 0); hist.nsamp = 0; hist.min = std::min(a, b); hist.max = std::max(a, b); - hist.step = (hist.max - hist.min) / hist.nbins; + hist.step = (hist.max - hist.min) / hist.numBins; if (hist.isUint32) - hist_init_bins(hist.nbins, INT_PTR(hist.ptr)); + hist_init_bins(hist.numBins, INT_PTR(hist.ptr)); else - hist_init_bins(hist.nbins, SHR_PTR(hist.ptr)); + hist_init_bins(hist.numBins, SHR_PTR(hist.ptr)); } static inline int -calc_bin(int nbins, float histMin, float histStep, float value) +calc_bin(int numBins, float histMin, float histStep, float value) { - return (histStep > 0.0f) ? std::min((int) ((value - histMin) / histStep), nbins - 1) : 0; + return (histStep > 0.0f) ? std::min((int) ((value - histMin) / histStep), numBins - 1) : 0; } template <typename T> static void histBinAddValue(HistogramEntry &hist, T *bins, float value) { - auto bin = calc_bin(hist.nbins, hist.min, hist.step, value); - if (bin >= 0 && bin < hist.nbins) bins[bin]++; + auto bin = calc_bin(hist.numBins, hist.min, hist.step, value); + if (bin >= 0 && bin < hist.numBins) bins[bin]++; } template <typename T> static void histBinSubValue(HistogramEntry &hist, T *bins, float value) { - auto bin = calc_bin(hist.nbins, hist.min, hist.step, value); - if (bin >= 0 && bin < hist.nbins && bins[bin] > 0) bins[bin]--; + auto bin = calc_bin(hist.numBins, hist.min, hist.step, value); + if (bin >= 0 && bin < hist.numBins && bins[bin] > 0) bins[bin]--; } static void @@ -88,9 +88,9 @@ histBin(HistogramEntry &hist) for (int i = 0; i < hist.nsamp; ++i) values[i] = fltptr[i]; if (hist.isUint32) - hist_init_bins(hist.nbins, INT_PTR(hist.ptr)); + hist_init_bins(hist.numBins, INT_PTR(hist.ptr)); else - hist_init_bins(hist.nbins, SHR_PTR(hist.ptr)); + hist_init_bins(hist.numBins, SHR_PTR(hist.ptr)); if (hist.isUint32) for (int i = 0; i < hist.nsamp; ++i) histBinAddValue(hist, INT_PTR(hist.ptr), values[i]); @@ -101,7 +101,7 @@ histBin(HistogramEntry &hist) static int histReset(HistogramEntry &hist) { - assert(hist.nbins > 0); + assert(hist.numBins > 0); if (hist.nsamp < hist.capacity) { @@ -110,9 +110,9 @@ histReset(HistogramEntry &hist) else { if (hist.isUint32) - for (int i = 0; i < hist.nbins; ++i) INT_PTR(hist.ptr)[i] = 0; + for (int i = 0; i < hist.numBins; ++i) INT_PTR(hist.ptr)[i] = 0; else - for (int i = 0; i < hist.nbins; ++i) SHR_PTR(hist.ptr)[i] = 0; + for (int i = 0; i < hist.numBins; ++i) SHR_PTR(hist.ptr)[i] = 0; } hist.nsamp = 0; @@ -132,7 +132,7 @@ histCheckValue(const HistogramEntry &hist, float &value) static int histAddValue(HistogramEntry &hist, float value) { - assert(hist.nbins > 0); + assert(hist.numBins > 0); // if (is_equal(hist.min, hist.max)) return 0; @@ -178,7 +178,7 @@ histRemoveValue(HistogramEntry &hist, float value) static int histSubValue(HistogramEntry &hist, float value) { - assert(hist.nbins > 0); + assert(hist.numBins > 0); if (is_equal(hist.min, hist.max)) return 0; @@ -211,7 +211,7 @@ histGetBin(int nbins, double s, const T *ptr) while (count < s); assert(i > 0); - assert(i - 1 < nbins); + assert((i - 1) < nbins); assert(ptr[i - 1] > 0); double t = (count - s) / ptr[i - 1]; @@ -223,7 +223,7 @@ static double histGetPercentile(const HistogramEntry &hist, double p) { assert(hist.nsamp > 0); - assert(hist.nbins > 0); + assert(hist.numBins > 0); assert(p >= 0.0); assert(p <= 100.0); @@ -233,12 +233,12 @@ histGetPercentile(const HistogramEntry &hist, double p) if (lprint && Options::cdoVerbose) { lprint = false; - cdo_print("Using percentile method: histogram with %d bins", hist.nbins); + cdo_print("Using percentile method: histogram with %d bins", hist.numBins); } double s = hist.nsamp * (p / 100.0); - auto bin = hist.isUint32 ? histGetBin(hist.nbins, s, INT_PTR(hist.ptr)) : histGetBin(hist.nbins, s, SHR_PTR(hist.ptr)); + auto bin = hist.isUint32 ? histGetBin(hist.numBins, s, INT_PTR(hist.ptr)) : histGetBin(hist.numBins, s, SHR_PTR(hist.ptr)); // assert(hist.step > 0.0f); @@ -248,37 +248,37 @@ histGetPercentile(const HistogramEntry &hist, double p) } void -HistogramSet::createVarLevels(int varID, int nlevels, size_t nhists) +HistogramSet::createVarLevels(int varID, int numLevels, size_t numHists) { - auto nbins = histGetEnvNBins(); + auto numBins = histGetEnvNumBins(); - assert(nbins > 0); - assert(nlevels > 0); - assert(nhists > 0); + assert(numBins > 0); + assert(numLevels > 0); + assert(numHists > 0); - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - this->var_nlevels[varID] = nlevels; - this->var_nhists[varID] = nhists; - this->histograms[varID].resize(nlevels); + this->var_numLevels[varID] = numLevels; + this->var_numHists[varID] = numHists; + this->histograms[varID].resize(numLevels); - for (int levelID = 0; levelID < nlevels; ++levelID) + for (int levelID = 0; levelID < numLevels; ++levelID) { - this->histograms[varID][levelID].resize(nhists); + this->histograms[varID][levelID].resize(numHists); auto &hists = this->histograms[varID][levelID]; - for (size_t histID = 0; histID < nhists; histID++) + for (size_t histID = 0; histID < numHists; histID++) { hists[histID].min = 0.0f; hists[histID].max = 0.0f; hists[histID].step = 0.0f; - hists[histID].nbins = nbins; + hists[histID].numBins = numBins; hists[histID].nsamp = 0; - auto isUint32 = (this->nsteps >= USHRT_MAX); + auto isUint32 = (this->numSteps >= USHRT_MAX); auto vsize = isUint32 ? sizeof(unsigned int) : sizeof(unsigned short); hists[histID].isUint32 = isUint32; - hists[histID].capacity = FLT_CAPACITY(nbins, vsize); - hists[histID].ptr = std::malloc(nbins * vsize); + hists[histID].capacity = FLT_CAPACITY(numBins, vsize); + hists[histID].ptr = std::malloc(numBins * vsize); if (hists[histID].ptr == nullptr) cdo_abort("Not enough memory (%s)", __func__); } } @@ -286,12 +286,12 @@ HistogramSet::createVarLevels(int varID, int nlevels, size_t nhists) template <typename T1, typename T2> static void -def_bounds(size_t nhists, std::vector<HistogramEntry> &hists, const Varray<T1> &v1, const Varray<T2> &v2, float mv1, float mv2) +def_bounds(size_t numHists, std::vector<HistogramEntry> &hists, Varray<T1> const &v1, Varray<T2> const &v2, float mv1, float mv2) { assert(!v1.empty()); assert(!v2.empty()); - for (size_t i = 0; i < nhists; ++i) + for (size_t i = 0; i < numHists; ++i) { float a = v1[i]; float b = v2[i]; @@ -306,26 +306,26 @@ def_bounds(size_t nhists, std::vector<HistogramEntry> &hists, const Varray<T1> & void HistogramSet::defVarLevelBounds(int varID, int levelID, Field const &field1, Field const &field2) { - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - auto nlevels = this->var_nlevels[varID]; - if (levelID < 0 || levelID >= nlevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); + auto numLevels = this->var_numLevels[varID]; + if (levelID < 0 || levelID >= numLevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); - auto nhists = this->var_nhists[varID]; - if (nhists != gridInqSize(field1.grid) || nhists != gridInqSize(field2.grid)) cdo_abort("Grids are different (%s)", __func__); + auto numHists = this->var_numHists[varID]; + if (numHists != gridInqSize(field1.grid) || numHists != gridInqSize(field2.grid)) cdo_abort("Grids are different (%s)", __func__); auto &hists = this->histograms[varID][levelID]; float mv1 = (float) field1.missval; float mv2 = (float) field2.missval; - auto func = [&](auto const &v1, auto const &v2) { def_bounds(nhists, hists, v1, v2, mv1, mv2); }; + auto func = [&](auto const &v1, auto const &v2) { def_bounds(numHists, hists, v1, v2, mv1, mv2); }; field_operation2(func, field1, field2); } template <typename T> static int -histAddVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray<T> const &v, size_t numMissVals, double missval) +histAddVarLevelValues(size_t numHists, std::vector<HistogramEntry> &hists, Varray<T> const &v, size_t numMissVals, double missval) { T mv = missval; assert(!v.empty()); @@ -333,12 +333,12 @@ histAddVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray< int nign = 0; if (numMissVals) { - for (size_t i = 0; i < nhists; ++i) + for (size_t i = 0; i < numHists; ++i) if (fp_is_not_equal(v[i], mv)) nign += histAddValue(hists[i], v[i]); } else { - for (size_t i = 0; i < nhists; ++i) nign += histAddValue(hists[i], v[i]); + for (size_t i = 0; i < numHists; ++i) nign += histAddValue(hists[i], v[i]); } return nign; @@ -346,7 +346,7 @@ histAddVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray< template <typename T> static int -histSubVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray<T> const &v, size_t numMissVals, double missval) +histSubVarLevelValues(size_t numHists, std::vector<HistogramEntry> &hists, Varray<T> const &v, size_t numMissVals, double missval) { T mv = missval; assert(!v.empty()); @@ -354,12 +354,12 @@ histSubVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray< int nign = 0; if (numMissVals) { - for (size_t i = 0; i < nhists; ++i) + for (size_t i = 0; i < numHists; ++i) if (fp_is_not_equal(v[i], mv)) nign += histSubValue(hists[i], v[i]); } else { - for (size_t i = 0; i < nhists; ++i) nign += histSubValue(hists[i], v[i]); + for (size_t i = 0; i < numHists; ++i) nign += histSubValue(hists[i], v[i]); } return nign; @@ -368,21 +368,21 @@ histSubVarLevelValues(size_t nhists, std::vector<HistogramEntry> &hists, Varray< int HistogramSet::addVarLevelValues(int varID, int levelID, Field const &field) { - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - auto nlevels = this->var_nlevels[varID]; - if (levelID < 0 || levelID >= nlevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); + auto numLevels = this->var_numLevels[varID]; + if (levelID < 0 || levelID >= numLevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); - auto nhists = this->var_nhists[varID]; - if (nhists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); + auto numHists = this->var_numHists[varID]; + if (numHists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); auto &hists = this->histograms[varID][levelID]; - auto func = [&](auto &v, auto numMissVals, double mv) { return histAddVarLevelValues(nhists, hists, v, numMissVals, mv); }; + auto func = [&](auto &v, auto numMissVals, double mv) { return histAddVarLevelValues(numHists, hists, v, numMissVals, mv); }; auto nign = field_operation(func, field, field.numMissVals, field.missval); if (nign) { - cdo_warning("%d out of %d grid values are out of bounds and have been ignored (%s)", nign, nhists, __func__); + cdo_warning("%d out of %d grid values are out of bounds and have been ignored (%s)", nign, numHists, __func__); return 1; } @@ -392,21 +392,21 @@ HistogramSet::addVarLevelValues(int varID, int levelID, Field const &field) int HistogramSet::subVarLevelValues(int varID, int levelID, Field const &field) { - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - auto nlevels = this->var_nlevels[varID]; - if (levelID < 0 || levelID >= nlevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); + auto numLevels = this->var_numLevels[varID]; + if (levelID < 0 || levelID >= numLevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); - auto nhists = this->var_nhists[varID]; - if (nhists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); + auto numHists = this->var_numHists[varID]; + if (numHists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); auto &hists = this->histograms[varID][levelID]; - auto func = [&](auto &v, auto numMissVals, double mv) { return histSubVarLevelValues(nhists, hists, v, numMissVals, mv); }; + auto func = [&](auto &v, auto numMissVals, double mv) { return histSubVarLevelValues(numHists, hists, v, numMissVals, mv); }; auto nign = field_operation(func, field, field.numMissVals, field.missval); if (nign) { - cdo_warning("%d out of %d grid values are out of bounds and have been ignored (%s)", nign, nhists, __func__); + cdo_warning("%d out of %d grid values are out of bounds and have been ignored (%s)", nign, numHists, __func__); return 1; } @@ -416,33 +416,33 @@ HistogramSet::subVarLevelValues(int varID, int levelID, Field const &field) void HistogramSet::reset(int varID, int levelID) { - assert(nvars > 0); + assert(numVars > 0); - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - auto nlevels = this->var_nlevels[varID]; - assert(nlevels > 0); + auto numLevels = this->var_numLevels[varID]; + assert(numLevels > 0); - if (levelID < 0 || levelID >= nlevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); + if (levelID < 0 || levelID >= numLevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); - auto nhists = this->var_nhists[varID]; - assert(nhists > 0); + auto numHists = this->var_numHists[varID]; + assert(numHists > 0); auto &hists = this->histograms[varID][levelID]; - for (size_t i = 0; i < nhists; ++i) histReset(hists[i]); + for (size_t i = 0; i < numHists; ++i) histReset(hists[i]); } */ template <typename T> static size_t -calcPercentile(size_t nhists, const std::vector<HistogramEntry> &hists, double p, Varray<T> &v, double missval) +calcPercentile(size_t numHists, const std::vector<HistogramEntry> &hists, double p, Varray<T> &v, double missval) { T mv = missval; assert(!v.empty()); size_t numMissVals = 0; - for (size_t i = 0; i < nhists; ++i) + for (size_t i = 0; i < numHists; ++i) { if (hists[i].nsamp) { v[i] = histGetPercentile(hists[i], p); } else @@ -458,16 +458,16 @@ calcPercentile(size_t nhists, const std::vector<HistogramEntry> &hists, double p void HistogramSet::getVarLevelPercentiles(Field &field, int varID, int levelID, double p) { - if (varID < 0 || varID >= nvars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); + if (varID < 0 || varID >= numVars) cdo_abort("Illegal argument: varID %d is undefined (%s)", varID, __func__); - auto nlevels = this->var_nlevels[varID]; - if (levelID < 0 || levelID >= nlevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); + auto numLevels = this->var_numLevels[varID]; + if (levelID < 0 || levelID >= numLevels) cdo_abort("Illegal argument: levelID %d is undefined (%s)", levelID, __func__); - auto nhists = this->var_nhists[varID]; - if (nhists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); + auto numHists = this->var_numHists[varID]; + if (numHists != gridInqSize(field.grid)) cdo_abort("Grids are different (%s)", __func__); auto const &hists = this->histograms[varID][levelID]; - auto func = [&](auto &v, double mv) { return calcPercentile(nhists, hists, p, v, mv); }; + auto func = [&](auto &v, double mv) { return calcPercentile(numHists, hists, p, v, mv); }; field.numMissVals = field_operation(func, field, field.missval); } diff --git a/src/percentiles_hist.h b/src/percentiles_hist.h index 4f9a444c5af13a960ed67287528e6c919e38dff3..14dc4802526b77773141bebd70144421e945800a 100644 --- a/src/percentiles_hist.h +++ b/src/percentiles_hist.h @@ -14,7 +14,7 @@ struct HistogramEntry float step = 0.0f; int nsamp = 0; int capacity = 0; - short nbins = 0; + int numBins = 0; bool isUint32 = false; }; @@ -27,57 +27,57 @@ HistogramSet // clang-format on { private: - int nvars = 0; - int nsteps = 0; - std::vector<int> var_nlevels; - std::vector<size_t> var_nhists; + int numVars = 0; + int numSteps = 0; + std::vector<int> var_numLevels; + std::vector<size_t> var_numHists; std::vector<std::vector<std::vector<HistogramEntry>>> histograms; void init() { - var_nlevels.resize(nvars, 0); - var_nhists.resize(nvars, 0); - histograms.resize(nvars); + var_numLevels.resize(numVars, 0); + var_numHists.resize(numVars, 0); + histograms.resize(numVars); } public: HistogramSet() {} - explicit HistogramSet(int _nvars) : nvars(_nvars) + explicit HistogramSet(int _numVars) : numVars(_numVars) { - assert(nvars > 0); + assert(numVars > 0); init(); } - HistogramSet(int _nvars, int _nsteps) : nvars(_nvars), nsteps(_nsteps) + HistogramSet(int _numVars, int _numSteps) : numVars(_numVars), numSteps(_numSteps) { - assert(nvars > 0); + assert(numVars > 0); init(); } void - create(int _nvars, int _nsteps = 0) + create(int _numVars, int _numSteps = 0) { - nvars = _nvars; - nsteps = _nsteps; - assert(nvars > 0); + numVars = _numVars; + numSteps = _numSteps; + assert(numVars > 0); init(); } ~HistogramSet() { - for (auto varID = nvars; varID-- > 0;) + for (auto varID = numVars; varID-- > 0;) { - auto nhists = this->var_nhists[varID]; - for (auto levelID = this->var_nlevels[varID]; levelID-- > 0;) + auto numHists = this->var_numHists[varID]; + for (auto levelID = this->var_numLevels[varID]; levelID-- > 0;) { - for (auto histID = nhists; histID-- > 0;) std::free(this->histograms[varID][levelID][histID].ptr); + for (auto histID = numHists; histID-- > 0;) std::free(this->histograms[varID][levelID][histID].ptr); } } } - void createVarLevels(int varID, int nlevels, size_t nhists); + void createVarLevels(int varID, int numLevels, size_t numHists); void defVarLevelBounds(int varID, int levelID, Field const &field1, Field const &field2); int addVarLevelValues(int varID, int levelID, Field const &field); int subVarLevelValues(int varID, int levelID, Field const &field);