NumCpp  2.10.1
A Templatized Header Only C++ Implementation of the Python NumPy Library
generateThreshold.hpp
Go to the documentation of this file.
1 
29 #pragma once
30 
31 #include <cmath>
32 #include <string>
33 
37 #include "NumCpp/Core/Types.hpp"
38 #include "NumCpp/NdArray.hpp"
40 
41 namespace nc::imageProcessing
42 {
43  //============================================================================
44  // Method Description:
53  template<typename dtype>
54  dtype generateThreshold(const NdArray<dtype>& inImageArray, double inRate)
55  {
57 
58  if (inRate < 0. || inRate > 1.)
59  {
60  THROW_INVALID_ARGUMENT_ERROR("input rate must be of the range [0, 1]");
61  }
62 
63  // first build a histogram
64  auto minValue = static_cast<int32>(std::floor(inImageArray.min().item()));
65  auto maxValue = static_cast<int32>(std::floor(inImageArray.max().item()));
66 
67  if (utils::essentiallyEqual(inRate, 0.))
68  {
69  return static_cast<dtype>(maxValue);
70  }
71 
72  if (utils::essentiallyEqual(inRate, 1.))
73  {
75  {
76  return static_cast<dtype>(minValue - 1);
77  }
78 
79  return dtype{ 0 };
80  }
81 
82  const auto histSize = static_cast<uint32>(maxValue - minValue + 1);
83 
84  NdArray<double> histogram(1, histSize);
85  histogram.zeros();
86  for (auto intensity : inImageArray)
87  {
88  const auto bin = static_cast<uint32>(static_cast<int32>(std::floor(intensity)) - minValue);
89  ++histogram[bin];
90  }
91 
92  // integrate the normalized histogram from right to left to make a survival function (1 - CDF)
93  const auto dNumPixels = static_cast<double>(inImageArray.size());
94  NdArray<double> survivalFunction(1, histSize + 1);
95  survivalFunction[-1] = 0.;
96  for (int32 i = histSize - 1; i > -1; --i)
97  {
98  double histValue = histogram[i] / dNumPixels;
99  survivalFunction[i] = survivalFunction[i + 1] + histValue;
100  }
101 
102  // binary search through the survival function to find the rate
103  uint32 indexLow = 0;
104  uint32 indexHigh = histSize - 1;
105  uint32 index = indexHigh / 2; // integer division
106 
107  constexpr bool keepGoing = true;
108  while (keepGoing)
109  {
110  const double value = survivalFunction[index];
111  if (value < inRate)
112  {
113  indexHigh = index;
114  }
115  else if (value > inRate)
116  {
117  indexLow = index;
118  }
119  else
120  {
121  const int32 thresh = static_cast<int32>(index) + minValue - 1;
123  {
124  return static_cast<dtype>(thresh);
125  }
126 
127  return thresh < 0 ? 0 : static_cast<dtype>(thresh);
128  }
129 
130  if (indexHigh - indexLow < 2)
131  {
132  return static_cast<dtype>(static_cast<int32>(indexHigh) + minValue - 1);
133  }
134 
135  index = indexLow + (indexHigh - indexLow) / 2;
136  }
137 
138  // shouldn't ever get here but stop the compiler from throwing a warning
139  return static_cast<dtype>(histSize - 1);
140  }
141 
142 } // namespace nc::imageProcessing
#define THROW_INVALID_ARGUMENT_ERROR(msg)
Definition: Error.hpp:37
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:39
Holds info about the dtype.
Definition: DtypeInfo.hpp:41
Holds 1D and 2D arrays, the main work horse of the NumCpp library.
Definition: NdArrayCore.hpp:138
self_type max(Axis inAxis=Axis::NONE) const
Definition: NdArrayCore.hpp:2964
size_type size() const noexcept
Definition: NdArrayCore.hpp:4415
value_type item() const
Definition: NdArrayCore.hpp:2945
self_type min(Axis inAxis=Axis::NONE) const
Definition: NdArrayCore.hpp:3008
Definition: applyThreshold.hpp:34
dtype generateThreshold(const NdArray< dtype > &inImageArray, double inRate)
Definition: generateThreshold.hpp:54
bool essentiallyEqual(dtype inValue1, dtype inValue2) noexcept
Definition: essentiallyEqual.hpp:48
dtype floor(dtype inValue) noexcept
Definition: floor.hpp:48
std::int32_t int32
Definition: Types.hpp:36
NdArray< uint32 > histogram(const NdArray< dtype > &inArray, const NdArray< double > &inBinEdges)
Definition: histogram.hpp:57
std::uint32_t uint32
Definition: Types.hpp:40