[u/mrichter/AliRoot.git] / STEER / AliMatrixSparse.h

#ifndef ALIMATRIXSPARSE_H
#define ALIMATRIXSPARSE_H

#include "AliMatrixSq.h"
#include "AliVectorSparse.h"


/**********************************************************************************************/
/* Sparse matrix class, used as a global matrix for AliMillePede2                             */
/*                                                                                            */ 
/* Author: ruben.shahoyan@cern.ch                                                             */
/*                                                                                            */ 
/**********************************************************************************************/


//
class AliMatrixSparse : public AliMatrixSq 
{
 public:
  AliMatrixSparse() : fVecs(0) {}
  AliMatrixSparse(Int_t size);
  AliMatrixSparse(const AliMatrixSparse &mat);
  virtual ~AliMatrixSparse() {Clear();}
  //
  AliVectorSparse* GetRow(Int_t ir)        const {return (ir<fNcols) ? fVecs[ir] : 0;}
  AliVectorSparse* GetRowAdd(Int_t ir);
  //
  virtual Int_t   GetSize()                            const {return fNrows;}
  virtual Int_t   GetNRows()                           const {return fNrows;}
  virtual Int_t   GetNCols()                           const {return fNcols;}
  //
  void Clear(Option_t* option="");
  void Reset()                            {for (int i=fNcols;i--;) GetRow(i)->Reset();}
  void Print(Option_t* option="")                      const;
  AliMatrixSparse& operator=(const AliMatrixSparse& src);
  Double_t&        operator()(Int_t row,Int_t col);
  Double_t         operator()(Int_t row,Int_t col)     const;
  void             SetToZero(Int_t row,Int_t col);
  Float_t          GetDensity()                        const;
  //
  Double_t         DiagElem(Int_t r)                   const;
  Double_t&        DiagElem(Int_t r);
  void             SortIndices(Bool_t valuesToo=kFALSE);
  //
  void MultiplyByVec(const TVectorD &vecIn, TVectorD &vecOut) const; 
  void MultiplyByVec(const Double_t* vecIn, Double_t* vecOut) const;
  //
  void AddToRow(Int_t r, Double_t *valc,Int_t *indc,Int_t n);
  //
 protected:
  //
  AliVectorSparse** fVecs;
  //
  ClassDef(AliMatrixSparse,0)
};

//___________________________________________________
inline void AliMatrixSparse::MultiplyByVec(const TVectorD &vecIn, TVectorD &vecOut) const 
{
  MultiplyByVec((Double_t*)vecIn.GetMatrixArray(),(Double_t*)vecOut.GetMatrixArray());
}

//___________________________________________________
inline void AliMatrixSparse::SetToZero(Int_t row,Int_t col)
{
  //  set existing element to 0
  if (IsSymmetric() && col>row) Swap(row,col); 
  AliVectorSparse* rowv = GetRow(row);
  if (rowv) rowv->SetToZero(col);
}

//___________________________________________________
inline Double_t AliMatrixSparse::operator()(Int_t row,Int_t col) const
{
  //  printf("M: find\n");
  if (IsSymmetric() && col>row) Swap(row,col); 
  AliVectorSparse* rowv = GetRow(row);
  if (!rowv) return 0;
  return rowv->FindIndex(col);
}

//___________________________________________________
inline Double_t& AliMatrixSparse::operator()(Int_t row,Int_t col)
{
  //  printf("M: findindexAdd\n");
  if (IsSymmetric() && col>row) Swap(row,col); 
  AliVectorSparse* rowv = GetRowAdd(row);
  if (col>=fNcols) fNcols = col+1;
  return rowv->FindIndexAdd(col);
}

//___________________________________________________
inline Double_t AliMatrixSparse::DiagElem(Int_t row) const
{
  AliVectorSparse* rowv = GetRow(row);
  if (!rowv) return 0;
  if (IsSymmetric()) return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ? rowv->GetLastElem() : 0.;
  else return rowv->FindIndex(row);
  //
}

//___________________________________________________
inline Double_t &AliMatrixSparse::DiagElem(Int_t row)
{
  AliVectorSparse* rowv = GetRowAdd(row);
  if (row>=fNcols) fNcols = row+1;
  if (IsSymmetric()) {
    return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ? 
		       rowv->GetLastElem() : rowv->FindIndexAdd(row);
  }
  else return rowv->FindIndexAdd(row);
  //
}


