[Insight-users] ItkHistogram errors and proposing changes to
BinMin/Max
Martin Styner
martin_styner@ieee.org
Mon, 09 Sep 2002 13:13:47 +0200
This is a multi-part message in MIME format.
--Boundary_(ID_eS91NU/zzxjmm+gQRCuRPg)
Content-type: text/plain; charset=us-ascii; format=flowed
Content-transfer-encoding: 7BIT
Hi
My earlier changes are not enough, one also needs to adjust some
'interval 'variables used in the code of itkHistogram.txx. I have
attached the changed files.
best wishes
Martin
Martin Styner wrote:
> Hi
>
> I found a (minor) error in itkHistogram.txx (mispelling of a member at
> line 341), here is the diff of that specific location:
> 341c341
> < m_TempMeasurmentVector[i] = this->GetBinMin(i, index[i]) ;
> ---
> > m_TempMeasurementVector[i] = this->GetBinMin(i, index[i]) ;
>
>
> Also, I think that the Min and Max values for the Bin's should not be of
> type float, but of the same types as the histogram values
> (Statistics::Histogram<..>::MeasurementType). On one hand, this would
> eliminate all the conversion froma and to the float type, also it would
> make more sense from a design-viewpoint. This would create following
> changes to itkHistogram.txx and itkHistogram.h (again as diffs):
>
> Index: Code/Numerics/Statistics/itkHistogram.h
> diff -r1.20 itkHistogram.h
> 176c176
> < const float min)
> ---
> > const MeasurementType min)
> 181c181
> < unsigned long nbin, float max)
> ---
> > unsigned long nbin, const MeasurementType max)
> 212c212
> < &measurement) const ;
> ---
> > &measurement) ;
> 216c216
> < &measurement) const;
> ---
> > &measurement) ;
> 219c219
> < MeasurementVectorType& GetHistogramMinFromIndex(const IndexType
> &index) const ;
> ---
> > MeasurementVectorType& GetHistogramMinFromIndex(const IndexType
> &index) ;
> 222c222
> < MeasurementVectorType& GetHistogramMaxFromIndex(const IndexType
> &index) const ;
> ---
> > MeasurementVectorType& GetHistogramMaxFromIndex(const IndexType
> &index) ;
> 387c387
> < std::vector< std::vector<float> > m_Min ;
> ---
> > std::vector< std::vector<MeasurementType> > m_Min ;
> 390c390
> < std::vector< std::vector<float> > m_Max ;
> ---
> > std::vector< std::vector<MeasurementType> > m_Max ;
>
> Index: Code/Numerics/Statistics/itkHistogram.txx
> diff -r1.14 itkHistogram.txx
> 308c308
> < ::GetHistogramMinFromValue(const MeasurementVectorType &measurement)
> const
> ---
> > ::GetHistogramMinFromValue(const MeasurementVectorType &measurement)
> 322c322
> < ::GetHistogramMaxFromValue(const MeasurementVectorType &measurement)
> const
> ---
> > ::GetHistogramMaxFromValue(const MeasurementVectorType &measurement)
> 337c337
> < ::GetHistogramMinFromIndex(const IndexType &index) const
> ---
> > ::GetHistogramMinFromIndex(const IndexType &index)
> 341c341
> < m_TempMeasurmentVector[i] = this->GetBinMin(i, index[i]) ;
> ---
> > m_TempMeasurementVector[i] = this->GetBinMin(i, index[i]) ;
> 343c343
> < return m_TempMeasurmentVector ;
> ---
> > return m_TempMeasurementVector ;
> 351c351
> < ::GetHistogramMaxFromIndex(const IndexType &index) const
> ---
> > ::GetHistogramMaxFromIndex(const IndexType &index)
>
> What do others think of the proposed changes
> best wishes
> Martin
>
--
Martin Styner, PhD. Ing. ETH
Group Head Medical Image Analysis for Orthopedics
M.E. Mueller Institute for Biomechanics
Center for Computed Assisted Surgery
University of Bern
Murtenstrasse 35
P.O.Box 30
CH - 3010 Bern
Switzerland
Tel office: ++41-31-632-8784 , FAX: ++41-31-632-4951
email: Martin.Styner@memot.unibe.ch, martin_styner@ieee.org
WWW: http://cranium.unibe.ch/~mstyner
--Boundary_(ID_eS91NU/zzxjmm+gQRCuRPg)
Content-type: text/plain; name=itkHistogram.h
Content-transfer-encoding: 7BIT
Content-disposition: inline; filename=itkHistogram.h
/*=========================================================================
Program: Insight Segmentation & Registration Toolkit
Module: $RCSfile: itkHistogram.h,v $
Language: C++
Date: $Date: 2002/09/03 16:59:01 $
Version: $Revision: 1.20 $
Copyright (c) 2002 Insight Consortium. All rights reserved.
See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notices for more information.
=========================================================================*/
#ifndef __itkHistogram_h
#define __itkHistogram_h
#include <vector>
#include "itkIndex.h"
#include "itkSize.h"
#include "itkFixedArray.h"
#include "itkSample.h"
#include "itkDenseFrequencyContainer.h"
#include "itkSparseFrequencyContainer.h"
namespace itk{
namespace Statistics{
/** \class Histogram
* \brief This class stores measurement vectors in the context of n-dimensional
* histogram.
*
* Users can set arbitrary value for each bins min max value for each
* dimension (variable interval). Each dimension of the histogram represents
* each dimension of measurement vectors. For example, an Image with each pixel
* having a gray level intensity value and a gradient magnitude can be imported
* to this class. Then the resulting Histogram has two dimensions, intensity
* and gradient magnitude.
*
* Before any operation, users have to call Initialize(SizeType) method to
* prepare the indexing mechanism and the internal frequency container.
* After this, users want to set range of each bin using
* SetBinMin(dimension, n) and SetBinMax(dimension, n) methods.
*
* The first two template arguments are same as those of
* Sample. The last one, "TFrequencyContainter", is
* the type of the internal frequency container. If you think your Histogram
* is dense, in other words, almost every bin is used, then use default.
* If you expect that a very little portion of bins will be used, replace it
* with SparseFrequencyContainer class
*
* Since this class is n-dimensional, it supports data access
* methods using ITK Index type in addition to the methods using
* "InstanceIdentifiers".
*
* \sa Sample, DenseFrequencyContainer,
* SparseFrequencyContainer
*/
template < class TMeasurement = float, unsigned int VMeasurementVectorSize = 1,
class TFrequencyContainer = DenseFrequencyContainer< float > >
class ITK_EXPORT Histogram
: public Sample < FixedArray< TMeasurement, VMeasurementVectorSize > >
{
public:
/** Standard typedefs */
typedef Histogram Self ;
typedef Sample< FixedArray< TMeasurement, VMeasurementVectorSize > > Superclass ;
typedef SmartPointer<Self> Pointer ;
/** Run-time type information (and related methods). */
itkTypeMacro(Histogram, Sample) ;
/** standard New() method support */
itkNewMacro(Self) ;
/** Dimension of a measurement vector */
enum { MeasurementVectorSize = VMeasurementVectorSize } ;
/** type of an element of a measurement vector */
typedef TMeasurement MeasurementType ;
/** Common sample class typedefs */
typedef typename Superclass::MeasurementVectorType MeasurementVectorType ;
typedef typename Superclass::InstanceIdentifier InstanceIdentifier ;
typedef MeasurementVectorType ValueType ;
/** frequency container typedef */
typedef TFrequencyContainer FrequencyContainerType ;
typedef typename FrequencyContainerType::Pointer FrequencyContainerPointer ;
/** Frequency value type from superclass */
typedef typename FrequencyContainerType::FrequencyType FrequencyType ;
/** Index typedef support. An index is used to access pixel values. */
typedef itk::Index< VMeasurementVectorSize > IndexType;
typedef typename IndexType::IndexValueType IndexValueType;
/** size array type */
typedef itk::Size< VMeasurementVectorSize > SizeType ;
typedef typename SizeType::SizeValueType SizeValueType ;
/** bin min max value storage types */
typedef std::vector< MeasurementType > BinMinVectorType ;
typedef std::vector< MeasurementType > BinMaxVectorType ;
typedef std::vector< BinMinVectorType > BinMinContainerType ;
typedef std::vector< BinMaxVectorType > BinMaxContainerType ;
/** generates the offset table.
* subclasses should call this method in their initialize() method
* the overide methods have prepare the frequency container for
* input and output. */
void Initialize(const SizeType &size) ;
/** Do the same thing as the above Initialize(SizeType) method do
* , and also creates equal size bins within the range given
* by lower and upper bound. If users want to assign bin's min and
* max values along each dimension use SetBinMin() and SetBinMax()
* functions*/
void Initialize(const SizeType &size, MeasurementVectorType lowerBound,
MeasurementVectorType upperBound) ;
/** returns the index of histogram corresponding to measurement value */
IndexType& GetIndex(const MeasurementVectorType &measurement) ;
/** returns the index that is uniquely labelled by an instance identifier
* The corresponding id is the offset of the index
* This method uses ImageBase::ComputeIndex() method */
IndexType& GetIndex(const InstanceIdentifier &id) ;
/** returns true if the given index is out of bound meaning one of index
* is not between [0, last index] */
bool IsIndexOutOfBound(const IndexType &index) const;
/** returns the instance identifier of the cell that is indexed by the
* index. The corresponding instance identifier is the offset of the index
* This method uses ImageBase::ComputeIndex() method */
InstanceIdentifier GetInstanceIdentifier(const IndexType &index) const ;
/** Get the total number of instances in this container */
unsigned int GetNumberOfInstances() const ;
/** Returns the number of instances (bins or cells) in this container */
unsigned int Size() const ;
/** Returns the number of bins along the "dimension" */
unsigned int Size(const unsigned int &dimension) const
{ return static_cast< int >(m_Size[dimension]) ; }
/** Method to get m_Size */
SizeType GetSize() const
{ return m_Size ; }
/** return the size of each dimension of the measurement vector container */
SizeValueType GetSize(const unsigned int dimension) const
{
return m_Size[dimension] ;
}
/** Method to get minimum value of n th bin of dimension d */
MeasurementType& GetBinMin(const unsigned int dimension,
const unsigned long nbin)
{ return m_Min[dimension][nbin] ; }
/** Method to get maximum value of n th bin of dimension d */
MeasurementType& GetBinMax(const unsigned int dimension,
const unsigned long nbin)
{ return m_Max[dimension][nbin] ; }
/** Method to set minimum value of n th bin of dimension d */
void SetBinMin(const unsigned int dimension, const unsigned long nbin,
const MeasurementType min)
{ m_Min[dimension][nbin] = min ; }
/** Method to set maximum value of n th bin of dimension d */
void SetBinMax(const unsigned int dimension,
unsigned long nbin, const MeasurementType max)
{ m_Max[dimension][nbin] = max ; }
/** Method to get the minimum of the bin corresponding to the gray level of
* dimension d. */
MeasurementType& GetBinMinFromValue(const unsigned int dimension,
const float value ) const ;
/** Method to get the maximum of the bin corresponding to the gray level of
* dimension d. */
MeasurementType& GetBinMaxFromValue(const unsigned int dimension,
const float value ) const ;
/** Method to get the minimum vector of a dimension */
BinMinVectorType& GetDimensionMins(const unsigned int dimension) const
{ return m_Min[dimension] ; }
/** Method to get the maximum vector of a dimension */
BinMaxVectorType& GetDimensionMaxs(const unsigned int dimension) const
{ return m_Max[dimension] ; }
/** Method to get the minimum vector */
BinMinContainerType& GetMins() const
{ return m_Min ; }
/** Method to get the maximum vector */
BinMaxContainerType& GetMaxs() const
{ return m_Max ; }
/** Method to get mins of each dimension for a measurement in the histogram */
MeasurementVectorType& GetHistogramMinFromValue(const MeasurementVectorType
&measurement) ;
/** Method to get maxs of each dimension for a measurement in the histogram */
MeasurementVectorType& GetHistogramMaxFromValue(const MeasurementVectorType
&measurement) ;
/** Method to get mins in the histogram by index */
MeasurementVectorType& GetHistogramMinFromIndex(const IndexType &index) ;
/** Method to get maxs in the histogram by index */
MeasurementVectorType& GetHistogramMaxFromIndex(const IndexType &index) ;
/** Method to get the frequency from histogram */
FrequencyType GetFrequency(const InstanceIdentifier &id) const
{ return m_FrequencyContainer->GetFrequency(id) ; }
/** returns frequency of a bin that is indexed by index */
FrequencyType GetFrequency(const IndexType &index) const ;
/** Method to set the frequency of histogram */
void SetFrequency(const FrequencyType value) ;
/** Method to set the frequency of histogram */
void SetFrequency(const InstanceIdentifier &id, const FrequencyType value)
{ m_FrequencyContainer->SetFrequency(id, value) ; }
/** Method to set the frequency of histogram */
void SetFrequency(const IndexType &index,
const FrequencyType value) ;
/** Method to set the frequency corresponding to gray levels measurement */
void SetFrequency(const MeasurementVectorType &measurement,
const FrequencyType value) ;
/** Method to increase the frequency by one. This function is convinent
* to create histogram. */
void IncreaseFrequency(const InstanceIdentifier &id,
const FrequencyType value)
{ m_FrequencyContainer->IncreaseFrequency(id, value) ; }
/** Method to increase the frequency by one. This function is convinent
* to create histogram. */
void IncreaseFrequency(const IndexType &index,
const FrequencyType value) ;
/** Method to increase the frequency by one. This function is convinent
* to create histogram. */
void IncreaseFrequency(const MeasurementVectorType &measurement,
const FrequencyType value) ;
/** Method to get measurement from the histogram using an instance identifier */
MeasurementVectorType& GetMeasurementVector(const InstanceIdentifier &id) ;
/** Method to get measurement from the histogram */
MeasurementVectorType& GetMeasurementVector(const IndexType &index) ;
/** Method to get measurement from the histogram */
MeasurementType& GetMeasurement(const unsigned long n,
const unsigned int dimension) const ;
/** returns the frequency of the'dimension' dimension's 'n'th element. */
FrequencyType GetFrequency(const unsigned long n,
const unsigned int dimension) const ;
/** returns the frequency of the 'dimension' dimension */
FrequencyType GetTotalFrequency(const unsigned int &dimension) const ;
/** returns 'p'th percentile value.
* Let assume n = the index of the bin where the p-th percentile value is,
* min = min value of the dimension of the bin,
* max = max value of the dimension of the bin,
* interval = max - min ,
* pp = cumlated proportion until n-1 bin ;
* and pb = frequency of the bin / total frequency of the dimension.
*
* If p is less than 0.5,
* the percentile value =
* min + ((p - pp ) / pb) * interval
* If p is greater than or equal to 0.5
* the percentile value =
* max - ((pp - p) / pb) * interval */
double Quantile(const unsigned int dimension, const double &p) ;
/** iterator support */
class Iterator ;
friend class Iterator ;
Iterator Begin()
{
Iterator iter(0, this) ;
return iter ;
}
Iterator End()
{
return Iterator(m_OffsetTable[MeasurementVectorSize], this) ;
}
class Iterator
{
public:
Iterator(){};
Iterator(Pointer histogram)
{
m_Id = 0 ;
m_Histogram = histogram;
}
Iterator(InstanceIdentifier id, Pointer histogram)
: m_Id(id), m_Histogram(histogram)
{}
FrequencyType GetFrequency() const
{
return m_Histogram->GetFrequency(m_Id) ;
}
void SetFrequency(const FrequencyType value)
{
m_Histogram->SetFrequency(m_Id, value);
}
InstanceIdentifier GetInstanceIdentifier() const
{ return m_Id ; }
MeasurementVectorType& GetMeasurementVector() const
{
return m_Histogram->GetMeasurementVector(m_Id) ;
}
Iterator& operator++()
{
++m_Id;
return *this;
}
bool operator!=(const Iterator& it)
{ return (m_Id != it.m_Id); }
bool operator==(const Iterator& it)
{ return (m_Id == it.m_Id); }
Iterator& operator=(const Iterator& it)
{
m_Id = it.m_Id;
m_Histogram = it.m_Histogram ;
return *this ;
}
Iterator(const Iterator& it)
{
m_Id = it.m_Id;
m_Histogram = it.m_Histogram ;
}
private:
// Iterator pointing DenseFrequencyContainer
InstanceIdentifier m_Id;
// Pointer of DenseFrequencyContainer
Pointer m_Histogram ;
} ; // end of iterator class
protected:
Histogram() ;
virtual ~Histogram() {}
void PrintSelf(std::ostream& os, Indent indent) const;
// The number of bins for each dimension
SizeType m_Size ;
// lower bound of each bin
std::vector< std::vector<MeasurementType> > m_Min ;
// upper bound of each bin
std::vector< std::vector<MeasurementType> > m_Max ;
private:
Histogram(const Self&); //purposely not implemented
void operator=(const Self&); //purposely not implemented
InstanceIdentifier m_OffsetTable[MeasurementVectorSize + 1] ;
FrequencyContainerPointer m_FrequencyContainer ;
unsigned int m_NumberOfInstances ;
MeasurementVectorType m_TempMeasurementVector ;
IndexType m_TempIndex ;
} ; // end of class
} // end of namespace Statistics
} // end of namespace itk
#ifndef ITK_MANUAL_INSTANTIATION
#include "itkHistogram.txx"
#endif
#endif
--Boundary_(ID_eS91NU/zzxjmm+gQRCuRPg)
Content-type: text/plain; name=itkHistogram.txx
Content-transfer-encoding: 7BIT
Content-disposition: inline; filename=itkHistogram.txx
/*=========================================================================
Program: Insight Segmentation & Registration Toolkit
Module: $RCSfile: itkHistogram.txx,v $
Language: C++
Date: $Date: 2002/09/03 16:59:01 $
Version: $Revision: 1.14 $
Copyright (c) 2002 Insight Consortium. All rights reserved.
See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notices for more information.
=========================================================================*/
#ifndef _itkHistogram_txx
#define _itkHistogram_txx
#include "itkHistogram.h"
#include "itkNumericTraits.h"
namespace itk{
namespace Statistics{
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::Histogram()
{
m_FrequencyContainer = FrequencyContainerType::New() ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
unsigned int
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::Size() const
{
unsigned int size = 1 ;
for (unsigned int i = 0 ; i < MeasurementVectorSize ; i++)
{
size *= m_Size[i] ;
}
return size ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
void
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::Initialize(const SizeType &size)
{
m_Size = size ;
// creates offset table which will be used for generation of
// instance identifiers.
InstanceIdentifier num = 1 ;
m_OffsetTable[0] = num ;
for (unsigned int i = 0 ; i < MeasurementVectorSize ; i++)
{
num *= m_Size[i] ;
m_OffsetTable[i + 1] = num ;
}
m_NumberOfInstances = num ;
// adjust the sizes of min max value containers
int dim;
m_Min.resize(MeasurementVectorSize);
for ( dim = 0; dim < MeasurementVectorSize; dim++)
{
m_Min[dim].resize(m_Size[dim]);
}
m_Max.resize(MeasurementVectorSize);
for ( dim = 0; dim < MeasurementVectorSize; dim++)
{
m_Max[dim].resize(m_Size[dim]);
}
// initialize the frequency container
m_FrequencyContainer->Initialize(m_OffsetTable[VMeasurementVectorSize]) ;
this->SetFrequency(0.0f) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
void
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::Initialize(const SizeType &size, MeasurementVectorType lowerBound,
MeasurementVectorType upperBound)
{
this->Initialize(size) ;
float interval ;
for ( int i = 0 ; i < MeasurementVectorSize ; i++)
{
interval = (float) (upperBound[i] - lowerBound[i]) /
static_cast< MeasurementType >(size[i]) ;
// Set the min vector and max vector
for (unsigned int j = 0; j < (size[i] - 1) ; j++)
{
this->SetBinMin(i, j, lowerBound[i] + ((float)j * interval)) ;
this->SetBinMax(i, j, lowerBound[i] + (((float)j + 1) * interval));
}
this->SetBinMin(i, size[i] - 1,
lowerBound[i] + (((float) size[i] - 1) * interval)) ;
this->SetBinMax(i, size[i] - 1,
upperBound[i]) ;
}
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
inline typename Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>::IndexType&
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetIndex(const MeasurementVectorType &measurement)
{
// now using somthing similar to binary search to find
// index.
int dim ;
int begin, mid, end ;
MeasurementType median ;
MeasurementType tempMeasurement ;
for (dim = 0 ; dim < MeasurementVectorSize ; dim++)
{
tempMeasurement = measurement[dim] ;
begin = 0 ;
if (tempMeasurement < m_Min[dim][begin])
{
// one of measurement is below the minimum
m_TempIndex[0] = (long) m_Size[0] ;
return m_TempIndex ;
}
end = m_Min[dim].size() - 1 ;
if (tempMeasurement >= m_Max[dim][end])
{
// one of measurement is above the maximum
m_TempIndex[0] = (long) m_Size[0] ;
return m_TempIndex ;
}
mid = (end + 1) / 2 ;
median = m_Min[dim][mid] ;
while(true)
{
if (tempMeasurement < median )
{
end = mid - 1 ;
}
else if (tempMeasurement > median)
{
if (tempMeasurement < m_Max[dim][mid])
{
m_TempIndex[dim] = mid ;
break ;
}
begin = mid + 1 ;
}
else
{
// measurement[dim] = m_Min[dim][med]
m_TempIndex[dim] = mid ;
break ;
}
mid = begin + (end - begin) / 2 ;
median = m_Min[dim][mid] ;
} // end of while
} // end of for()
return m_TempIndex;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
inline typename Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>::IndexType&
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetIndex(const InstanceIdentifier &id)
{
InstanceIdentifier id2 = id ;
for (int i = MeasurementVectorSize - 1 ; i > 0 ; i--)
{
m_TempIndex[i] = static_cast<IndexValueType>(id2 / m_OffsetTable[i]);
id2 -= (m_TempIndex[i] * m_OffsetTable[i]);
}
m_TempIndex[0] = static_cast<IndexValueType>(id2);
return m_TempIndex;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline bool
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::IsIndexOutOfBound(const IndexType &index) const
{
for (int dim = 0 ; dim < MeasurementVectorSize ; dim++)
{
if (index[dim] < 0 || index[dim] >= static_cast<IndexValueType>(m_Size[dim]))
{
return true ;
}
}
return false ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram<TMeasurement, VMeasurementVectorSize,
TFrequencyContainer>::InstanceIdentifier
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetInstanceIdentifier(const IndexType &index) const
{
InstanceIdentifier id = 0 ;
for (int i= MeasurementVectorSize - 1 ; i > 0 ; i-- )
{
id += index[i] * m_OffsetTable[i];
}
id += index[0] ;
return id ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
unsigned int
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetNumberOfInstances() const
{
return m_NumberOfInstances ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram<TMeasurement, VMeasurementVectorSize,
TFrequencyContainer>::MeasurementType&
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetBinMinFromValue(const unsigned int dimension, const float value ) const
{
// If the value is lower than any of min value in the Histogram,
// it returns the lowest min value
if ( value <= this->m_Min[dimension][0] )
{
return this->m_Min[dimension][0];
}
// If the value is higher than any of min value in the Histogram,
// it returns the highest min value
if ( value >= m_Min[dimension][m_Size[dimension]-1] )
{
return m_Min[dimension][this->m_Size[dimension]-1];
}
for ( int i=0; i < this->m_Size[dimension]; i++ )
{
if ( (value >= this->m_Min[dimension][i])
&& (value < this->m_Max[dimension][i]) )
{
return this->m_Min[dimension][i];
}
}
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetBinMaxFromValue(const unsigned int dimension, const float value ) const
{
// If the value is lower than any of max value in the Histogram,
// it returns the lowest max value
if ( value <= this->m_Max[dimension][0] )
{
return this->m_Max[dimension][0];
}
// If the value is higher than any of max value in the Histogram,
// it returns the highest max value
if ( value >= m_Max[dimension][m_Size[dimension]-1] )
{
return m_Max[dimension][this->m_Size[dimension]-1];
}
for ( int i = 0 ; i < this->m_Size[dimension]; i++ )
{
if ( (value >= this->m_Min[dimension][i])
&& (value < this->m_Max[dimension][i]) )
{
return this->m_Max[dimension][i];
}
}
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetHistogramMinFromValue(const MeasurementVectorType &measurement)
{
for ( int i = 0; i < MeasurementVectorSize; i++ )
{
m_TempMeasurementVector[i] = this->GetDimensionMinByValue(i,measurement[i]);
}
return m_TempMeasurementVector ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram<TMeasurement, VMeasurementVectorSize, TFrequencyContainer>
::GetHistogramMaxFromValue(const MeasurementVectorType &measurement)
{
for ( int i=0; i < MeasurementVectorSize; i++ )
{
m_TempMeasurementVector[i] = this->GetDimensionMaxByValue(i,measurement[i]);
}
return m_TempMeasurementVector ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetHistogramMinFromIndex(const IndexType &index)
{
for ( int i=0; i < MeasurementVectorSize; i++ )
{
m_TempMeasurementVector[i] = this->GetBinMin(i, index[i]) ;
}
return m_TempMeasurementVector ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetHistogramMaxFromIndex(const IndexType &index)
{
for ( int i=0; i < MeasurementVectorSize; i++ )
{
m_TempMeasurementVector[i] = this->GetBinMax(i, index[i]) ;
}
return m_TempMeasurementVector ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetMeasurementVector(const IndexType &index)
{
for ( int i = 0; i < MeasurementVectorSize; i++)
{
m_TempMeasurementVector[i] = (m_Min[i][index[i]] + m_Max[i][index[i]])/2;
}
return m_TempMeasurementVector ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementVectorType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetMeasurementVector(const InstanceIdentifier &id)
{
return GetMeasurementVector(GetIndex(id)) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::SetFrequency(const FrequencyType value)
{
typename Self::Iterator iter = this->Begin() ;
typename Self::Iterator end = this->End() ;
while ( iter != end )
{
iter.SetFrequency(value) ;
++iter ;
}
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::SetFrequency(const IndexType &index, const FrequencyType value)
{
this->SetFrequency(GetInstanceIdentifier(index), value) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::SetFrequency(const MeasurementVectorType &measurement, const FrequencyType value)
{
this->SetFrequency(GetInstanceIdentifier(GetIndex(measurement)), value) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::IncreaseFrequency(const IndexType &index, const FrequencyType value)
{
this->IncreaseFrequency(GetInstanceIdentifier(index), value) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::IncreaseFrequency(const MeasurementVectorType &measurement, const FrequencyType value)
{
this->IncreaseFrequency(GetInstanceIdentifier(GetIndex(measurement)), value) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::FrequencyType
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetFrequency(const IndexType &index) const
{
return ( GetFrequency(GetInstanceIdentifier(index)) ) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer>
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::MeasurementType&
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetMeasurement(const unsigned long n, const unsigned int dimension) const
{
return static_cast< MeasurementType >((m_Min[dimension][n] +
m_Max[dimension][n]) / 2) ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::FrequencyType
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetFrequency(const unsigned long n, const unsigned int dimension) const
{
InstanceIdentifier nextOffset = m_OffsetTable[dimension + 1] ;
InstanceIdentifier current = m_OffsetTable[dimension] * n ;
InstanceIdentifier includeLength = m_OffsetTable[dimension] ;
InstanceIdentifier include ;
InstanceIdentifier includeEnd ;
InstanceIdentifier last = m_OffsetTable[VMeasurementVectorSize] ;
FrequencyType frequency = 0 ;
while (current < last)
{
include = current ;
includeEnd = include + includeLength ;
while(include < includeEnd)
{
frequency += GetFrequency(include) ;
include++ ;
}
current += nextOffset ;
}
return frequency ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
inline typename Histogram< TMeasurement, VMeasurementVectorSize,
TFrequencyContainer >::FrequencyType
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::GetTotalFrequency(const unsigned int &dimension) const
{
return m_FrequencyContainer->GetTotalFrequency() ;
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
double
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::Quantile(const unsigned int dimension, const double &p)
{
InstanceIdentifier n ;
const unsigned int size = this->Size(dimension) ;
double p_n_prev ;
double p_n ;
double f_n ;
double cumulated = 0 ;
double totalFrequency = double(GetTotalFrequency(dimension)) ;
double binProportion ;
double min, max, interval ;
if ( p < 0.5 )
{
n = 0 ;
p_n_prev = NumericTraits< double >::Zero ;
p_n = NumericTraits< double >::Zero ;
do
{
f_n = GetFrequency(n, dimension) ;
cumulated += f_n ;
p_n_prev = p_n ;
p_n = cumulated / totalFrequency ;
n++ ;
}
while( n < size && p_n < p) ;
binProportion = f_n / totalFrequency ;
min = double(GetBinMin(dimension, n - 1)) ;
max = double(GetBinMax(dimension, n - 1)) ;
interval = max - min ;
return min + ((p - p_n_prev) / binProportion) * interval ;
}
else
{
n = size - 1 ;
InstanceIdentifier m = NumericTraits< InstanceIdentifier >::Zero;
p_n_prev = NumericTraits< double >::One ;
p_n = NumericTraits< double >::One ;
do
{
f_n = GetFrequency(n, dimension) ;
cumulated += f_n ;
p_n_prev = p_n ;
p_n = NumericTraits< double >::One - cumulated / totalFrequency ;
n--;
m++;
}
while( m < size && p_n > p);
binProportion = f_n / totalFrequency ;
double min = double(GetBinMin(dimension, n + 1)) ;
double max = double(GetBinMax(dimension, n + 1)) ;
double interval = max - min ;
return max - ((p_n_prev - p) / binProportion) * interval ;
}
}
template< class TMeasurement, unsigned int VMeasurementVectorSize,
class TFrequencyContainer >
void
Histogram< TMeasurement, VMeasurementVectorSize, TFrequencyContainer >
::PrintSelf(std::ostream& os, Indent indent) const
{
Superclass::PrintSelf(os,indent);
os << indent << "OffsetTable: " << *m_OffsetTable << std::endl;
os << indent << "FrequencyContainerPointer: " << m_FrequencyContainer
<< std::endl;
}
} // end of namespace Statistics
} // end of namespace itk
#endif
--Boundary_(ID_eS91NU/zzxjmm+gQRCuRPg)--