]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliTrackPointArray.cxx
Partial fix for bug #59809: putting HLT trigger decision in the ESD, not yet in stand...
[u/mrichter/AliRoot.git] / STEER / AliTrackPointArray.cxx
index 92b2b15db0ff3628f9bc1a09c3807ec98ef92356..5572817e12d02680e32ebdc183bf97f2f634b6cf 100644 (file)
 //////////////////////////////////////////////////////////////////////////////
 
 #include <TMath.h>
+#include <TMatrixD.h>
 #include <TMatrixDSym.h>
+#include <TGeoMatrix.h>
+#include <TMatrixDSymEigen.h>
 
 #include "AliTrackPointArray.h"
 
@@ -32,10 +35,16 @@ ClassImp(AliTrackPointArray)
 //______________________________________________________________________________
 AliTrackPointArray::AliTrackPointArray() :
   TObject(),
+  fSorted(kFALSE),
   fNPoints(0),
   fX(0),
   fY(0),
   fZ(0),
+  fCharge(0),
+  fDriftTime(0),
+  fChargeRatio(0),
+  fClusterType(0),
+  fIsExtra(0),
   fSize(0),
   fCov(0),
   fVolumeID(0)
@@ -45,25 +54,50 @@ AliTrackPointArray::AliTrackPointArray() :
 //______________________________________________________________________________
 AliTrackPointArray::AliTrackPointArray(Int_t npoints):
   TObject(),
+  fSorted(kFALSE),
   fNPoints(npoints),
   fX(new Float_t[npoints]),
   fY(new Float_t[npoints]),
   fZ(new Float_t[npoints]),
+  fCharge(new Float_t[npoints]),
+  fDriftTime(new Float_t[npoints]),
+  fChargeRatio(new Float_t[npoints]),
+  fClusterType(new Int_t[npoints]),
+  fIsExtra(new Bool_t[npoints]),
   fSize(6*npoints),
   fCov(new Float_t[fSize]),
   fVolumeID(new UShort_t[npoints])
 {
   // Constructor
   //
+  for (Int_t ip=0; ip<npoints;ip++){
+    fX[ip]=0;
+    fY[ip]=0;
+    fZ[ip]=0;
+    fCharge[ip]=0;
+    fDriftTime[ip]=0;
+    fChargeRatio[ip]=0;
+    fClusterType[ip]=0;
+    fIsExtra[ip]=kFALSE;
+    fVolumeID[ip]=0;
+    for (Int_t icov=0;icov<6; icov++)
+      fCov[6*ip+icov]=0;
+  }
 }
 
 //______________________________________________________________________________
 AliTrackPointArray::AliTrackPointArray(const AliTrackPointArray &array):
   TObject(array),
+  fSorted(array.fSorted),
   fNPoints(array.fNPoints),
   fX(new Float_t[fNPoints]),
   fY(new Float_t[fNPoints]),
   fZ(new Float_t[fNPoints]),
+  fCharge(new Float_t[fNPoints]),
+  fDriftTime(new Float_t[fNPoints]),
+  fChargeRatio(new Float_t[fNPoints]),
+  fClusterType(new Int_t[fNPoints]),
+  fIsExtra(new Bool_t[fNPoints]),
   fSize(array.fSize),
   fCov(new Float_t[fSize]),
   fVolumeID(new UShort_t[fNPoints])
@@ -73,6 +107,31 @@ AliTrackPointArray::AliTrackPointArray(const AliTrackPointArray &array):
   memcpy(fX,array.fX,fNPoints*sizeof(Float_t));
   memcpy(fY,array.fY,fNPoints*sizeof(Float_t));
   memcpy(fZ,array.fZ,fNPoints*sizeof(Float_t));
