ITK  6.0.0
Insight Toolkit
itkConnectedImageNeighborhoodShape.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * https://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18
19#ifndef itkConnectedImageNeighborhoodShape_h
20#define itkConnectedImageNeighborhoodShape_h
21
22#include "itkIntTypes.h" // For uintmax_t
23#include "itkMath.h"
24#include "itkOffset.h"
25
26#include <array>
27#include <cassert>
28#include <limits>
29
30namespace itk
31{
32
86template <unsigned int VImageDimension>
88{
89public:
90 static constexpr unsigned int ImageDimension = VImageDimension;
104 constexpr explicit ConnectedImageNeighborhoodShape(const size_t maximumCityblockDistance,
105 const bool includeCenterPixel) noexcept
106 : m_MaximumCityblockDistance{ maximumCityblockDistance }
107 , m_IncludeCenterPixel{ includeCenterPixel }
108 , m_NumberOfOffsets{ CalculateNumberOfOffsets(maximumCityblockDistance, includeCenterPixel) }
109 {}
114 constexpr size_t
115 GetNumberOfOffsets() const noexcept
116 {
117 return m_NumberOfOffsets;
118 }
119
120
122 void
123 FillOffsets(Offset<ImageDimension> * const offsets) const noexcept
124 {
125 if (m_NumberOfOffsets > 0)
126 {
127 assert(offsets != nullptr);
129 std::fill_n(offset.begin(), ImageDimension, -1);
132 size_t i = 0;
133
134 while (i < m_NumberOfOffsets)
135 {
136 const size_t numberOfNonZeroOffsetValues =
137 ImageDimension - static_cast<size_t>(std::count(offset.begin(), offset.end(), 0));
138
139 if ((m_IncludeCenterPixel || (numberOfNonZeroOffsetValues > 0)) &&
140 (numberOfNonZeroOffsetValues <= m_MaximumCityblockDistance))
141 {
142 offsets[i] = offset;
143 ++i;
144 }
145
146 // Go to the next offset:
147 for (unsigned int direction = 0; direction < ImageDimension; ++direction)
148 {
149 auto & offsetValue = offset[direction];
150
151 ++offsetValue;
152
153 if (offsetValue <= 1)
154 {
155 break;
156 }
157 offsetValue = -1;
158 }
159 }
160 }
161 }
162
163private:
164 // The maximum city-block distance (Manhattan distance) between the center
165 // pixel and each connected neighbor pixel.
167
168 // Specifies whether or not the center pixel (offset zero) should be included
169 // with the offsets for this shape.
171
172 // The number of offsets needed to represent this shape.
174
175
176 // Calculates a + b. Numeric overflow triggers a compilation error in
177 // "constexpr context" and a debug assert failure at run-time.
178 static constexpr uintmax_t
179 CalculateSum(const uintmax_t a, const uintmax_t b) noexcept
180 {
181 return ((a + b) >= a) && ((a + b) >= b) ? (a + b) : (assert(!"CalculateSum overflow!"), 0);
182 }
183
184
185 // Calculates 2 ^ n. Numeric overflow triggers a compilation error in
186 // "constexpr context" and a debug assert failure at run-time.
187 static constexpr uintmax_t
188 CalculatePowerOfTwo(const size_t n) noexcept
189 {
190 return (n < std::numeric_limits<uintmax_t>::digits) ? (uintmax_t{ 1 } << n)
191 : (assert(!"CalculatePowerOfTwo overflow!"), 0);
192 }
193
194
195 // Calculates the binomial coefficient, 'n' over 'k'.
196 // Inspired by the 'binom' function from Walter, June 23, 2017:
197 // https://stackoverflow.com/questions/44718971/calculate-binomial-coffeficient-very-reliable/44719165#44719165
198 // Optimized for small values of 'k' (k <= n/2).
199 static constexpr uintmax_t
200 CalculateBinomialCoefficient(const uintmax_t n, const uintmax_t k) noexcept
201 {
202 return (k > n) ? (assert(!"Out of range!"), 0)
203 : (k == 0) ? 1 : Math::UnsignedProduct(n, CalculateBinomialCoefficient(n - 1, k - 1)) / k;
204 }
205
206
207 // Calculates the number of m-dimensional hypercubes on the boundary of an
208 // n-cube, as described at https://en.wikipedia.org/wiki/Hypercube#Elements
209 // (Which refers to H.S.M. Coxeter, Regular polytopes, 3rd ed., 1973, p.120.)
210 static constexpr size_t
211 CalculateNumberOfHypercubesOnBoundaryOfCube(const size_t m, const size_t n) noexcept
212 {
213 // Calculate 2^(n-m) * BinomialCoefficient(n, m)
215 (((2 * m) < n) ?
216 // Calculate either the binomial coefficient of (n, m) or (n, n - m).
217 // Mathematically, both should yield the same number, but the
218 // calculation is optimized for a smaller second argument.
220 : CalculateBinomialCoefficient(n, n - m)));
221 }
222
223
224 // Iterates recursively from i = ImageDimension-1 down to m (inclusive).
225 static constexpr size_t
226 CalculateSumOfNumberOfHypercubesOnBoundaryOfCube(const size_t i, const size_t m) noexcept
227 {
229 ((i <= m) ? 0 : CalculateSumOfNumberOfHypercubesOnBoundaryOfCube(i - 1, m)));
230 }
231
232
234 static constexpr size_t
235 CalculateNumberOfConnectedNeighbors(const size_t maximumCityblockDistance) noexcept
236 {
237 return (((maximumCityblockDistance == 0) || (ImageDimension == 0))
238 ? 0
239 : ((maximumCityblockDistance >= ImageDimension)
242 (ImageDimension - maximumCityblockDistance))));
243 }
244
245
247 static constexpr size_t
248 CalculateNumberOfOffsets(const size_t maximumCityblockDistance, const bool includeCenterPixel) noexcept
249 {
250 return (includeCenterPixel ? 1 : 0) + CalculateNumberOfConnectedNeighbors(maximumCityblockDistance);
251 }
252};
257template <unsigned int VImageDimension, size_t VMaximumCityblockDistance, bool VIncludeCenterPixel>
258auto
260{
261 constexpr ConnectedImageNeighborhoodShape<VImageDimension> shape{ VMaximumCityblockDistance, VIncludeCenterPixel };
262 std::array<Offset<VImageDimension>, shape.GetNumberOfOffsets()> offsets;
263 shape.FillOffsets(offsets.data());
264 return offsets;
265}
268} // namespace itk
269
270#endif
static constexpr vcl_size_t CalculateNumberOfConnectedNeighbors(const vcl_size_t maximumCityblockDistance) noexcept
static constexpr vcl_size_t CalculateNumberOfOffsets(const vcl_size_t maximumCityblockDistance, const bool includeCenterPixel) noexcept
static constexpr uintmax_t CalculatePowerOfTwo(const vcl_size_t n) noexcept
constexpr vcl_size_t GetNumberOfOffsets() const noexcept
static constexpr vcl_size_t CalculateSumOfNumberOfHypercubesOnBoundaryOfCube(const vcl_size_t i, const vcl_size_t m) noexcept
static constexpr uintmax_t CalculateSum(const uintmax_t a, const uintmax_t b) noexcept
static constexpr vcl_size_t CalculateNumberOfHypercubesOnBoundaryOfCube(const vcl_size_t m, const vcl_size_t n) noexcept
constexpr ConnectedImageNeighborhoodShape(const vcl_size_t maximumCityblockDistance, const bool includeCenterPixel) noexcept
static constexpr uintmax_t CalculateBinomialCoefficient(const uintmax_t n, const uintmax_t k) noexcept
void FillOffsets(Offset< ImageDimension > *const offsets) const noexcept
constexpr TReturnType UnsignedPower(const uintmax_t base, const uintmax_t exponent) noexcept
Definition: itkMath.h:797
constexpr TReturnType UnsignedProduct(const uintmax_t a, const uintmax_t b) noexcept
Definition: itkMath.h:775
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
auto GenerateConnectedImageNeighborhoodShapeOffsets() noexcept
constexpr iterator begin()
Definition: itkOffset.h:316
constexpr iterator end()
Definition: itkOffset.h:334