Fit spectral components to a Vector of data from an image. More...

#include <ImageFit1D.h>

Public Types
enum	AbcissaType { PIXEL, IM_NATIVE, VELOCITY, N_TYPES }
Public Member Functions
	ImageFit1D ()
	Default Constructor.
	ImageFit1D (const ImageInterface< T > &image, uInt axis=0)
	Constructor.
	ImageFit1D (const ImageInterface< T > &image, const ImageInterface< T > &weights, uInt axis=0)
	Constructor with fitting weights image.
	~ImageFit1D ()
	Destructor.
	ImageFit1D (const ImageFit1D &other)
	Copy constructor.
ImageFit1D &	operator= (const ImageFit1D &other)
	Assignment operator.
void	setImage (const ImageInterface< T > &im, const ImageInterface< T > &weights, uInt axis)
	Set Image(s) and axis.
void	setImage (const ImageInterface< T > &im, uInt axis)
Bool	setData (const IPosition &pos, const ImageFit1D< T >::AbcissaType type, const Bool doAbs=True)
	Set the data to be fit.
Bool	setData (const ImageRegion &region, const ImageFit1D< T >::AbcissaType type, Bool doAbs=True)
void	setElements (const SpectralList &list)
	Set a SpectralList of SpectralElements to fit for.
void	addElement (const SpectralElement &el)
	Add new SpectralElement(s) to the SpectralList (can be empty) of SpectralElements to be fit for.
void	addElements (const SpectralList &list)
Bool	setGaussianElements (uInt nGauss)
	Set a SpectralList of Gaussian SpectralElements to fit for.
void	clearList ()
	Clear the SpectralList of elements to be fit for.
Bool	fit ()
	Do the fit and return convergence status.
Double	getChiSquared () const
	Get Chi Squared of fit.
Double	getNumberIterations () const
	Get number of iterations for last fit.
const SpectralList &	getList (Bool fit=True) const
	Recover the list of elements.
Vector< T >	getEstimate (Int which=-1) const
	Recover vectors for the estimate, fit and residual.
Vector< T >	getFit (Int which=-1) const
Vector< T >	getResidual (Int which=-1, Bool fit=True) const
Vector< Bool >	getTotalMask () const
	Get Total Mask (data and range mask)
Bool	succeeded () const
	did the fit succeed? should only be called after fit().
Bool	converged () const
	did the fit converge? should only be called after fit().
String	errorMessage () const
	Recover the error message.
void	invalidate ()
	flag the solution as invalid based on external criteria.
Bool	isValid () const
	is the solution valid? If False, some external logic has called invalidate()
Static Public Member Functions
static Bool	setAbcissaState (String &errMsg, ImageFit1D< T >::AbcissaType &type, CoordinateSystem &cSys, const String &xUnit, const String &doppler, uInt pixelAxis)
	Helper function.
Private Types
typedef NumericTraits< T > ::PrecisionType	FitterType
	In the future I will be able to template the fitter on T.
Private Member Functions
void	check () const
	Functions.
void	checkType () const
void	copy (const ImageFit1D< T > &other)
Bool	makeAbcissa (Vector< Double > &x, ImageFit1D< T >::AbcissaType type, Bool doAbs)
void	setWeightsImage (const ImageInterface< T > &im)
void	_resetFitter ()
	reset the fitter, for example if we've done a fit and want to move to the next position in the image
Private Attributes
std::auto_ptr< ImageInterface < T > >	itsImagePtr
std::auto_ptr< ImageInterface < T > >	itsWeightPtr
uInt	itsAxis
ProfileFit1D< FitterType >	itsFitter
CoordinateSystem	itsCS
String	itsError
Bool	_converged
Bool	_success
Bool	_isValid

Detailed Description

template<class T>
class casa::ImageFit1D< T >

Fit spectral components to a Vector of data from an image.

Intended use:

Public interface

Review Status

Test programs:: tImageFit1D

Prerequisite

Synopsis

Fit lists (held in class SpectralList) of SpectralElements to a Vector of data from the image. Each SpectralElement can be one from a variety of types. The values of the parameters for each SpectralElement provide the initial starting guesses for the fitting process.

