STEER/AliMatrixSparse.h

   1 #ifndef ALIMATRIXSPARSE_H
   2 #define ALIMATRIXSPARSE_H
   3
   4 #include "AliMatrixSq.h"
   5
   6
   7 ///////////////////////////////////////////////////////////////////////////////////////
   8 class AliVectorSparse {
   9  public:
  10   AliVectorSparse();
  11   AliVectorSparse(const AliVectorSparse& src);
  12   virtual ~AliVectorSparse() {Clear();}
  13   virtual void Print(Option_t* option="")                 const;
  14   //
  15   Int_t     GetNElems()                                   const {return fNElems;}
  16   UShort_t *GetIndices()                                  const {return fIndex;}
  17   Double_t *GetElems()                                    const {return fElems;}
  18   UShort_t& GetIndex(Int_t i)                                   {return fIndex[i];}
  19   Double_t& GetElem(Int_t i)                              const {return fElems[i];}
  20   void      Clear();
  21   void      Reset()                                             {memset(fElems,0,fNElems*sizeof(Double_t));}
  22   void      ReSize(Int_t sz,Bool_t copy=kFALSE);
  23   void      SortIndices(Bool_t valuesToo=kFALSE);
  24   //
  25   AliVectorSparse& operator=(const AliVectorSparse& src);
  26   //
  27   virtual Double_t         operator()(Int_t ind)         const;
  28   virtual Double_t&        operator()(Int_t ind);
  29   virtual void             SetToZero(Int_t ind);
  30   Double_t                 FindIndex(Int_t ind)          const;
  31   Double_t&                FindIndexAdd(Int_t ind);
  32   //
  33   Int_t     GetLastIndex()                               const {return fIndex[fNElems-1];}
  34   Double_t  GetLastElem()                                const {return fElems[fNElems-1];}
  35   Double_t &GetLastElem()                                      {return fElems[fNElems-1];}
  36   //
  37  protected:
  38   Int_t            fNElems;   //
  39   UShort_t*        fIndex;    // Index of stored elems
  40   Double_t*        fElems;    // pointer on elements
  41 };
  42
  43
  44 //___________________________________________________
  45 inline Double_t AliVectorSparse::operator()(Int_t ind) const
  46 {
  47   return FindIndex(ind);
  48 }
  49
  50 //___________________________________________________
  51 inline Double_t& AliVectorSparse::operator()(Int_t ind)
  52 {
  53   return FindIndexAdd(ind);
  54 }
  55
  56 //////////////////////////////////////////////////////////////////////////////////////
  57
  58 /////////////////////////////////////////////////////////////////////////////////////////////
  59 // Sparse matrix class
  60 //
  61 class AliMatrixSparse : public AliMatrixSq
  62 {
  63  public:
  64   AliMatrixSparse() : fVecs(0) {}
  65   AliMatrixSparse(Int_t size);
  66   AliMatrixSparse(const AliMatrixSparse &mat);
  67   virtual ~AliMatrixSparse() {Clear();}
  68   //
  69   AliVectorSparse* GetRow(Int_t ir)        const {return (ir<fNcols) ? fVecs[ir] : 0;}
  70   AliVectorSparse* GetRowAdd(Int_t ir);
  71
  72   void Clear(Option_t* option="");
  73   void Reset()                            {for (int i=fNcols;i--;) GetRow(i)->Reset();}
  74   void Print(Option_t* option="")                      const;
  75   AliMatrixSparse& operator=(const AliMatrixSparse& src);
  76   Double_t&        operator()(Int_t row,Int_t col);
  77   Double_t         operator()(Int_t row,Int_t col)     const;
  78   void             SetToZero(Int_t row,Int_t col);
  79   Float_t          GetDensity()                        const;
  80   //
  81   Double_t         DiagElem(Int_t r)                   const;
  82   Double_t&        DiagElem(Int_t r);
  83   void             SortIndices(Bool_t valuesToo=kFALSE);
  84   //
  85   void MultiplyByVec(Double_t* vecIn, Double_t* vecOut) const;
  86   void MultiplyByVec(TVectorD &vecIn, TVectorD &vecOut) const {MultiplyByVec((Double_t*)vecIn.GetMatrixArray(),(Double_t*)vecOut.GetMatrixArray());}
  87   //
  88  protected:
  89   //
  90   AliVectorSparse** fVecs;
  91   //
  92   ClassDef(AliMatrixSparse,0)
  93 };
  94
  95 //___________________________________________________
  96 inline void AliMatrixSparse::SetToZero(Int_t row,Int_t col)
  97 {
  98   //  set existing element to 0
  99   if (IsSymmetric() && col>row) Swap(row,col);
 100   AliVectorSparse* rowv = GetRow(row);
 101   if (rowv) rowv->SetToZero(col);
 102 }
 103
 104 //___________________________________________________
 105 inline Double_t AliMatrixSparse::operator()(Int_t row,Int_t col) const
 106 {
 107   //  printf("M: find\n");
 108   if (IsSymmetric() && col>row) Swap(row,col);
 109   AliVectorSparse* rowv = GetRow(row);
 110   if (!rowv) return 0;
 111   return rowv->FindIndex(col);
 112 }
 113
 114 //___________________________________________________
 115 inline Double_t& AliMatrixSparse::operator()(Int_t row,Int_t col)
 116 {
 117   //  printf("M: findindexAdd\n");
 118   if (IsSymmetric() && col>row) Swap(row,col);
 119   AliVectorSparse* rowv = GetRowAdd(row);
 120   return rowv->FindIndexAdd(col);
 121 }
 122
 123 //___________________________________________________
 124 inline Double_t AliMatrixSparse::DiagElem(Int_t row) const
 125 {
 126   AliVectorSparse* rowv = GetRow(row);
 127   if (!rowv) return 0;
 128   if (IsSymmetric()) return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ? rowv->GetLastElem() : 0.;
 129   else return rowv->FindIndex(row);
 130   //
 131 }
 132
 133 //___________________________________________________
 134 inline Double_t &AliMatrixSparse::DiagElem(Int_t row)
 135 {
 136   AliVectorSparse* rowv = GetRowAdd(row);
 137   if (IsSymmetric()) return (rowv->GetNElems()>0 && rowv->GetLastIndex()==row) ?
 138                        rowv->GetLastElem() : rowv->FindIndexAdd(row);
 139   else return rowv->FindIndexAdd(row);
 140   //
 141 }
 142
 143 #endif
 144