Diana Software
QMatrixCArray.hh
Go to the documentation of this file.
1 
2 #ifndef _Q_MATRIX_C_ARRAY_HH_
3 #define _Q_MATRIX_C_ARRAY_HH_
4 
5 #include "QObject.hh"
6 #include "QVectorC.hh"
7 #include "QMatrixC.hh"
8 #include <vector>
9 #include <map>
10 
36 
37 public:
38 
41 
42  QMatrixCArray(const QMatrixCArray &in);
43 
44  virtual ~QMatrixCArray();
46  virtual void Clear();
47 
49  void Initialize(double Re,double Im=0);
50 
52  void Initialize(const QComplex &Val);
63  void Invert();
64 
74  void CholeskyDecompose() ;
75 
89  virtual void CholeskySolveX(const Diana::QVectorC &Y,Diana::QVectorC &X)
90  const;
91 
92 protected:
93 
94  QMatrixCArray(UInt_t nMatrix,UInt_t nRow,UInt_t nCol);
95 
96  const QMatrixCArray &operator=(const QMatrixCArray &in);
97 
99  void Resize(UInt_t nMatrix,UInt_t nRow,UInt_t nCol);
100 
102  UInt_t GetSize() const { return fSize; }
104  UInt_t GetNMatrix() const { return fNMatrix; }
105 
107  UInt_t GetNRow() const { return fNRow; }
108 
110  UInt_t GetNCol() const { return fNCol; }
111 
119  QMatrixC GetMatrix(UInt_t idx);
120 
130  const QMatrixC GetMatrix(UInt_t idx) const;
131 
140  Diana::QVectorC GetVector(UInt_t idx1,UInt_t idx2);
141 
152  const Diana::QVectorC GetVector(UInt_t idx1,UInt_t idx2) const;
153 
154 
155  Bool_t fDataOwner;
156 
158  UInt_t fSize;
159 
161  UInt_t fNMatrix;
162 
164  UInt_t fNRow;
165 
167  UInt_t fNCol;
168 
170  Double_t *fData; //[fSize]
171 
173 };
175 #endif
double Im(const Diana::QComplex &z)
Function to get the imag part.
Definition: QComplex.cc:262
double Re(const Diana::QComplex &z)
Function to get the real part.
Definition: QComplex.cc:260
#define Q_END_NAMESPACE
Definition: QDiana.hh:22
#define Q_BEGIN_NAMESPACE
Definition: QDiana.hh:20
An array of complex matrices that have been laid out in memory.
UInt_t fNMatrix
Number of identical matrices.
UInt_t fSize
Number of double values in fData.
Diana::QVectorC GetVector(UInt_t idx1, UInt_t idx2)
Return the vector of values in the (idx1,idx2) position.
Double_t * fData
Pointer to the head of the data.
void Resize(UInt_t nMatrix, UInt_t nRow, UInt_t nCol)
Resize, destroys all data in memory.
UInt_t fNRow
Number of rows in each matrix.
void CholeskyDecompose()
Cholesky decompose the matrices in place.
ClassDef(QMatrixCArray, 2)
UInt_t fNCol
Number of columns in each matrix.
UInt_t GetNRow() const
Return the number of rows in each matrix.
QMatrixCArray()
Default ctor.
void Initialize(double Re, double Im=0)
Initialize the values of the matrix.
void Invert()
Invert the matrices in place.
UInt_t GetSize() const
Return the total number of doubles in memory.
QMatrixC GetMatrix(UInt_t idx)
Return the vector of values in the (idx1,idx2) position.
virtual ~QMatrixCArray()
UInt_t GetNMatrix() const
Return the number of matrices stored.
UInt_t GetNCol() const
Return the number of columns in each matrix.
virtual void CholeskySolveX(const Diana::QVectorC &Y, Diana::QVectorC &X) const
Solve the equation AX=Y for X.
virtual void Clear()
Clear the data in memory.
const QMatrixCArray & operator=(const QMatrixCArray &in)
Interface for complex matrices in Diana analysis.
Definition: QMatrixC.hh:27