#endif
Commit	Line	Data
	1	#ifndef ALIMATRIXSPARSE_H
	2	#define ALIMATRIXSPARSE_H
	3
	4	#include "AliMatrixSq.h"
	5	#include "AliVectorSparse.h"
	6
	7
	8	/**********************************************************************************************/
	9	/* Sparse matrix class, used as a global matrix for AliMillePede2 */
	10	/* */
	11	/* Author: ruben.shahoyan@cern.ch */
	12	/* */
	13	/**********************************************************************************************/
	14
	15
	16	//
	17	class AliMatrixSparse : public AliMatrixSq
	18	{
	19	public:
	20	AliMatrixSparse() : fVecs(0) {}
	21	AliMatrixSparse(Int_t size);
	22	AliMatrixSparse(const AliMatrixSparse &mat);
	23	virtual ~AliMatrixSparse() {Clear();}
	24	//
	25	AliVectorSparse* GetRow(Int_t ir) const {return (ir<fNcols) ? fVecs[ir] : 0;}
	26	AliVectorSparse* GetRowAdd(Int_t ir);
	27	//
	28	virtual Int_t GetSize() const {return fNrows;}
	29	virtual Int_t GetNRows() const {return fNrows;}
	30	virtual Int_t GetNCols() const {return fNcols;}
	31	//
	32	void Clear(Option_t* option="");
	33	void Reset() {for (int i=fNcols;i--;) GetRow(i)->Reset();}
	34	void Print(Option_t* option="") const;
	35	AliMatrixSparse& operator=(const AliMatrixSparse& src);
	36	Double_t& operator()(Int_t row,Int_t col);
	37	Double_t operator()(Int_t row,Int_t col) const;
	38	void SetToZero(Int_t row,Int_t col);
	39	Float_t GetDensity() const;
	40	//
	41	Double_t DiagElem(Int_t r) const;
	42	Double_t& DiagElem(Int_t r);
	43	void SortIndices(Bool_t valuesToo=kFALSE);
	44	//
	45	void MultiplyByVec(const TVectorD &vecIn, TVectorD &vecOut) const;
	46	void MultiplyByVec(const Double_t* vecIn, Double_t* vecOut) const;
	47	//
	48	void AddToRow(Int_t r, Double_t valc,Int_t indc,Int_t n);
	49	//
	50	protected:
	51	//
	52	AliVectorSparse** fVecs;
	53	//
	54	ClassDef(AliMatrixSparse,0)
	55	};
	56
	57	//___________________________________________________
	58	inline void AliMatrixSparse::MultiplyByVec(const TVectorD &vecIn, TVectorD &vecOut) const
	59	{
	60	MultiplyByVec((Double_t)vecIn.GetMatrixArray(),(Double_t)vecOut.GetMatrixArray());
	61	}
	62
	63	//___________________________________________________
	64	inline void AliMatrixSparse::SetToZero(Int_t row,Int_t col)
	65	{
	66	// set existing element to 0
	67	if (IsSymmetric() && col>row) Swap(row,col);
	68	AliVectorSparse* rowv = GetRow(row);
	69	if (rowv) rowv->SetToZero(col);
	70	}
	71
	72	//___________________________________________________
	73	inline Double_t AliMatrixSparse::operator()(Int_t row,Int_t col) const
	74	{
	75	// printf("M: find\n");
	76	if (IsSymmetric() && col>row) Swap(row,col);
	77	AliVectorSparse* rowv = GetRow(row);
	78	if (!rowv) return 0;
	79	return rowv->FindIndex(col);
	80	}
	81
	82	//___________________________________________________
	83	inline Double_t& AliMatrixSparse::operator()(Int_t row,Int_t col)
	84	{
	85	// printf("M: findindexAdd\n");
	86	if (IsSymmetric() && col>row) Swap(row,col);
	87	AliVectorSparse* rowv = GetRowAdd(row);
	88	if (col>=fNcols) fNcols = col+1;
	89	return rowv->FindIndexAdd(col);
	90	}
	91
	92	//___________________________________________________
	93	inline Double_t AliMatrixSparse::DiagElem(Int_t row) const
	94	{
	95	AliVectorSparse* rowv = GetRow(row);
	96	if (!rowv) return 0;
	97	if (IsSymmetric()) return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ? rowv->GetLastElem() : 0.;
	98	else return rowv->FindIndex(row);
	99	//
	100	}
	101
	102	//___________________________________________________
	103	inline Double_t &AliMatrixSparse::DiagElem(Int_t row)
	104	{
	105	AliVectorSparse* rowv = GetRowAdd(row);
	106	if (row>=fNcols) fNcols = row+1;
	107	if (IsSymmetric()) {
	108	return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ?
	109	rowv->GetLastElem() : rowv->FindIndexAdd(row);
	110	}
	111	else return rowv->FindIndexAdd(row);
	112	//
	113	}
	114
	115
	116	#endif
	117