You specify the domain in which the fit is to be done via the enum AbcissaType. The CoordinateSystem in the image is used to convert the pixel coordinates to the desired abcissa. You can change the units of the CoordinateSystem if you want to fit in different units. If you set an estimate yourself (function setElements or addElement) it is the callers responsibility that the elements are in the correct abcissa domain. Function setGaussianElements will automatically make an estimate in the correct domain.

Also, a SpectralElement object holds a mask indicating whether a parameter should be held fixed or solved for. After the fitting is done, a new SpectralList holding SpectralElements with the fitted parameters is created.

For all the functions that return a status Bool, True is good. If False is returned, an error message can be recovered with function errorMessage, You should not proceed if False is returned.

Exceptions will be thrown if you do not set the Image and axis via the constructor or setImage function.

Example

    PagedImage<Float> im("myimage");
    Int axis = 2;
    ImageFit1D<Float> fitter(image, axis);
    IPosition pos(in.ndim(),0);
    fitter.setData(pos, ImageFit1D<Float>::IM_NATIVE);     // Fit in native coordinate space
    fitter.setGaussianElements(3);                      // FIt 3 Gaussians
    if (fitter.fit()) {
       cerr << fitter.getList() << endl;                // Print result
    }

To Do

Add constraints

Definition at line 117 of file ImageFit1D.h.

Member Typedef Documentation

template<class T>

typedef NumericTraits<T>::PrecisionType casa::ImageFit1D< T >::FitterType [private]

In the future I will be able to template the fitter on T.

For now it must be Double.

Definition at line 264 of file ImageFit1D.h.

Member Enumeration Documentation

template<class T>

enum casa::ImageFit1D::AbcissaType

Enumerator:

PIXEL
IM_NATIVE
VELOCITY
N_TYPES

Definition at line 121 of file ImageFit1D.h.

Constructor & Destructor Documentation

template<class T>

casa::ImageFit1D< T >::ImageFit1D ( )

Default Constructor.

No image or axis set, so fit() cannot be called until other parameters are set, via eg setImage.

template<class T>

casa::ImageFit1D< T >::ImageFit1D	(	const ImageInterface< T > &	image,
		uInt	axis = `0`
	)

Constructor.

Fitting weights are assumed all unity.

template<class T>

casa::ImageFit1D< T >::ImageFit1D	(	const ImageInterface< T > &	image,
		const ImageInterface< T > &	weights,
		uInt	axis = `0`
	)

Constructor with fitting weights image.

The data and weights images must be the same shape.

template<class T>

casa::ImageFit1D< T >::~ImageFit1D ( )

Destructor.

template<class T>

casa::ImageFit1D< T >::ImageFit1D ( const ImageFit1D< T > & other )

Copy constructor.

Uses reference semantics.

Member Function Documentation

template<class T>

void casa::ImageFit1D< T >::_resetFitter ( ) [private]

reset the fitter, for example if we've done a fit and want to move to the next position in the image

template<class T>

void casa::ImageFit1D< T >::addElement ( const SpectralElement & el ) [inline]

Add new SpectralElement(s) to the SpectralList (can be empty) of SpectralElements to be fit for.

You must have already called setData to call this function.

Definition at line 187 of file ImageFit1D.h.

template<class T>

void casa::ImageFit1D< T >::addElements ( const SpectralList & list ) [inline]

Definition at line 188 of file ImageFit1D.h.

template<class T>

void casa::ImageFit1D< T >::check ( ) const [private]

Functions.

template<class T>

void casa::ImageFit1D< T >::checkType ( ) const [private]

template<class T>

void casa::ImageFit1D< T >::clearList ( ) [inline]

Clear the SpectralList of elements to be fit for.

Definition at line 200 of file ImageFit1D.h.

template<class T>

Bool casa::ImageFit1D< T >::converged ( ) const

did the fit converge? should only be called after fit().

template<class T>

void casa::ImageFit1D< T >::copy ( const ImageFit1D< T > & other ) [private]

template<class T>

String casa::ImageFit1D< T >::errorMessage ( ) const [inline]

Recover the error message.

Definition at line 240 of file ImageFit1D.h.

template<class T>

Bool casa::ImageFit1D< T >::fit ( )

Do the fit and return convergence status.

Errors in the fitting process will generate an AipsError exception and you should catch these yourself.