+  if (array.fCharge) {
+    memcpy(fCharge,array.fCharge,fNPoints*sizeof(Float_t));
+  } else {
+    memset(fCharge, 0, fNPoints*sizeof(Float_t));
+  }
+  if (array.fDriftTime) {
+    memcpy(fDriftTime,array.fDriftTime,fNPoints*sizeof(Float_t));
+  } else {
+    memset(fDriftTime, 0, fNPoints*sizeof(Float_t));
+  }
+  if (array.fChargeRatio) {
+    memcpy(fChargeRatio,array.fChargeRatio,fNPoints*sizeof(Float_t));
+  } else {
+    memset(fChargeRatio, 0, fNPoints*sizeof(Float_t));
+  }
+  if (array.fClusterType) {
+    memcpy(fClusterType,array.fClusterType,fNPoints*sizeof(Int_t));
+  } else {
+    memset(fClusterType, 0, fNPoints*sizeof(Int_t));
+  }
+  if (array.fIsExtra) {
+    memcpy(fIsExtra,array.fIsExtra,fNPoints*sizeof(Bool_t));
+  } else {
+    memset(fIsExtra, 0, fNPoints*sizeof(Bool_t));
+  }
   memcpy(fVolumeID,array.fVolumeID,fNPoints*sizeof(UShort_t));
   memcpy(fCov,array.fCov,fSize*sizeof(Float_t));
 }
@@ -85,16 +144,37 @@ AliTrackPointArray &AliTrackPointArray::operator =(const AliTrackPointArray& arr
   if(this==&array) return *this;
   ((TObject *)this)->operator=(array);
 
+  fSorted = array.fSorted;
   fNPoints = array.fNPoints;
   fSize = array.fSize;
+  delete [] fX;
   fX = new Float_t[fNPoints];
+  delete [] fY;
   fY = new Float_t[fNPoints];
+  delete [] fZ;
   fZ = new Float_t[fNPoints];
+  delete [] fCharge;
+  fCharge = new Float_t[fNPoints];
+  delete [] fDriftTime;
+  fDriftTime = new Float_t[fNPoints];
+  delete [] fChargeRatio;
+  fChargeRatio = new Float_t[fNPoints];
+  delete [] fClusterType;
+  fClusterType = new Int_t[fNPoints];
+  delete [] fIsExtra;
+  fIsExtra = new Bool_t[fNPoints];
+  delete [] fVolumeID;
   fVolumeID = new UShort_t[fNPoints];
+  delete [] fCov;
   fCov = new Float_t[fSize];
   memcpy(fX,array.fX,fNPoints*sizeof(Float_t));
   memcpy(fY,array.fY,fNPoints*sizeof(Float_t));
   memcpy(fZ,array.fZ,fNPoints*sizeof(Float_t));
+  memcpy(fCharge,array.fCharge,fNPoints*sizeof(Float_t));
+  memcpy(fDriftTime,array.fDriftTime,fNPoints*sizeof(Float_t));
+  memcpy(fChargeRatio,array.fChargeRatio,fNPoints*sizeof(Float_t));
+  memcpy(fClusterType,array.fClusterType,fNPoints*sizeof(Int_t));
+  memcpy(fIsExtra,array.fIsExtra,fNPoints*sizeof(Bool_t));
   memcpy(fVolumeID,array.fVolumeID,fNPoints*sizeof(UShort_t));
   memcpy(fCov,array.fCov,fSize*sizeof(Float_t));
 
@@ -109,6 +189,11 @@ AliTrackPointArray::~AliTrackPointArray()
   delete [] fX;
   delete [] fY;
   delete [] fZ;
+  delete [] fCharge;
+  delete [] fDriftTime;
+  delete [] fChargeRatio;
+  delete [] fClusterType;
+  delete [] fIsExtra;
   delete [] fVolumeID;
   delete [] fCov;
 }
@@ -123,11 +208,17 @@ Bool_t AliTrackPointArray::AddPoint(Int_t i, const AliTrackPoint *p)
   fX[i] = p->GetX();
   fY[i] = p->GetY();
   fZ[i] = p->GetZ();
+  fCharge[i] = p->GetCharge();
+  fDriftTime[i] = p->GetDriftTime();
+  fChargeRatio[i]=p->GetChargeRatio();
+  fClusterType[i]=p->GetClusterType();
+  fIsExtra[i] = p->IsExtra();
   fVolumeID[i] = p->GetVolumeID();
   memcpy(&fCov[6*i],p->GetCov(),6*sizeof(Float_t));
   return kTRUE;
 }
 
