ITK  4.3.0
Insight Segmentation and Registration Toolkit
itkHistogramMatchingImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
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  * http://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 __itkHistogramMatchingImageFilter_h
19 #define __itkHistogramMatchingImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
22 #include "itkHistogram.h"
23 #include "vnl/vnl_matrix.h"
24 
25 namespace itk
26 {
65 /* THistogramMeasurement -- The precision level for which to do
66  HistogramMeasurmenets */
67 template< class TInputImage, class TOutputImage, class THistogramMeasurement = typename TInputImage::PixelType >
69  public ImageToImageFilter< TInputImage, TOutputImage >
70 {
71 public:
77 
79  itkNewMacro(Self);
80 
83 
85  itkStaticConstMacro(ImageDimension, unsigned int,
86  TInputImage::ImageDimension);
87  itkStaticConstMacro(OutputImageDimension, unsigned int,
88  TOutputImage::ImageDimension);
90 
92  typedef typename TOutputImage::RegionType OutputImageRegionType;
93 
96  typedef typename Superclass::InputImagePointer InputImagePointer;
97  typedef typename Superclass::InputImageConstPointer InputImageConstPointer;
98  typedef typename Superclass::OutputImageType OutputImageType;
99  typedef typename Superclass::OutputImagePointer OutputImagePointer;
100 
102  typedef typename InputImageType::PixelType InputPixelType;
103  typedef typename OutputImageType::PixelType OutputPixelType;
104 
108 
110  void SetSourceImage(const InputImageType *source)
111  { this->SetInput(source); }
112  const InputImageType * GetSourceImage(void)
113  { return this->GetInput(); }
115 
117  void SetReferenceImage(const InputImageType *reference);
118 
119  const InputImageType * GetReferenceImage(void);
120 
122  itkSetMacro(NumberOfHistogramLevels, SizeValueType);
123  itkGetConstMacro(NumberOfHistogramLevels, SizeValueType);
125 
127  itkSetMacro(NumberOfMatchPoints, SizeValueType);
128  itkGetConstMacro(NumberOfMatchPoints, SizeValueType);
130 
136  itkSetMacro(ThresholdAtMeanIntensity, bool);
137  itkGetConstMacro(ThresholdAtMeanIntensity, bool);
138  itkBooleanMacro(ThresholdAtMeanIntensity);
140 
142  virtual void GenerateInputRequestedRegion();
143 
147  itkGetObjectMacro(SourceHistogram, HistogramType);
148  itkGetObjectMacro(ReferenceHistogram, HistogramType);
149  itkGetObjectMacro(OutputHistogram, HistogramType);
151 
152 #ifdef ITK_USE_CONCEPT_CHECKING
153 
154  itkConceptMacro( IntConvertibleToInputCheck,
156  itkConceptMacro( SameDimensionCheck,
158  itkConceptMacro( DoubleConvertibleToInputCheck,
160  itkConceptMacro( DoubleConvertibleToOutputCheck,
162  itkConceptMacro( InputConvertibleToDoubleCheck,
164  itkConceptMacro( OutputConvertibleToDoubleCheck,
166  itkConceptMacro( SameTypeCheck,
168 
170 #endif
171 
172 protected:
175  void PrintSelf(std::ostream & os, Indent indent) const;
176 
177  void BeforeThreadedGenerateData();
178 
179  void AfterThreadedGenerateData();
180 
181  void ThreadedGenerateData(const OutputImageRegionType & outputRegionForThread,
182  ThreadIdType threadId);
183 
189  virtual void VerifyInputInformation() {}
190 
192  void ComputeMinMaxMean(const InputImageType *image,
193  THistogramMeasurement & minValue,
194  THistogramMeasurement & maxValue,
195  THistogramMeasurement & meanValue);
196 
198  void ConstructHistogram(const InputImageType *image,
199  HistogramType *histogram, const THistogramMeasurement minValue,
200  const THistogramMeasurement maxValue);
201 
202 private:
203  HistogramMatchingImageFilter(const Self &); //purposely not implemented
204  void operator=(const Self &); //purposely not implemented
205 
209 
213 
214  THistogramMeasurement m_SourceMinValue;
215  THistogramMeasurement m_SourceMaxValue;
216  THistogramMeasurement m_SourceMeanValue;
217  THistogramMeasurement m_ReferenceMinValue;
218  THistogramMeasurement m_ReferenceMaxValue;
219  THistogramMeasurement m_ReferenceMeanValue;
220  THistogramMeasurement m_OutputMinValue;
221  THistogramMeasurement m_OutputMaxValue;
222  THistogramMeasurement m_OutputMeanValue;
223 
227 
230 
235 };
236 } // end namespace itk
237 
238 #ifndef ITK_MANUAL_INSTANTIATION
239 #include "itkHistogramMatchingImageFilter.hxx"
240 #endif
241 
242 #endif
243