Fixes

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

diff --git a/STEER/AliMatrixSparse.h b/STEER/AliMatrixSparse.h
index 20874a58eaafac06dc48365b971507e29370aca0..12898b4a4e4951fa2682627356c390a6501d2678 100644 (file)
--- a/STEER/AliMatrixSparse.h
+++ b/STEER/AliMatrixSparse.h
@@ -2,61 +2,17 @@
 #define ALIMATRIXSPARSE_H
 
 #include "AliMatrixSq.h"
+#include "AliVectorSparse.h"
 
 
-///////////////////////////////////////////////////////////////////////////////////////
-class AliVectorSparse {
- public:
-  AliVectorSparse();
-  AliVectorSparse(const AliVectorSparse& src);
-  virtual ~AliVectorSparse() {Clear();}
-  virtual void Print(Option_t* option="")                 const;
-  //
-  Int_t     GetNElems()                                   const {return fNElems;}
-  UShort_t *GetIndices()                                  const {return fIndex;}
-  Double_t *GetElems()                                    const {return fElems;}
-  UShort_t& GetIndex(Int_t i)                                   {return fIndex[i];}
-  Double_t& GetElem(Int_t i)                              const {return fElems[i];}
-  void      Clear();
-  void      Reset()                                             {memset(fElems,0,fNElems*sizeof(Double_t));}
-  void      ReSize(Int_t sz,Bool_t copy=kFALSE);
-  void      SortIndices(Bool_t valuesToo=kFALSE);
-  //
-  AliVectorSparse& operator=(const AliVectorSparse& src);
-  //
-  virtual Double_t         operator()(Int_t ind)         const;
-  virtual Double_t&        operator()(Int_t ind);
-  virtual void             SetToZero(Int_t ind);
-  Double_t                 FindIndex(Int_t ind)          const;
-  Double_t&                FindIndexAdd(Int_t ind);
-  //
-  Int_t     GetLastIndex()                               const {return fIndex[fNElems-1];}
-  Double_t  GetLastElem()                                const {return fElems[fNElems-1];}
-  Double_t &GetLastElem()                                      {return fElems[fNElems-1];}
-  //
- protected:
-  Int_t            fNElems;   // 
-  UShort_t*        fIndex;    // Index of stored elems
-  Double_t*        fElems;    // pointer on elements
-};
-
-
-//___________________________________________________
-inline Double_t AliVectorSparse::operator()(Int_t ind) const
-{
-  return FindIndex(ind);
-}
+/**********************************************************************************************/
+/* Sparse matrix class, used as a global matrix for AliMillePede2                             */
+/*                                                                                            */ 
+/* Author: ruben.shahoyan@cern.ch                                                             */
+/*                                                                                            */ 
+/**********************************************************************************************/
 
-//___________________________________________________
-inline Double_t& AliVectorSparse::operator()(Int_t ind)
-{
-  return FindIndexAdd(ind);
-}
-
-//////////////////////////////////////////////////////////////////////////////////////
 
-/////////////////////////////////////////////////////////////////////////////////////////////
-// Sparse matrix class
 //
 class AliMatrixSparse : public AliMatrixSq 
 {
@@ -68,7 +24,11 @@ class AliMatrixSparse : public AliMatrixSq
   //
   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;
@@ -82,8 +42,10 @@ class AliMatrixSparse : public AliMatrixSq
   Double_t&        DiagElem(Int_t r);
   void             SortIndices(Bool_t valuesToo=kFALSE);
   //
-  void MultiplyByVec(Double_t* vecIn, Double_t* vecOut) const;
-  void MultiplyByVec(TVectorD &vecIn, TVectorD &vecOut) const {MultiplyByVec((Double_t*)vecIn.GetMatrixArray(),(Double_t*)vecOut.GetMatrixArray());}
+  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:
   //
@@ -92,6 +54,12 @@ class AliMatrixSparse : public AliMatrixSq
   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)
 {
@@ -117,6 +85,7 @@ 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);
 }
 
@@ -134,11 +103,15 @@ inline Double_t AliMatrixSparse::DiagElem(Int_t row) const
 inline Double_t &AliMatrixSparse::DiagElem(Int_t row)
 {
   AliVectorSparse* rowv = GetRowAdd(row);
-  if (IsSymmetric()) return (rowv->GetNElems()>0 && rowv->GetLastIndex()==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