ITK  6.0.0
Insight Toolkit
itkConstNeighborhoodIterator.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#ifndef itkConstNeighborhoodIterator_h
19#define itkConstNeighborhoodIterator_h
20
21#include <vector>
22#include <cstring>
23#include <iostream>
24#include "itkImage.h"
25#include "itkNeighborhood.h"
26#include "itkMacro.h"
28
29namespace itk
30{
50template <typename TImage, typename TBoundaryCondition = ZeroFluxNeumannBoundaryCondition<TImage>>
51class ITK_TEMPLATE_EXPORT ConstNeighborhoodIterator
52 : public Neighborhood<typename TImage::InternalPixelType *, TImage::ImageDimension>
53{
54public:
56 using InternalPixelType = typename TImage::InternalPixelType;
57 using PixelType = typename TImage::PixelType;
58
60 using DimensionValueType = unsigned int;
61
63 static constexpr DimensionValueType Dimension = TImage::ImageDimension;
64
68
70 using typename Superclass::OffsetType;
71 using typename Superclass::RadiusType;
72 using typename Superclass::SizeType;
73 using typename Superclass::Iterator;
74 using typename Superclass::ConstIterator;
75
77 using ImageType = TImage;
81
82
84 using BoundaryConditionType = TBoundaryCondition;
86 using OutputImageType = typename BoundaryConditionType::OutputImageType;
87
90
94 using NeighborhoodAccessorFunctorType = typename ImageType::NeighborhoodAccessorFunctorType;
95
99
101 ConstNeighborhoodIterator() = default;
102
104 ~ConstNeighborhoodIterator() override = default;
105
108
111 ConstNeighborhoodIterator(const SizeType & radius, const ImageType * ptr, const RegionType & region)
112 {
113 this->Initialize(radius, ptr, region);
114 for (DimensionValueType i = 0; i < Dimension; ++i)
115 {
116 m_InBounds[i] = false;
117 }
118 this->ResetBoundaryCondition();
119 m_NeighborhoodAccessorFunctor = ptr->GetNeighborhoodAccessor();
120 m_NeighborhoodAccessorFunctor.SetBegin(ptr->GetBufferPointer());
121 }
126 operator=(const Self & orig);
127
129 void
130 PrintSelf(std::ostream &, Indent) const override;
131
135 ComputeInternalIndex(const NeighborIndexType n) const;
136
139 GetBound() const
140 {
141 return m_Bound;
142 }
143
147 GetBound(NeighborIndexType n) const
148 {
149 return m_Bound[n];
150 }
151
154 GetCenterPointer() const
155 {
156 return (this->operator[]((this->Size()) >> 1));
157 }
158
162 GetCenterPixel() const
163 {
164 return m_NeighborhoodAccessorFunctor.Get(this->GetCenterPointer());
165 }
166
168 const ImageType *
169 GetImagePointer() const
170 {
171 return m_ConstImage;
172 }
173
177 GetIndex() const
178 {
179 return m_Loop;
180 }
181
183 GetFastIndexPlusOffset(const OffsetType & o) const
184 {
185 return m_Loop + o;
186 }
187
191 GetNeighborhood() const;
192
195 GetPixel(const NeighborIndexType i) const
196 {
197 if (!m_NeedToUseBoundaryCondition || this->InBounds())
198 {
199 return (m_NeighborhoodAccessorFunctor.Get(this->operator[](i)));
200 }
201
202 OffsetType internalIndex;
203 OffsetType offset;
204
205 return this->IndexInBounds(i, internalIndex, offset)
206 ? m_NeighborhoodAccessorFunctor.Get(this->operator[](i))
207 : m_NeighborhoodAccessorFunctor.BoundaryCondition(internalIndex, offset, this, m_BoundaryCondition);
208 }
209
216 GetPixel(NeighborIndexType i, bool & IsInBounds) const;
217
221 GetPixel(const OffsetType & o) const
222 {
223 bool inbounds;
224
225 return (this->GetPixel(this->GetNeighborhoodIndex(o), inbounds));
226 }
227
234 GetPixel(const OffsetType & o, bool & IsInBounds) const
235 {
236 return (this->GetPixel(this->GetNeighborhoodIndex(o), IsInBounds));
237 }
238
243 GetNext(const unsigned int axis, NeighborIndexType i) const
244 {
245 return (this->GetPixel(this->GetCenterNeighborhoodIndex() + (i * this->GetStride(axis))));
246 }
247
252 GetNext(const unsigned int axis) const
253 {
254 return (this->GetPixel(this->GetCenterNeighborhoodIndex() + this->GetStride(axis)));
255 }
256
261 GetPrevious(const unsigned int axis, NeighborIndexType i) const
262 {
263 return (this->GetPixel(this->GetCenterNeighborhoodIndex() - (i * this->GetStride(axis))));
264 }
265
270 GetPrevious(const unsigned int axis) const
271 {
272 return (this->GetPixel(this->GetCenterNeighborhoodIndex() - this->GetStride(axis)));
273 }
274
278 GetIndex(const OffsetType & o) const
279 {
280 return (this->GetIndex() + o);
281 }
282
286 GetIndex(NeighborIndexType i) const
287 {
288 return (this->GetIndex() + this->GetOffset(i));
289 }
290
293 GetRegion() const
294 {
295 return m_Region;
296 }
297
301 GetBeginIndex() const
302 {
303 return m_BeginIndex;
304 }
305
309 GetBoundingBoxAsImageRegion() const;
310
313 GetWrapOffset() const
314 {
315 return m_WrapOffset;
316 }
317
324 GetWrapOffset(NeighborIndexType n) const
325 {
326 return m_WrapOffset[n];
327 }
328
330 void
331 GoToBegin();
332
335 void
336 GoToEnd();
337
340 void
341 Initialize(const SizeType & radius, const ImageType * ptr, const RegionType & region);
342
345 bool
346 IsAtBegin() const
347 {
348 return (this->GetCenterPointer() == m_Begin);
349 }
350
353 bool
354 IsAtEnd() const
355 {
356 if (this->GetCenterPointer() > m_End)
357 {
358 ExceptionObject e(__FILE__, __LINE__);
359 std::ostringstream msg;
360 msg << "In method IsAtEnd, CenterPointer = " << this->GetCenterPointer() << " is greater than End = " << m_End
361 << std::endl
362 << " " << *this;
363 e.SetDescription(msg.str().c_str());
364 throw e;
365 }
366 return (this->GetCenterPointer() == m_End);
367 }
375 operator++();
376
382 operator--();
383
387 bool
388 operator==(const Self & it) const
389 {
390 return it.GetCenterPointer() == this->GetCenterPointer();
391 }
393 ITK_UNEQUAL_OPERATOR_MEMBER_FUNCTION(Self);
394
398 bool
399 operator<(const Self & it) const
400 {
401 return this->GetCenterPointer() < it.GetCenterPointer();
402 }
403
407 bool
408 operator<=(const Self & it) const
409 {
410 return this->GetCenterPointer() <= it.GetCenterPointer();
411 }
412
416 bool
417 operator>(const Self & it) const
418 {
419 return this->GetCenterPointer() > it.GetCenterPointer();
420 }
421
425 bool
426 operator>=(const Self & it) const
427 {
428 return this->GetCenterPointer() >= it.GetCenterPointer();
429 }
430
435 void
436 SetLocation(const IndexType & position)
437 {
438 this->SetLoop(position);
439 this->SetPixelPointers(position);
440 }
447 operator+=(const OffsetType &);
448
453 operator-=(const OffsetType &);
454
457 operator-(const Self & b) const
458 {
459 return m_Loop - b.m_Loop;
460 }
461
465 bool
466 InBounds() const;
467
479 bool
480 IndexInBounds(const NeighborIndexType n, OffsetType & internalIndex, OffsetType & offset) const;
481
484 bool
485 IndexInBounds(const NeighborIndexType n) const;
486
492 void
493 OverrideBoundaryCondition(const ImageBoundaryConditionPointerType i)
494 {
495 m_BoundaryCondition = i;
496 }
497
500 void
501 ResetBoundaryCondition()
502 {
503 m_BoundaryCondition = &m_InternalBoundaryCondition;
504 }
505
507 void
508 SetBoundaryCondition(const TBoundaryCondition & c)
509 {
510 m_InternalBoundaryCondition = c;
511 }
512
515 GetBoundaryCondition() const
516 {
517 return m_BoundaryCondition;
518 }
519
521 void
522 NeedToUseBoundaryConditionOn()
523 {
524 this->SetNeedToUseBoundaryCondition(true);
525 }
526
527 void
528 NeedToUseBoundaryConditionOff()
529 {
530 this->SetNeedToUseBoundaryCondition(false);
531 }
532
533 void
534 SetNeedToUseBoundaryCondition(bool b)
535 {
536 m_NeedToUseBoundaryCondition = b;
537 }
538
539 bool
540 GetNeedToUseBoundaryCondition() const
541 {
542 return m_NeedToUseBoundaryCondition;
543 }
544
546 void
547 SetRegion(const RegionType & region);
548
549protected:
552 void
553 SetLoop(const IndexType & p)
554 {
555 m_Loop = p;
556 m_IsInBoundsValid = false;
557 }
558
562 void
563 SetBound(const SizeType &);
564
569 void
570 SetPixelPointers(const IndexType &);
571
574 void
575 SetBeginIndex(const IndexType & start)
576 {
577 m_BeginIndex = start;
578 }
579
582 void
583 SetEndIndex();
584
587 IndexType m_BeginIndex{ { 0 } };
588
590 IndexType m_Bound{ { 0 } };
591
593 const InternalPixelType * m_Begin{ nullptr };
594
596 typename ImageType::ConstWeakPointer m_ConstImage{};
597
599 const InternalPixelType * m_End{ nullptr };
600
603 IndexType m_EndIndex{ { 0 } };
604
606 IndexType m_Loop{ { 0 } };
607
609 RegionType m_Region{};
610
616 OffsetType m_WrapOffset{ { 0 } };
617
619 TBoundaryCondition m_InternalBoundaryCondition{};
620
625 ImageBoundaryConditionPointerType m_BoundaryCondition{ &m_InternalBoundaryCondition };
626
629 mutable bool m_InBounds[Dimension]{ false };
630
632 mutable bool m_IsInBounds{ false };
633
637 mutable bool m_IsInBoundsValid{ false };
638
640 IndexType m_InnerBoundsLow{};
641
643 IndexType m_InnerBoundsHigh{};
644
646 bool m_NeedToUseBoundaryCondition{ false };
647
649 NeighborhoodAccessorFunctorType m_NeighborhoodAccessorFunctor{};
650};
651
652template <typename TImage>
656{
658 ret += ind;
659 return ret;
660}
661
662template <typename TImage>
666{
667 return (it + ind);
668}
669
670template <typename TImage>
674{
676 ret -= ind;
677 return ret;
678}
679} // namespace itk
680
681#ifndef ITK_MANUAL_INSTANTIATION
682# include "itkConstNeighborhoodIterator.hxx"
683#endif
684
685#endif
Pixel-wise addition of two images.
Const version of NeighborhoodIterator, defining iteration of a local N-dimensional neighborhood of pi...
const InternalPixelType * GetCenterPointer() const
typename TImage::RegionType RegionType
typename ImageType::NeighborhoodAccessorFunctorType NeighborhoodAccessorFunctorType
typename TImage::InternalPixelType InternalPixelType
typename NeighborhoodType::NeighborIndexType NeighborIndexType
typename BoundaryConditionType::OutputImageType OutputImageType
Control indentation during Print() invocation.
Definition: itkIndent.h:50
A light-weight container object for storing an N-dimensional neighborhood of values.
typename AllocatorType::iterator Iterator
SizeValueType NeighborIndexType
typename AllocatorType::const_iterator ConstIterator
static constexpr double e
Definition: itkMath.h:56
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
bool operator>(const Index< VDimension > &one, const Index< VDimension > &two)
Definition: itkIndex.h:564
bool operator<=(const Index< VDimension > &one, const Index< VDimension > &two)
Definition: itkIndex.h:571
long IndexValueType
Definition: itkIntTypes.h:93
ConstNeighborhoodIterator< TImage > operator-(const ConstNeighborhoodIterator< TImage > &it, const typename ConstNeighborhoodIterator< TImage >::OffsetType &ind)
bool operator==(const Index< VDimension > &one, const Index< VDimension > &two)
Definition: itkIndex.h:543
ConstNeighborhoodIterator< TImage > operator+(const ConstNeighborhoodIterator< TImage > &it, const typename ConstNeighborhoodIterator< TImage >::OffsetType &ind)
bool operator<(const Index< VDimension > &one, const Index< VDimension > &two)
Definition: itkIndex.h:557
bool operator>=(const Index< VDimension > &one, const Index< VDimension > &two)
Definition: itkIndex.h:578
long OffsetValueType
Definition: itkIntTypes.h:97