[Insight-users] Submission: PCA decomposition calculator
Zachary Pincus
zpincus at stanford.edu
Mon, 3 May 2004 00:15:28 -0700
Hello,
I just updated/rebuilt my ITK tree from the CVS and built the test code
on my machine. Unfortunately, everything works fine on my box, so I'm
in a bad position to help debug the error.
I'm using GCC 3.3 (using c++ as the compiler) under OS X 10.3.3 on a G4
(7450).
Here is the output I get from the test (after renaming the test routine
"main" so I could build it separately.)
> [atlotl:ITK/Projects/PCADecomposer] zpincus%
> ./itkImagePCADecompositionCalculatorTest
> ImagePCADecompositionCalculator (0x2e009b0)
> RTTI typeinfo:
> N3itk31ImagePCADecompositionCalculatorINS_5ImageIdLj2EEES2_EE
> Reference Count: 1
> Modified Time: 25
> Debug: Off
> Observers:
> none
> Projection: 3.5574 -2.3118
> Image: 0x2e00800
> The basis of projection is:
> [-0.3853 0.5929 0.5929 -0.3853 ]
> [-0.5929 -0.3853 -0.3853 -0.5929 ]
> The projection of [0 2 2 0] is [-1.5412 -2.3716]
> this should be approx [-1.5412 -2.3716]
> The projection of [0 3 3 0] is [3.5574 -2.3118]
> this should be approx [3.5574 -2.3119]
>
> Test Passed
What platform(s) does the test fail on? What's the output? (The
dashboard seems to be down, else I'd just check there.)
Thanks,
Zach
On May 2, 2004, at 6:53 PM, Luis Ibanez wrote:
>
> Hi Zach,
>
> Thanks for contributing this filter to ITK.
>
> The files have been checked in. Some minor
> changes in coding Style were made.
>
> The test is executing but reporting a failure.
> If you have a chance could you please take a look
> at the values reported at the end ?
>
> It seems that the output doesn't corresponds to
> the range of values reported by Matlab.
>
>
> Thanks a lot,
>
>
> Luis
>
>
>
> ---------------------
> Zachary Pincus wrote:
>
>> I have written a calculator class that projects images down to a
>> space spanned by a specified orthonormal basis set and returns that
>> projection. (Also test code too!)
>> This task is generically useful, but it is especially used in image
>> processing for computing the "PCA decomposition coefficients" of
>> images, given a set of some number of principal components. (As a
>> reminder, the principal components are an orthonormal basis set
>> spanning some subspace of the original image space, the
>> "decomposition coefficients" vector is the projection of the image
>> down to this subspace.)
>> People tend to use PCA decompositions for image classification tasks,
>> because a decomposition captures how an image varies along
>> "important" dimensions, and disregards noise. (cf. the "eigenfaces"
>> literature for face recognition.)
>> Hopefully this code will prove useful and merits inclusion in ITK.
>> Zach Pincus
>> Department of Biochemistry and Program in Biomedical Informatics
>> Stanford University School of Medicine
>> ----------------------------------------------------------------------
>> --
>> #ifndef __itkImagePCADecompositionCalculator_h
>> #define __itkImagePCADecompositionCalculator_h
>> #include "itkObject.h"
>> #include "itkImagePCAShapeModelEstimator.h"
>> #include "vnl/vnl_vector.h"
>> #include "vnl/vnl_matrix.h"
>> namespace itk
>> {
>> /** \class ImagePCADecompositionCalculator
>> * \brief Decomposes an image into directions along basis components.
>> * * This calculator computes the projection of an image into a
>> subspace specified
>> * by some orthonormal basis.
>> * Typically, this basis will be the principal components of an image
>> data set,
>> * as calculated by an ImagePCAShapeModelEstimator. The output of the
>> calculator
>> * is a vnl_vector containing the coefficients along each dimension
>> of the
>> * provided basis set.
>> * To use this calculator, first set each basis image with the
>> SetBasisImage * method. In the PCA case, the basis images are the
>> outputs of the * ImagePCAShapeModelEstimator (except the zeroth
>> output, which is the average * image).
>> * SetBasisFromModel is a convenience method to set all of this
>> information from
>> * a given ImagePCAShapeModelEstimator instance.
>> *
>> * This class is templated over the input image type and the type of
>> images
>> * used to describe the basis.
>> *
>> * \warning This method assumes that the input image consists of
>> scalar pixel
>> * types.
>> *
>> * \warning All images (input and basis) must be the same size. This
>> is not * checked at runtime. *
>> * \ingroup Operators
>> */
>> template <class TInputImage, class TBasisImage = Image<double,
>> TInputImage::ImageDimension>
>> >//::itk::GetImageDimension<TInputImage>::ImageDimension> >
>> class ITK_EXPORT ImagePCADecompositionCalculator : public Object {
>> public:
>> /** Standard class typedefs. */
>> typedef ImagePCADecompositionCalculator Self;
>> typedef Object Superclass;
>> typedef SmartPointer<Self> Pointer;
>> typedef SmartPointer<const Self> ConstPointer;
>> /** Method for creation through the object factory. */
>> itkNewMacro(Self);
>> /** Run-time type information (and related methods). */
>> itkTypeMacro(ImagePCADecompositionCalculator, Object);
>> /** Type definitions for the input images. */
>> typedef TInputImage InputImageType;
>> typedef TBasisImage BasisImageType;
>> /** Pointer types for the image. */
>> typedef typename TInputImage::Pointer InputImagePointer;
>> typedef typename TBasisImage::Pointer BasisImagePointer;
>> /** Const Pointer type for the image. */
>> typedef typename TInputImage::ConstPointer InputImageConstPointer;
>> /** Basis image pixel type: this is also the type of the optput
>> vector */
>> typedef typename TBasisImage::PixelType BasisPixelType;
>> /** Input Image dimension */
>> itkStaticConstMacro(InputImageDimension, unsigned int,
>> TInputImage::ImageDimension); /** Basis Image
>> dimension */
>> itkStaticConstMacro(BasisImageDimension, unsigned int,
>> TBasisImage::ImageDimension); /** Vector
>> of basis image pointers. */
>> typedef std::vector< BasisImagePointer > BasisImagePointerVector;
>> /** Type definitions for internal vectors and matrices */
>> typedef vnl_matrix<BasisPixelType> BasisMatrixType;
>> typedef vnl_vector<BasisPixelType> BasisVectorType;
>> /** Set the input image. */
>> itkSetConstObjectMacro(Image,InputImageType);
>> /** Set the basis images. */
>> void SetBasisImages(const BasisImagePointerVector _arg); /**
>> Get the basis images. */
>> BasisImagePointerVector& GetBasisImages(void) {return
>> m_BasisImages;}
>> /** Type definition of a compatible ImagePCAShapeModelEstimator */
>> typedef typename ImagePCAShapeModelEstimator<TInputImage,
>> TBasisImage>::Pointer ModelPointerType;
>> /** Set the basis images from a ImagePCAShapeModelEstimator */
>> void SetBasisFromModel(ModelPointerType model);
>> /** Compute the PCA decomposition of the input image. */
>> void Compute(void);
>> /** Return the projection of the image. */
>> itkGetMacro(Projection,BasisVectorType);
>> protected:
>> ImagePCADecompositionCalculator();
>> virtual ~ImagePCADecompositionCalculator() {};
>> void PrintSelf(std::ostream& os, Indent indent) const;
>> void CalculateBasisMatrix(void);
>> void CalculateImageAsVector(void);
>> private:
>> typedef typename BasisImageType::SizeType BasisSizeType;
>> ImagePCADecompositionCalculator(const Self&); //purposely not
>> implemented
>> void operator=(const Self&); //purposely not implemented
>> BasisVectorType m_Projection;
>> BasisVectorType m_ImageAsVector;
>> BasisImagePointerVector m_BasisImages;
>> BasisSizeType m_Size;
>> InputImageConstPointer m_Image;
>> BasisMatrixType m_BasisMatrix;
>> bool m_BasisMatrixCalculated;
>> unsigned long m_NumPixels;
>> };
>> } // end namespace itk
>> #ifndef ITK_MANUAL_INSTANTIATION
>> #include "itkImagePCADecompositionCalculator.txx"
>> #endif
>> #endif
>
>
>
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users
>