AliESDHeader: AliTriggerConfiguration and more trigger scalers added

[u/mrichter/AliRoot.git] / STEER / AliSymMatrix.cxx

diff --git a/STEER/AliSymMatrix.cxx b/STEER/AliSymMatrix.cxx
index 184c511f5b6ea02880e780f541e126c04b751d74..2c000cce0830047b6fdd7bd1f3244953ea231d37 100644 (file)
--- a/STEER/AliSymMatrix.cxx
+++ b/STEER/AliSymMatrix.cxx
@@ -14,6 +14,7 @@
 #include <stdio.h>
 #include <iostream>
 #include <float.h>
+#include <string.h>
 //
 #include <TClass.h>
 #include <TMath.h>
@@ -32,6 +33,7 @@ Int_t         AliSymMatrix::fgCopyCnt = 0;
 AliSymMatrix::AliSymMatrix() 
 : fElems(0),fElemsAdd(0)
 {
+  // default constructor
   fSymmetric = kTRUE;
   fgCopyCnt++;
 }
@@ -40,7 +42,7 @@ AliSymMatrix::AliSymMatrix()
 AliSymMatrix::AliSymMatrix(Int_t size)
   : AliMatrixSq(),fElems(0),fElemsAdd(0)
 {
-  //
+  //constructor for matrix with defined size
   fNrows = 0;
   fNrowIndex = fNcols = fRowLwb = size;
   fElems     = new Double_t[fNcols*(fNcols+1)/2];
@@ -54,6 +56,7 @@ AliSymMatrix::AliSymMatrix(Int_t size)
 AliSymMatrix::AliSymMatrix(const AliSymMatrix &src) 
   : AliMatrixSq(src),fElems(0),fElemsAdd(0)
 {
+  // copy constructor
   fNrowIndex = fNcols = src.GetSize();
   fNrows = 0;
   fRowLwb = src.GetSizeUsed();
@@ -88,7 +91,7 @@ AliSymMatrix::~AliSymMatrix()
 //___________________________________________________________
 AliSymMatrix&  AliSymMatrix::operator=(const AliSymMatrix& src)
 {
-  //
+  // assignment operator
   if (this != &src) {
     TObject::operator=(src);
     if (GetSizeBooked()!=src.GetSizeBooked() && GetSizeAdded()!=src.GetSizeAdded()) {
@@ -105,12 +108,12 @@ AliSymMatrix&  AliSymMatrix::operator=(const AliSymMatrix& src)
       int nmainel = src.GetSizeBooked()*(src.GetSizeBooked()+1);
       memcpy(fElems,src.fElems,nmainel*sizeof(Double_t));
       if (src.GetSizeAdded()) { // transfer extra rows to main matrix
-       Double_t *pnt = fElems + nmainel*sizeof(Double_t);
+       Double_t *pnt = fElems + nmainel;//*sizeof(Double_t);
        int ncl = src.GetSizeBooked() + 1;
        for (int ir=0;ir<src.GetSizeAdded();ir++) {
          ncl += ir; 
          memcpy(pnt,src.fElemsAdd[ir],ncl*sizeof(Double_t));
-         pnt += ncl*sizeof(Double_t);
+         pnt += ncl;//*sizeof(Double_t);
        }
       }
       //
@@ -131,7 +134,7 @@ AliSymMatrix&  AliSymMatrix::operator=(const AliSymMatrix& src)
 //___________________________________________________________
 AliSymMatrix& AliSymMatrix::operator+=(const AliSymMatrix& src)
 {
-  //
+  // add operator
   if (GetSizeUsed() != src.GetSizeUsed()) {
     AliError("Matrix sizes are different");
     return *this;
@@ -143,6 +146,7 @@ AliSymMatrix& AliSymMatrix::operator+=(const AliSymMatrix& src)
 //___________________________________________________________
 void AliSymMatrix::Clear(Option_t*)
 {
+  // clear dynamic part
   if (fElems) {delete[] fElems; fElems = 0;}
   //  
   if (fElemsAdd) {
@@ -159,13 +163,14 @@ Float_t AliSymMatrix::GetDensity() const
 {
   // get fraction of non-zero elements
   Int_t nel = 0;
-  for (int i=GetSizeUsed();i--;) for (int j=i+1;j--;) if (TMath::Abs(GetEl(i,j))>DBL_MIN) nel++;
+  for (int i=GetSizeUsed();i--;) for (int j=i+1;j--;) if (!IsZero(GetEl(i,j))) nel++;
   return 2.*nel/( (GetSizeUsed()+1)*GetSizeUsed() );
 }
 
 //___________________________________________________________
 void AliSymMatrix::Print(Option_t* option) const
 {
+  // print itself
   printf("Symmetric Matrix: Size = %d (%d rows added dynamically), %d used\n",GetSize(),GetSizeAdded(),GetSizeUsed());
   TString opt = option; opt.ToLower();
   if (opt.IsNull()) return;
@@ -178,7 +183,7 @@ void AliSymMatrix::Print(Option_t* option) const
 }
 
 //___________________________________________________________
-void AliSymMatrix::MultiplyByVec(Double_t *vecIn,Double_t *vecOut) const
+void AliSymMatrix::MultiplyByVec(const Double_t *vecIn,Double_t *vecOut) const
 {
   // fill vecOut by matrix*vecIn
   // vector should be of the same size as the matrix
@@ -206,7 +211,7 @@ AliSymMatrix* AliSymMatrix::DecomposeChol()
       fgBuffer = new AliSymMatrix(*this);
     }
     catch(bad_alloc&) {
-      printf("Failed to allocate memory for Choleski decompostions\n");
+      AliInfo("Failed to allocate memory for Choleski decompostions");
       return 0;
     }
   }
@@ -222,8 +227,9 @@ AliSymMatrix* AliSymMatrix::DecomposeChol()
       for (int k=i-1;k>=0;k--) if (rowi[k]&&rowj[k]) sum -= rowi[k]*rowj[k];
       if (i == j) {
        if (sum <= 0.0) { // not positive-definite
-         printf("The matrix is not positive definite [%e]\n"
-                "Choleski decomposition is not possible\n",sum);
+         AliInfo(Form("The matrix is not positive definite [%e]\n"
+                      "Choleski decomposition is not possible",sum));
+         Print("l");
          return 0;
        }
        rowi[i] = TMath::Sqrt(sum);
@@ -241,7 +247,7 @@ Bool_t AliSymMatrix::InvertChol()
   //
   AliSymMatrix* mchol = DecomposeChol();
   if (!mchol) {
-    printf("Failed to invert the matrix\n");
+    AliInfo("Failed to invert the matrix");
     return kFALSE;
   }
   //
@@ -304,7 +310,7 @@ Bool_t AliSymMatrix::SolveChol(Double_t *b, Bool_t invert)
   //
   AliSymMatrix *pmchol = DecomposeChol();
   if (!pmchol) {
-    printf("SolveChol failed\n");
+    AliInfo("SolveChol failed");
     //    Print("l");
     return kFALSE;
   }
@@ -352,6 +358,7 @@ Bool_t AliSymMatrix::SolveChol(TVectorD &brhs, TVectorD &bsol,Bool_t invert)
 //___________________________________________________________
 void AliSymMatrix::AddRows(int nrows)
 {
+  // add empty rows
   if (nrows<1) return;
   Double_t **pnew = new Double_t*[nrows+fNrows];
   for (int ir=0;ir<fNrows;ir++) pnew[ir] = fElemsAdd[ir]; // copy old extra rows
@@ -419,10 +426,11 @@ void AliSymMatrix::AddToRow(Int_t r, Double_t *valc,Int_t *indc,Int_t n)
 //___________________________________________________________
 Double_t* AliSymMatrix::GetRow(Int_t r)
 {
+  // get pointer on the row
   if (r>=GetSize()) {
     int nn = GetSize();
     AddRows(r-GetSize()+1); 
-    printf("create %d of %d\n",r, nn);
+    AliDebug(2,Form("create %d of %d\n",r, nn));
     return &((fElemsAdd[r-GetSizeBooked()])[0]);
   }
   else return &fElems[GetIndex(r,0)];
@@ -451,7 +459,7 @@ int AliSymMatrix::SolveSpmInv(double *vecB, Bool_t stabilize)
     for (Int_t i=nGlo; i--;) rowMax[i] = colMax[i] = 0.0;
     for (Int_t i=nGlo; i--;) for (Int_t j=i+1;j--;) { 
        double vl = TMath::Abs(Query(i,j));
-       if (vl<DBL_MIN) continue;
+       if (IsZero(vl)) continue;
        if (vl > rowMax[i]) rowMax[i] = vl; // Max elemt of row i
        if (vl > colMax[j]) colMax[j] = vl; // Max elemt of column j
        if (i==j) continue;
@@ -460,8 +468,8 @@ int AliSymMatrix::SolveSpmInv(double *vecB, Bool_t stabilize)
       }
     //
     for (Int_t i=nGlo; i--;) {
-      if (TMath::Abs(rowMax[i])>DBL_MIN) rowMax[i] = 1./rowMax[i]; // Max elemt of row i
-      if (TMath::Abs(colMax[i])>DBL_MIN) colMax[i] = 1./colMax[i]; // Max elemt of column i
+      if (!IsZero(rowMax[i])) rowMax[i] = 1./rowMax[i]; // Max elemt of row i
+      if (!IsZero(colMax[i])) colMax[i] = 1./colMax[i]; // Max elemt of column i
     }
     //
   }
@@ -474,7 +482,7 @@ int AliSymMatrix::SolveSpmInv(double *vecB, Bool_t stabilize)
       fgBuffer = new AliSymMatrix(*this);
     }
     catch(bad_alloc&) {
-      printf("Failed to allocate memory for matrix inversion buffer\n");
+      AliError("Failed to allocate memory for matrix inversion buffer");
       return 0;
     }
   }
@@ -483,11 +491,11 @@ int AliSymMatrix::SolveSpmInv(double *vecB, Bool_t stabilize)
   if (stabilize) for (int i=0;i<nGlo; i++) { // Small loop for matrix equilibration (gives a better conditioning) 
       for (int j=0;j<=i; j++) {
        double vl = Query(i,j);
-       if (TMath::Abs(vl)>DBL_MIN) SetEl(i,j, TMath::Sqrt(rowMax[i])*vl*TMath::Sqrt(colMax[j]) ); // Equilibrate the V matrix
+       if (!IsZero(vl)) SetEl(i,j, TMath::Sqrt(rowMax[i])*vl*TMath::Sqrt(colMax[j]) ); // Equilibrate the V matrix
       }
       for (int j=i+1;j<nGlo;j++) {
        double vl = Query(j,i);
-       if (TMath::Abs(vl)>DBL_MIN) fgBuffer->SetEl(j,i,TMath::Sqrt(rowMax[i])*vl*TMath::Sqrt(colMax[j]) ); // Equilibrate the V matrix
+       if (!IsZero(vl)) fgBuffer->SetEl(j,i,TMath::Sqrt(rowMax[i])*vl*TMath::Sqrt(colMax[j]) ); // Equilibrate the V matrix
       }
     }
   //
@@ -541,11 +549,11 @@ int AliSymMatrix::SolveSpmInv(double *vecB, Bool_t stabilize)
     }
   }
   //
-  for (Int_t i=0; i<nGlo; i++) for (Int_t j=0; j<nGlo; j++) {
-      double vl = TMath::Sqrt(colMax[i])*TMath::Sqrt(rowMax[j]); // Correct matrix V
-      if (i>=j) (*this)(i,j) *= vl;
-      else      (*fgBuffer)(j,i) *= vl;
-    }
+  if (stabilize) for (Int_t i=0; i<nGlo; i++) for (Int_t j=0; j<nGlo; j++) {
+       double vl = TMath::Sqrt(colMax[i])*TMath::Sqrt(rowMax[j]); // Correct matrix V
+       if (i>=j) (*this)(i,j) *= vl;
+       else      (*fgBuffer)(j,i) *= vl;
+      }
   //
   for (Int_t j=0; j<nGlo; j++) {
     rowMax[j] = 0.0;