+
 //______________________________________________________________________________
 Bool_t AliTrackPointArray::GetPoint(AliTrackPoint &p, Int_t i) const
 {
@@ -136,6 +227,11 @@ Bool_t AliTrackPointArray::GetPoint(AliTrackPoint &p, Int_t i) const
   if (i >= fNPoints) return kFALSE;
   p.SetXYZ(fX[i],fY[i],fZ[i],&fCov[6*i]);
   p.SetVolumeID(fVolumeID[i]);
+  p.SetCharge(fCharge ? fCharge[i] : 0);
+  p.SetDriftTime(fDriftTime ? fDriftTime[i] : 0);
+  p.SetChargeRatio(fChargeRatio ? fChargeRatio[i] : 0);
+  p.SetClusterType(fClusterType ? fClusterType[i] : 0);
+  p.SetExtra(fIsExtra ? fIsExtra[i] : kFALSE);
   return kTRUE;
 }
 
@@ -152,6 +248,29 @@ Bool_t AliTrackPointArray::HasVolumeID(UShort_t volid) const
   return check;
 }
 
+//______________________________________________________________________________
+void AliTrackPointArray::Sort(Bool_t down)
+{
+  // Sort the array by the values of Y-coordinate of the track points.
+  // The order is given by "down".
+  // Optimized more for maintenance rather than for speed.
+  if (fSorted) return;
+
+  Int_t *index=new Int_t[fNPoints];
+  AliTrackPointArray a(*this);
+  TMath::Sort(fNPoints,a.GetY(),index,down);
+  AliTrackPoint p;
+  for (Int_t i = 0; i < fNPoints; i++) {
+    a.GetPoint(p,index[i]);
+    AddPoint(i,&p);
+  }
+
+  delete[] index;
+  fSorted=kTRUE;
+}
+
 ClassImp(AliTrackPoint)
 
 //______________________________________________________________________________
@@ -160,6 +279,11 @@ AliTrackPoint::AliTrackPoint() :
   fX(0),
   fY(0),
   fZ(0),
+  fCharge(0),
+  fDriftTime(0),
+  fChargeRatio(0),
+  fClusterType(0),
+  fIsExtra(kFALSE),
   fVolumeID(0)
 {
   // Default constructor
@@ -169,31 +293,49 @@ AliTrackPoint::AliTrackPoint() :
 
 
 //______________________________________________________________________________
-AliTrackPoint::AliTrackPoint(Float_t x, Float_t y, Float_t z, const Float_t *cov, UShort_t volid) :
+AliTrackPoint::AliTrackPoint(Float_t x, Float_t y, Float_t z, const Float_t *cov, UShort_t volid, Float_t charge, Float_t drifttime,Float_t chargeratio, Int_t clutyp) :
   TObject(),
   fX(0),
   fY(0),
   fZ(0),
+  fCharge(0),
+  fDriftTime(0),
+  fChargeRatio(0),
+  fClusterType(0),
+  fIsExtra(kFALSE),
   fVolumeID(0)
 {
   // Constructor
   //
   SetXYZ(x,y,z,cov);
+  SetCharge(charge);
+  SetDriftTime(drifttime);
+  SetChargeRatio(chargeratio);
+  SetClusterType(clutyp);
   SetVolumeID(volid);
 }
 
 //______________________________________________________________________________
-AliTrackPoint::AliTrackPoint(const Float_t *xyz, const Float_t *cov, UShort_t volid) :
+AliTrackPoint::AliTrackPoint(const Float_t *xyz, const Float_t *cov, UShort_t volid, Float_t charge, Float_t drifttime,Float_t chargeratio, Int_t clutyp) :
   TObject(),
   fX(0),
   fY(0),
   fZ(0),
+  fCharge(0),
+  fDriftTime(0),
+  fChargeRatio(0),
+  fClusterType(0),
+  fIsExtra(kFALSE),
   fVolumeID(0)
 {
   // Constructor
   //
   SetXYZ(xyz[0],xyz[1],xyz[2],cov);
+  SetCharge(charge);
+  SetDriftTime(drifttime);
+  SetChargeRatio(chargeratio);
   SetVolumeID(volid);
+  SetClusterType(clutyp);
 }
 
 //______________________________________________________________________________
@@ -202,11 +344,21 @@ AliTrackPoint::AliTrackPoint(const AliTrackPoint &p):
   fX(0),
   fY(0),
   fZ(0),
+  fCharge(0),
+  fDriftTime(0),
+  fChargeRatio(0),
+  fClusterType(0),
+  fIsExtra(kFALSE),
   fVolumeID(0)
 {
   // Copy constructor
   //
   SetXYZ(p.fX,p.fY,p.fZ,&(p.fCov[0]));
+  SetCharge(p.fCharge);
+  SetDriftTime(p.fDriftTime);
+  SetChargeRatio(p.fChargeRatio);
+  SetClusterType(p.fClusterType);
+  SetExtra(p.fIsExtra);
   SetVolumeID(p.fVolumeID);
 }
 
@@ -219,6 +371,11 @@ AliTrackPoint &AliTrackPoint::operator =(const AliTrackPoint& p)
   ((TObject *)this)->operator=(p);
 
   SetXYZ(p.fX,p.fY,p.fZ,&(p.fCov[0]));
+  SetCharge(p.fCharge);
+  SetDriftTime(p.fDriftTime);
+  SetChargeRatio(p.fChargeRatio);
+  SetClusterType(p.fClusterType);
+  SetExtra(p.fIsExtra);
   SetVolumeID(p.fVolumeID);
 
   return *this;
@@ -299,6 +456,64 @@ Float_t AliTrackPoint::GetResidual(const AliTrackPoint &p, Bool_t weighted) cons
   return res;
 }
 
