ITK  4.1.0
Insight Segmentation and Registration Toolkit
itkStatisticsImageFilter.h
Go to the documentation of this file.
00001 /*=========================================================================
00002  *
00003  *  Copyright Insight Software Consortium
00004  *
00005  *  Licensed under the Apache License, Version 2.0 (the "License");
00006  *  you may not use this file except in compliance with the License.
00007  *  You may obtain a copy of the License at
00008  *
00009  *         http://www.apache.org/licenses/LICENSE-2.0.txt
00010  *
00011  *  Unless required by applicable law or agreed to in writing, software
00012  *  distributed under the License is distributed on an "AS IS" BASIS,
00013  *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
00014  *  See the License for the specific language governing permissions and
00015  *  limitations under the License.
00016  *
00017  *=========================================================================*/
00018 #ifndef __itkStatisticsImageFilter_h
00019 #define __itkStatisticsImageFilter_h
00020 
00021 #include "itkImageToImageFilter.h"
00022 #include "itkNumericTraits.h"
00023 #include "itkArray.h"
00024 #include "itkSimpleDataObjectDecorator.h"
00025 
00026 namespace itk
00027 {
00048 template< class TInputImage >
00049 class ITK_EXPORT StatisticsImageFilter:
00050   public ImageToImageFilter< TInputImage, TInputImage >
00051 {
00052 public:
00054   typedef StatisticsImageFilter                          Self;
00055   typedef ImageToImageFilter< TInputImage, TInputImage > Superclass;
00056   typedef SmartPointer< Self >                           Pointer;
00057   typedef SmartPointer< const Self >                     ConstPointer;
00058 
00060   itkNewMacro(Self);
00061 
00063   itkTypeMacro(StatisticsImageFilter, ImageToImageFilter);
00064 
00066   typedef typename TInputImage::Pointer InputImagePointer;
00067 
00068   typedef typename TInputImage::RegionType RegionType;
00069   typedef typename TInputImage::SizeType   SizeType;
00070   typedef typename TInputImage::IndexType  IndexType;
00071   typedef typename TInputImage::PixelType  PixelType;
00072 
00074   itkStaticConstMacro(ImageDimension, unsigned int,
00075                       TInputImage::ImageDimension);
00076 
00078   typedef typename NumericTraits< PixelType >::RealType RealType;
00079 
00081   typedef typename DataObject::Pointer DataObjectPointer;
00082 
00084   typedef SimpleDataObjectDecorator< RealType >  RealObjectType;
00085   typedef SimpleDataObjectDecorator< PixelType > PixelObjectType;
00086 
00088   PixelType GetMinimum() const
00089   { return this->GetMinimumOutput()->Get(); }
00090   PixelObjectType * GetMinimumOutput();
00092 
00093   const PixelObjectType * GetMinimumOutput() const;
00094 
00096   PixelType GetMaximum() const
00097   { return this->GetMaximumOutput()->Get(); }
00098   PixelObjectType * GetMaximumOutput();
00100 
00101   const PixelObjectType * GetMaximumOutput() const;
00102 
00104   RealType GetMean() const
00105   { return this->GetMeanOutput()->Get(); }
00106   RealObjectType * GetMeanOutput();
00108 
00109   const RealObjectType * GetMeanOutput() const;
00110 
00112   RealType GetSigma() const
00113   { return this->GetSigmaOutput()->Get(); }
00114   RealObjectType * GetSigmaOutput();
00116 
00117   const RealObjectType * GetSigmaOutput() const;
00118 
00120   RealType GetVariance() const
00121   { return this->GetVarianceOutput()->Get(); }
00122   RealObjectType * GetVarianceOutput();
00124 
00125   const RealObjectType * GetVarianceOutput() const;
00126 
00128   RealType GetSum() const
00129   { return this->GetSumOutput()->Get(); }
00130   RealObjectType * GetSumOutput();
00132 
00133   const RealObjectType * GetSumOutput() const;
00134 
00137   typedef ProcessObject::DataObjectPointerArraySizeType DataObjectPointerArraySizeType;
00138   using Superclass::MakeOutput;
00139   virtual DataObjectPointer MakeOutput(DataObjectPointerArraySizeType idx);
00140 
00141 #ifdef ITK_USE_CONCEPT_CHECKING
00142 
00143   itkConceptMacro( InputHasNumericTraitsCheck,
00144                    ( Concept::HasNumericTraits< PixelType > ) );
00145 
00147 #endif
00148 protected:
00149   StatisticsImageFilter();
00150   ~StatisticsImageFilter(){}
00151   void PrintSelf(std::ostream & os, Indent indent) const;
00153 
00157   void AllocateOutputs();
00158 
00160   void BeforeThreadedGenerateData();
00161 
00164   void AfterThreadedGenerateData();
00165 
00167   void  ThreadedGenerateData(const RegionType &
00168                              outputRegionForThread,
00169                              ThreadIdType threadId);
00170 
00171   // Override since the filter needs all the data for the algorithm
00172   void GenerateInputRequestedRegion();
00173 
00174   // Override since the filter produces all of its output
00175   void EnlargeOutputRequestedRegion(DataObject *data);
00176 
00177 private:
00178   StatisticsImageFilter(const Self &); //purposely not implemented
00179   void operator=(const Self &);        //purposely not implemented
00180 
00181   Array< RealType >       m_ThreadSum;
00182   Array< RealType >       m_SumOfSquares;
00183   Array< SizeValueType >  m_Count;
00184   Array< PixelType >      m_ThreadMin;
00185   Array< PixelType >      m_ThreadMax;
00186 }; // end of class
00187 } // end namespace itk
00188 
00189 #ifndef ITK_MANUAL_INSTANTIATION
00190 #include "itkStatisticsImageFilter.hxx"
00191 #endif
00192 
00193 #endif
00194