Referenced by casa::ImageFit1D< Float >::getList().

template<class T>

Double casa::ImageFit1D< T >::getChiSquared ( ) const [inline]

Get Chi Squared of fit.

Definition at line 208 of file ImageFit1D.h.

template<class T>

Vector<T> casa::ImageFit1D< T >::getEstimate ( Int which = -1 ) const

Recover vectors for the estimate, fit and residual.

If you don't specify which element, all elements are included If the Vectors are returned with zero length, it means an error condition exists (e.g. asking for fit before you do one). In this case an error message can be recovered with function errorMessage.

template<class T>

Vector<T> casa::ImageFit1D< T >::getFit ( Int which = -1 ) const

template<class T>

const SpectralList& casa::ImageFit1D< T >::getList ( Bool fit = True ) const [inline]

Recover the list of elements.

You can get the elements as initially estimated (fit=False), or after fitting (fit=True). In the latter case, the SpectralElements hold the parameters and errors of the fit.

Definition at line 217 of file ImageFit1D.h.

template<class T>

Double casa::ImageFit1D< T >::getNumberIterations ( ) const [inline]

Get number of iterations for last fit.

Definition at line 211 of file ImageFit1D.h.

template<class T>

Vector<T> casa::ImageFit1D< T >::getResidual	(	Int	which = `-1`,
		Bool	fit = `True`
	)		const

template<class T>

Vector<Bool> casa::ImageFit1D< T >::getTotalMask ( ) const [inline]

Get Total Mask (data and range mask)

Definition at line 231 of file ImageFit1D.h.

template<class T>

void casa::ImageFit1D< T >::invalidate ( )

flag the solution as invalid based on external criteria.

template<class T>

Bool casa::ImageFit1D< T >::isValid ( ) const

is the solution valid? If False, some external logic has called invalidate()

template<class T>

Bool casa::ImageFit1D< T >::makeAbcissa	(	Vector< Double > &	x,
		ImageFit1D< T >::AbcissaType	type,
		Bool	doAbs
	)		`[private]`

template<class T>

ImageFit1D& casa::ImageFit1D< T >::operator= ( const ImageFit1D< T > & other )

Assignment operator.

Uses reference semantics.

template<class T>

static Bool casa::ImageFit1D< T >::setAbcissaState	(	String &	errMsg,
		ImageFit1D< T >::AbcissaType &	type,
		CoordinateSystem &	cSys,
		const String &	xUnit,
		const String &	doppler,
		uInt	pixelAxis
	)		`[static]`

Helper function.

Sets up the CoordinateSystem to reflect the choice of abcissa unit and the doppler (if the axis is spectral).

template<class T>

Bool casa::ImageFit1D< T >::setData	(	const IPosition &	pos,
		const ImageFit1D< T >::AbcissaType	type,
		const Bool	doAbs = `True`
	)

Set the data to be fit.

All non-profile axes data are averaged. For the profile axis, the full spectrum is taken. The abcissa world values are computed when you call these functions. The domain of the abcissa values is controlled by AbcissaType and doAbs (absolute coordinates). The CoordinateSystem in the image is used to convert from pixels to world values.

template<class T>

Bool casa::ImageFit1D< T >::setData	(	const ImageRegion &	region,
		const ImageFit1D< T >::AbcissaType	type,
		Bool	doAbs = `True`
	)

template<class T>

void casa::ImageFit1D< T >::setElements ( const SpectralList & list ) [inline]

Set a SpectralList of SpectralElements to fit for.

These elements must be in the correct abcissa domain set in function setData. You must have already called setData to call this function. The SpectralElements in the list hold the initial estimates. They also contain the information about whether specific parameters are to be held fixed or allowed to vary in the fitting process. You can recover the list of elements with function getList.

Definition at line 181 of file ImageFit1D.h.

template<class T>

Bool casa::ImageFit1D< T >::setGaussianElements ( uInt nGauss )

Set a SpectralList of Gaussian SpectralElements to fit for.

The initial estimates for the Gaussians will be automatically determined in the correct abcissa domain. All of the parameters created by this function will be solved for by default. You can recover the list of elements with function getList. Status is returned, if False, error message can be recovered with errorMessage

template<class T>