+//_____________________________________________________________________________
+Bool_t AliTrackPoint::GetPCA(const AliTrackPoint &p, AliTrackPoint &out) const
+{
+  //
+  // Get the intersection point between this point and
+  // the point "p" belongs to.
+  // The result is stored as a point 'out'
+  // return kFALSE in case of failure.
+  out.SetXYZ(0,0,0);
+
+  TMatrixD t(3,1);
+  t(0,0)=GetX();
+  t(1,0)=GetY();
+  t(2,0)=GetZ();
+  TMatrixDSym tC(3);
+  {
+  const Float_t *cv=GetCov();
+  tC(0,0)=cv[0]; tC(0,1)=cv[1]; tC(0,2)=cv[2];
+  tC(1,0)=cv[1]; tC(1,1)=cv[3]; tC(1,2)=cv[4];
+  tC(2,0)=cv[2]; tC(2,1)=cv[4]; tC(2,2)=cv[5];
+  }
+
+  TMatrixD m(3,1);
+  m(0,0)=p.GetX();
+  m(1,0)=p.GetY();
+  m(2,0)=p.GetZ();
+  TMatrixDSym mC(3);
+  {
+  const Float_t *cv=p.GetCov();
+  mC(0,0)=cv[0]; mC(0,1)=cv[1]; mC(0,2)=cv[2];
+  mC(1,0)=cv[1]; mC(1,1)=cv[3]; mC(1,2)=cv[4];
+  mC(2,0)=cv[2]; mC(2,1)=cv[4]; mC(2,2)=cv[5];
+  }
+
+  TMatrixDSym tmW(tC);
+  tmW+=mC;
+  tmW.Invert();
+  if (!tmW.IsValid()) return kFALSE; 
+
+  TMatrixD mW(tC,TMatrixD::kMult,tmW);
+  TMatrixD tW(mC,TMatrixD::kMult,tmW);
+
+  TMatrixD mi(mW,TMatrixD::kMult,m);
+  TMatrixD ti(tW,TMatrixD::kMult,t);
+  ti+=mi;
+
+  TMatrixD iC(tC,TMatrixD::kMult,tmW);
+  iC*=mC;
+
+  out.SetXYZ(ti(0,0),ti(1,0),ti(2,0));
+  UShort_t id=p.GetVolumeID();
+  out.SetVolumeID(id);
+
+  return kTRUE;
+}
+
 //______________________________________________________________________________
 Float_t AliTrackPoint::GetAngle() const
 {
@@ -320,6 +535,32 @@ Float_t AliTrackPoint::GetAngle() const
 
 }
 
