NumCpp  2.7.0
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 
34 #include "NumCpp/Core/Types.hpp"
35 #include "NumCpp/NdArray.hpp"
37 
38 #include <cmath>
39 #include <string>
40 
41 namespace nc
42 {
43  namespace imageProcessing
44  {
45  //============================================================================
46  // Method Description:
55  template<typename dtype>
56  dtype generateThreshold(const NdArray<dtype>& inImageArray, double inRate)
57  {
59 
60  if (inRate < 0.0 || inRate > 1.0)
61  {
62  THROW_INVALID_ARGUMENT_ERROR("input rate must be of the range [0, 1]");
63  }
64 
65  // first build a histogram
66  auto minValue = static_cast<int32>(std::floor(inImageArray.min().item()));
67  auto maxValue = static_cast<int32>(std::floor(inImageArray.max().item()));
68 
69  if (utils::essentiallyEqual(inRate, 0.0))
70  {
71  return static_cast<dtype>(maxValue);
72  }
73 
74  if (utils::essentiallyEqual(inRate, 1.0))
75  {
77  {
78  return static_cast<dtype>(minValue - 1);
79  }
80 
81  return dtype{ 0 };
82  }
83 
84  const auto histSize = static_cast<uint32>(maxValue - minValue + 1);
85 
86  NdArray<double> histogram(1, histSize);
87  histogram.zeros();
88  for (auto intensity : inImageArray)
89  {
90  const auto bin = static_cast<uint32>(static_cast<int32>(std::floor(intensity)) - minValue);
91  ++histogram[bin];
92  }
93 
94  // integrate the normalized histogram from right to left to make a survival function (1 - CDF)
95  const auto dNumPixels = static_cast<double>(inImageArray.size());
96  NdArray<double> survivalFunction(1, histSize + 1);
97  survivalFunction[-1] = 0.0;
98  for (int32 i = histSize - 1; i > -1; --i)
99  {
100  double histValue = histogram[i] / dNumPixels;
101  survivalFunction[i] = survivalFunction[i + 1] + histValue;
102  }
103 
104  // binary search through the survival function to find the rate
105  uint32 indexLow = 0;
106  uint32 indexHigh = histSize - 1;
107  uint32 index = indexHigh / 2; // integer division
108 
109  constexpr bool keepGoing = true;
110  while (keepGoing)
111  {
112  const double value = survivalFunction[index];
113  if (value < inRate)
114  {
115  indexHigh = index;
116  }
117  else if (value > inRate)
118  {
119  indexLow = index;
120  }
121  else
122  {
123  const int32 thresh = static_cast<int32>(index) + minValue - 1;
125  {
126  return static_cast<dtype>(thresh);
127  }
128 
129  return thresh < 0 ? 0 : static_cast<dtype>(thresh);
130  }
131 
132  if (indexHigh - indexLow < 2)
133  {
134  return static_cast<dtype>(static_cast<int32>(indexHigh) + minValue - 1);
135  }
136 
137  index = indexLow + (indexHigh - indexLow) / 2;
138  }
139 
140  // shouldn't ever get here but stop the compiler from throwing a warning
141  return static_cast<dtype>(histSize - 1);
142  }
143 
144  } // namespace imageProcessing
145 } // namespace nc
#define THROW_INVALID_ARGUMENT_ERROR(msg)
Definition: Error.hpp:36
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:37
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:72
size_type size() const noexcept
Definition: NdArrayCore.hpp:4296
NdArray< dtype > min(Axis inAxis=Axis::NONE) const
Definition: NdArrayCore.hpp:3005
NdArray< dtype > max(Axis inAxis=Axis::NONE) const
Definition: NdArrayCore.hpp:2950
dtype generateThreshold(const NdArray< dtype > &inImageArray, double inRate)
Definition: generateThreshold.hpp:56
bool essentiallyEqual(dtype inValue1, dtype inValue2) noexcept
Definition: essentiallyEqual.hpp:52
Definition: Coordinate.hpp:45
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