ITK  6.0.0
Insight Toolkit
itkMRFImageFilter.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 itkMRFImageFilter_h
19#define itkMRFImageFilter_h
20
21#include "vnl/vnl_vector.h"
22#include "vnl/vnl_matrix.h"
23
25
27
30#include "itkSize.h"
31#include "ITKMarkovRandomFieldsClassifiersExport.h"
32
33namespace itk
34{
40{
41public:
47 enum class MRFStop : uint8_t
48 {
49 MaximumNumberOfIterations = 1,
50 ErrorTolerance
51 };
52};
53
54// Define how to print enumeration
55extern ITKMarkovRandomFieldsClassifiers_EXPORT std::ostream &
56 operator<<(std::ostream & out, const MRFImageFilterEnums::MRFStop value);
148template <typename TInputImage, typename TClassifiedImage>
149class ITK_TEMPLATE_EXPORT MRFImageFilter : public ImageToImageFilter<TInputImage, TClassifiedImage>
150{
151public:
152 ITK_DISALLOW_COPY_AND_MOVE(MRFImageFilter);
160 using typename Superclass::OutputImagePointer;
161
163 itkNewMacro(Self);
164
166 itkOverrideGetNameOfClassMacro(MRFImageFilter);
167
169 using InputImageType = TInputImage;
172
174 using InputImagePixelType = typename TInputImage::PixelType;
175
178
182
184 static constexpr unsigned int InputImageDimension = TInputImage::ImageDimension;
185
188
190 using TrainingImagePixelType = typename TClassifiedImage::PixelType;
191
195
198 using LabelledImagePixelType = typename TClassifiedImage::PixelType;
199
203
207
209 using LabelledImageOffsetType = typename TClassifiedImage::OffsetType;
210
213
215 static constexpr unsigned int ClassifiedImageDimension = TClassifiedImage::ImageDimension;
216
219
222
225
228
230
232
234
235 using InputImageFaceListIterator = typename InputImageFaceListType::iterator;
236
239
241
243
245
246 using LabelledImageFaceListIterator = typename LabelledImageFaceListType::iterator;
247
249 void
250 SetClassifier(typename ClassifierType::Pointer ptrToClassifier);
251
253 itkSetMacro(NumberOfClasses, unsigned int);
254 itkGetConstMacro(NumberOfClasses, unsigned int);
259 itkSetMacro(MaximumNumberOfIterations, unsigned int);
260 itkGetConstMacro(MaximumNumberOfIterations, unsigned int);
265 itkSetMacro(ErrorTolerance, double);
266 itkGetConstMacro(ErrorTolerance, double);
271 itkSetMacro(SmoothingFactor, double);
272 itkGetConstMacro(SmoothingFactor, double);
276 void
278
283 void
285
286 void
288
292 {
294
295 for (int i = 0; i < InputImageDimension; ++i)
296 {
297 radius[i] = m_InputImageNeighborhoodRadius[i];
298 }
299 return radius;
300 }
301
307 virtual void
308 SetMRFNeighborhoodWeight(std::vector<double> betaMatrix);
309
310 virtual std::vector<double>
312 {
313 return m_MRFNeighborhoodWeight;
314 }
315
317#if !defined(ITK_LEGACY_REMOVE)
319 // We need to expose the enum values at the class level
320 // for backwards compatibility
321 static constexpr MRFStopEnum MaximumNumberOfIterations = MRFStopEnum::MaximumNumberOfIterations;
322 static constexpr MRFStopEnum ErrorTolerance = MRFStopEnum::ErrorTolerance;
323#endif
324
327 itkGetConstReferenceMacro(StopCondition, MRFStopEnum);
328
329 /* Get macro for number of iterations */
330 itkGetConstReferenceMacro(NumberOfIterations, unsigned int);
331
332#ifdef ITK_USE_CONCEPT_CHECKING
333 // Begin concept checking
334 itkConceptMacro(UnsignedIntConvertibleToClassifiedCheck,
336 itkConceptMacro(ClassifiedConvertibleToUnsignedIntCheck,
338 itkConceptMacro(ClassifiedConvertibleToIntCheck, (Concept::Convertible<LabelledImagePixelType, int>));
339 itkConceptMacro(IntConvertibleToClassifiedCheck, (Concept::Convertible<int, LabelledImagePixelType>));
341 // End concept checking
342#endif
343
344protected:
346 ~MRFImageFilter() override = default;
347 void
348 PrintSelf(std::ostream & os, Indent indent) const override;
349
351 void
353
358 virtual void
360
362 virtual void
364
370
373
375 virtual void
378 LabelStatusImageNeighborhoodIterator & labelStatusIter);
379
380 void
381 GenerateData() override;
382
383 void
385
386 void
388
389 void
391
392private:
394
396
398
400
401 using LabelStatusImageFaceListIterator = typename LabelStatusImageFaceListType::iterator;
402
403 InputImageNeighborhoodRadiusType m_InputImageNeighborhoodRadius{};
404 LabelledImageNeighborhoodRadiusType m_LabelledImageNeighborhoodRadius{};
405 LabelStatusImageNeighborhoodRadiusType m_LabelStatusImageNeighborhoodRadius{};
406
407 unsigned int m_NumberOfClasses{ 0 };
408 unsigned int m_MaximumNumberOfIterations{ 50 };
409 unsigned int m_KernelSize{};
410
411 int m_ErrorCounter{ 0 };
412 int m_NeighborhoodSize{ 27 };
413 int m_TotalNumberOfValidPixelsInOutputImage{ 1 };
414 int m_TotalNumberOfPixelsInInputImage{ 1 };
415 double m_ErrorTolerance{ 0.2 };
416 double m_SmoothingFactor{ 1 };
417 double * m_ClassProbability{ nullptr }; // Class likelihood
418 unsigned int m_NumberOfIterations{ 0 };
419 MRFStopEnum m_StopCondition{ MRFStopEnum::MaximumNumberOfIterations };
420
421 LabelStatusImagePointer m_LabelStatusImage{};
422
423 std::vector<double> m_MRFNeighborhoodWeight{};
424 std::vector<double> m_NeighborInfluence{};
425 std::vector<double> m_MahalanobisDistance{};
426 std::vector<double> m_DummyVector{};
427
429 typename ClassifierType::Pointer m_ClassifierPtr{};
430
434 virtual void
436
438 void
440}; // class MRFImageFilter
441} // namespace itk
442
443#ifndef ITK_MANUAL_INSTANTIATION
444# include "itkMRFImageFilter.hxx"
445#endif
446
447#endif
Base class for all data objects in ITK.
Base class for the ImageClassifierBase object.
A multi-dimensional iterator templated over image type that walks a region of pixels.
A multi-dimensional iterator templated over image type that walks a region of pixels.
Base class for filters that take an image as input and produce an image as output.
Templated n-dimensional image class.
Definition: itkImage.h:89
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Contains all enum classes in MRFImageFilter class;.
Implementation of a labeler object that uses Markov Random Fields to classify pixels in an image data...
virtual void SetMRFNeighborhoodWeight(std::vector< double > betaMatrix)
typename TClassifiedImage::OffsetType LabelledImageOffsetType
typename LabelledImageIndexType::IndexValueType IndexValueType
void GenerateData() override
virtual void DoNeighborhoodOperation(const InputImageNeighborhoodIterator &imageIter, LabelledImageNeighborhoodIterator &labelledIter, LabelStatusImageNeighborhoodIterator &labelStatusIter)
void SetNeighborhoodRadius(const SizeValueType)
typename InputImageNeighborhoodIterator::RadiusType InputImageNeighborhoodRadiusType
void SetClassifier(typename ClassifierType::Pointer ptrToClassifier)
typename TInputImage::RegionType InputImageRegionType
typename TClassifiedImage::Pointer TrainingImagePointer
void SetNeighborhoodRadius(const NeighborhoodRadiusType &)
void GenerateOutputInformation() override
typename InputImageFacesCalculator::FaceListType InputImageFaceListType
typename TInputImage::SizeType SizeType
virtual void MinimizeFunctional()
const NeighborhoodRadiusType GetNeighborhoodRadius() const
typename TInputImage::ConstPointer InputImageConstPointer
typename InputImageFaceListType::iterator InputImageFaceListIterator
typename LabelStatusImageType::RegionType LabelStatusRegionType
virtual void SetDefaultMRFNeighborhoodWeight()
typename TInputImage::SizeType InputImageSizeType
void GenerateInputRequestedRegion() override
typename TClassifiedImage::PixelType LabelledImagePixelType
typename TInputImage::SizeType NeighborhoodRadiusType
void PrintSelf(std::ostream &os, Indent indent) const override
typename LabelStatusImageFacesCalculator::FaceListType LabelStatusImageFaceListType
void SetNeighborhoodRadius(const SizeValueType *radiusArray)
typename TClassifiedImage::IndexType LabelledImageIndexType
~MRFImageFilter() override=default
typename TInputImage::Pointer InputImagePointer
typename TClassifiedImage::Pointer LabelledImagePointer
typename TClassifiedImage::RegionType LabelledImageRegionType
typename TClassifiedImage::PixelType TrainingImagePixelType
typename LabelledImageFaceListType::iterator LabelledImageFaceListIterator
void EnlargeOutputRequestedRegion(DataObject *) override
typename TInputImage::PixelType InputImagePixelType
virtual void ApplyMRFImageFilter()
typename LabelStatusImageFaceListType::iterator LabelStatusImageFaceListIterator
typename LabelStatusImageNeighborhoodIterator::RadiusType LabelStatusImageNeighborhoodRadiusType
typename LabelStatusImageType::IndexType LabelStatusIndexType
typename LabelledImageFacesCalculator::FaceListType LabelledImageFaceListType
typename LabelStatusImageType::Pointer LabelStatusImagePointer
typename LabelledImageNeighborhoodIterator::RadiusType LabelledImageNeighborhoodRadiusType
virtual std::vector< double > GetMRFNeighborhoodWeight()
Defines iteration of a local N-dimensional neighborhood of pixels across an itk::Image.
SmartPointer< const Self > ConstPointer
SmartPointer< Self > Pointer
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
ITKCommon_EXPORT std::ostream & operator<<(std::ostream &out, typename AnatomicalOrientation::CoordinateEnum value)
long IndexValueType
Definition: itkIntTypes.h:93
unsigned long SizeValueType
Definition: itkIntTypes.h:86
Splits an image into a main region and several "face" regions which are used to handle computations o...