+//______________________________________________________________________________
+Bool_t AliTrackPoint::GetRotMatrix(TGeoRotation& rot) const
+{
+  // Returns the orientation of the
+  // sensitive layer (using cluster
+  // covariance matrix).
+  // Assumes that cluster has errors only in the layer's plane.
+  // Return value is kTRUE in case of success.
+
+  TMatrixDSym mcov(3);
+  {
+    const Float_t *cov=GetCov();
+    mcov(0,0)=cov[0]; mcov(0,1)=cov[1]; mcov(0,2)=cov[2];
+    mcov(1,0)=cov[1]; mcov(1,1)=cov[3]; mcov(1,2)=cov[4];
+    mcov(2,0)=cov[2]; mcov(2,1)=cov[4]; mcov(2,2)=cov[5];
+  }
+
+  TMatrixDSymEigen eigen(mcov);
+  TMatrixD eigenMatrix = eigen.GetEigenVectors();
+
+  rot.SetMatrix(eigenMatrix.GetMatrixArray());
+
+  return kTRUE;
+}
+
+
 //_____________________________________________________________________________
 AliTrackPoint& AliTrackPoint::Rotate(Float_t alpha) const
 {
@@ -383,5 +624,54 @@ void AliTrackPoint::Print(Option_t *) const
   printf("X = %12.6f    Tx = %12.6f%12.6f%12.6f\n", fX, fCov[0], fCov[1], fCov[2]);
   printf("Y = %12.6f    Ty = %12.6f%12.6f%12.6f\n", fY, fCov[1], fCov[3], fCov[4]);
   printf("Z = %12.6f    Tz = %12.6f%12.6f%12.6f\n", fZ, fCov[2], fCov[4], fCov[5]);
+  printf("Charge = %f\n", fCharge);
+  printf("Drift Time = %f\n", fDriftTime);
+  printf("Charge Ratio  = %f\n", fChargeRatio);
+  printf("Cluster Type  = %d\n", fClusterType);
+  if(fIsExtra) printf("This is an extra point\n");
 
 }
+
+
+//________________________________
+void AliTrackPoint::SetAlignCovMatrix(const TMatrixDSym alignparmtrx){
+  // Add the uncertainty on the cluster position due to alignment
+  // (using the 6x6 AliAlignObj Cov. Matrix alignparmtrx) to the already
+  // present Cov. Matrix 
+
+  TMatrixDSym cov(3);
+  TMatrixD coval(3,3);
+  TMatrixD jacob(3,6);
+  Float_t newcov[6];
+
+  cov(0,0)=fCov[0];
+  cov(1,0)=cov(0,1)=fCov[1];
+  cov(2,0)=cov(0,2)=fCov[2];
+  cov(1,1)=fCov[3];
+  cov(2,1)=cov(1,2)=fCov[4];
+  cov(2,2)=fCov[5];
+  
+  jacob(0,0) = 1;      jacob(1,0) = 0;       jacob(2,0) = 0;
+  jacob(0,1) = 0;      jacob(1,1) = 1;       jacob(2,1) = 0;
+  jacob(0,2) = 0;      jacob(1,2) = 0;       jacob(2,2) = 1;
+  jacob(0,3) = 0;      jacob(1,3) =-fZ;      jacob(2,3) = fY;
+  jacob(0,4) = fZ;     jacob(1,4) = 0;       jacob(2,4) =-fX;
+  jacob(0,5) = -fY;    jacob(1,5) = fX;      jacob(2,5) = 0;
+  
+  TMatrixD jacobT=jacob.T();jacob.T();
+  
+  coval=jacob*alignparmtrx*jacobT+cov;
+
+
+  newcov[0]=coval(0,0);
+  newcov[1]=coval(1,0);
+  newcov[2]=coval(2,0);
+  newcov[3]=coval(1,1);
+  newcov[4]=coval(2,1);
+  newcov[5]=coval(2,2);
+  SetXYZ(fX,fY,fZ,newcov);
+
+}
+
+