ITS/UPGRADE/AliITSUClusterLines.cxx

   1 #include "AliITSUClusterLines.h"
   2 #include <TMath.h>
   3
   4 ClassImp(AliITSUClusterLines);
   5
   6 //////////////////////////////////////////////////////////////////////
   7 // This class is used by the AliITSUVertexer to compute and store   //
   8 // information about vertex candidates.                             //
   9 // The cluster is seeded starting from two lines, then it is        //
  10 // possible to attach other lines. Whenever a new line is attached  //
  11 // the weight and coefficient matrices are computed.                //
  12 // Origin puccio@to.infn.it  Feb. 20 2014                           //
  13 //////////////////////////////////////////////////////////////////////
  14
  15
  16 //____________________________________________________________________________________________________________
  17 AliITSUClusterLines::AliITSUClusterLines() : TObject(),
  18                                              fA(),
  19                                              fB(),
  20                                              fLabels(),
  21                                              fV() {
  22   // Default Constructor
  23
  24   for(Int_t i=0;i<9;++i) fW[i] = 0.;
  25 }
  26
  27 //____________________________________________________________________________________________________________
  28 AliITSUClusterLines::AliITSUClusterLines(UInt_t first, AliStrLine *line1, UInt_t second, AliStrLine *line2,Bool_t weight) : TObject(),
  29                                                                                                                             fA(),
  30                                                                                                                             fB(),
  31                                                                                                                             fLabels(),
  32                                                                                                                             fV() {
  33   // Standard constructor
  34   fLabels.push_back(first);
  35   fLabels.push_back(second);
  36   Double_t wmat1[9],wmat2[9],p1[3],p2[3],cd1[3],cd2[3],sq1[3]={1.,1.,1.},sq2[3]={1.,1.,1.};
  37   line1->GetWMatrix(wmat1);
  38   line2->GetWMatrix(wmat2);
  39   for(Int_t i=0;i<9;++i) fW[i]=wmat1[i]+wmat2[i];
  40   line1->GetP0(p1);
  41   line2->GetP0(p2);
  42   line1->GetCd(cd1);
  43   line2->GetCd(cd2);
  44   if(weight) {
  45     line1->GetSigma2P0(sq1);
  46     line2->GetSigma2P0(sq2);
  47   }
  48
  49   Double_t det1=cd1[2]*cd1[2]*sq1[0]*sq1[1]+cd1[1]*cd1[1]*sq1[0]*sq1[2]+cd1[0]*cd1[0]*sq1[1]*sq1[2];
  50   Double_t det2=cd2[2]*cd2[2]*sq2[0]*sq2[1]+cd2[1]*cd2[1]*sq2[0]*sq2[2]+cd2[0]*cd2[0]*sq2[1]*sq2[2];
  51
  52   fA[0]=(cd1[2]*cd1[2]*sq1[1]+cd1[1]*cd1[1]*sq1[2])/det1+(cd2[2]*cd2[2]*sq2[1]+cd2[1]*cd2[1]*sq2[2])/det2;
  53   fA[1]=-cd1[0]*cd1[1]*sq1[2]/det1-cd2[0]*cd2[1]*sq2[2]/det2;
  54   fA[2]=-cd1[0]*cd1[2]*sq1[1]/det1-cd2[0]*cd2[2]*sq2[1]/det2;
  55   fA[3]=(cd1[2]*cd1[2]*sq1[0]+cd1[0]*cd1[0]*sq1[2])/det1+(cd2[2]*cd2[2]*sq2[0]+cd2[0]*cd2[0]*sq2[2])/det2;
  56   fA[4]=-cd1[1]*cd1[2]*sq1[0]/det1-cd2[1]*cd2[2]*sq2[0]/det2;
  57   fA[5]=(cd1[1]*cd1[1]*sq1[0]+cd1[0]*cd1[0]*sq1[1])/det1+(cd2[1]*cd2[1]*sq2[0]+cd2[0]*cd2[0]*sq2[1])/det2;
  58
  59   fB[0]=(cd1[1]*sq1[2]*(-cd1[1]*p1[0]+cd1[0]*p1[1])+cd1[2]*sq1[1]*(-cd1[2]*p1[0]+cd1[0]*p1[2]))/det1;
  60   fB[0]+=(cd2[1]*sq2[2]*(-cd2[1]*p2[0]+cd2[0]*p2[1])+cd2[2]*sq2[1]*(-cd2[2]*p2[0]+cd2[0]*p2[2]))/det2;
  61   fB[1]=(cd1[0]*sq1[2]*(-cd1[0]*p1[1]+cd1[1]*p1[0])+cd1[2]*sq1[0]*(-cd1[2]*p1[1]+cd1[1]*p1[2]))/det1;
  62   fB[1]+=(cd2[0]*sq2[2]*(-cd2[0]*p2[1]+cd2[1]*p2[0])+cd2[2]*sq2[0]*(-cd2[2]*p2[1]+cd2[1]*p2[2]))/det2;
  63   fB[2]=(cd1[0]*sq1[1]*(-cd1[0]*p1[2]+cd1[2]*p1[0])+cd1[1]*sq1[0]*(-cd1[1]*p1[2]+cd1[2]*p1[1]))/det1;
  64   fB[2]+=(cd2[0]*sq2[1]*(-cd2[0]*p2[2]+cd2[2]*p2[0])+cd2[1]*sq2[0]*(-cd2[1]*p2[2]+cd2[2]*p2[1]))/det2;
  65
  66   this->ComputeClusterCentroid();
  67
  68 }
  69
  70 //____________________________________________________________________________________________________________
  71 AliITSUClusterLines::~AliITSUClusterLines() {
  72   // Destructor
  73 }
  74
  75 //____________________________________________________________________________________________________________
  76 void AliITSUClusterLines::Add(UInt_t label, AliStrLine *line, Bool_t weight) {
  77   // Add a line to the cluster. It changes the weight matrix of the cluster and its parameters
  78   fLabels.push_back(label);
  79   Double_t wmat[9],p[3],cd[3],sq[3]={1.,1.,1.};
  80   line->GetWMatrix(wmat);
  81   line->GetP0(p);
  82   line->GetCd(cd);
  83
  84   for(Int_t i=0;i<9;++i) fW[i]+=wmat[i];
  85   if(weight) line->GetSigma2P0(sq);
  86
  87   Double_t det=cd[2]*cd[2]*sq[0]*sq[1]+cd[1]*cd[1]*sq[0]*sq[2]+cd[0]*cd[0]*sq[1]*sq[2];
  88   fA[0]+=(cd[2]*cd[2]*sq[1]+cd[1]*cd[1]*sq[2])/det;
  89   fA[1]+=-cd[0]*cd[1]*sq[2]/det;
  90   fA[2]+=-cd[0]*cd[2]*sq[1]/det;
  91   fA[3]+=(cd[2]*cd[2]*sq[0]+cd[0]*cd[0]*sq[2])/det;
  92   fA[4]+=-cd[1]*cd[2]*sq[0]/det;
  93   fA[5]+=(cd[1]*cd[1]*sq[0]+cd[0]*cd[0]*sq[1])/det;
  94
  95   fB[0]+=(cd[1]*sq[2]*(-cd[1]*p[0]+cd[0]*p[1])+cd[2]*sq[1]*(-cd[2]*p[0]+cd[0]*p[2]))/det;
  96   fB[1]+=(cd[0]*sq[2]*(-cd[0]*p[1]+cd[1]*p[0])+cd[2]*sq[0]*(-cd[2]*p[1]+cd[1]*p[2]))/det;
  97   fB[2]+=(cd[0]*sq[1]*(-cd[0]*p[2]+cd[2]*p[0])+cd[1]*sq[0]*(-cd[1]*p[2]+cd[2]*p[1]))/det;
  98
  99   this->ComputeClusterCentroid();
 100 }
 101
 102 //____________________________________________________________________________________________________________
 103 Int_t AliITSUClusterLines::Compare(const TObject* obj) const {
 104   // Comparison criteria between two clusters
 105   const AliITSUClusterLines *cl=(const AliITSUClusterLines*)obj;
 106   if(fLabels.size()<cl->GetSize()) return 1;
 107   if(fLabels.size()>cl->GetSize()) return -1;
 108   return 0;
 109 }
 110
 111 //____________________________________________________________________________________________________________
 112 void AliITSUClusterLines::ComputeClusterCentroid() {
 113   // Calculation of the centroid
 114   Double_t *a=fA;
 115   Double_t *b=fB;
 116   Double_t *v=fV;
 117
 118   Double_t det=a[0]*(a[3]*a[5]-a[4]*a[4])-a[1]*(a[1]*a[5]-a[4]*a[2])+a[2]*(a[1]*a[4]-a[2]*a[3]);
 119
 120   if(det==0) {
 121     cout << "Could not invert weight matrix" << endl;
 122     return;
 123   }
 124   v[0]=-(b[0]*(a[3]*a[5]-a[4]*a[4])-a[1]*(b[1]*a[5]-a[4]*b[2])+a[2]*(b[1]*a[4]-b[2]*a[3]))/det;
 125   v[1]=-(a[0]*(b[1]*a[5]-b[2]*a[4])-b[0]*(a[1]*a[5]-a[4]*a[2])+a[2]*(a[1]*b[2]-a[2]*b[1]))/det;
 126   v[2]=-(a[0]*(a[3]*b[2]-b[1]*a[4])-a[1]*(a[1]*b[2]-b[1]*a[2])+b[0]*(a[1]*a[4]-a[2]*a[3]))/det;
 127 }
 128
 129 //____________________________________________________________________________________________________________
 130 void AliITSUClusterLines::GetCovMatrix(Float_t cov[6]) {
 131   // Returns the covariance matrix (single precision)
 132   Double_t *w=fW;
 133   Double_t den=w[0]*(w[4]*w[8]-w[5]*w[7])-w[1]*(w[3]*w[8]-w[5]*w[6])+w[2]*(w[3]*w[7]-w[4]*w[6]);
 134   if(den==0) {
 135     cout << "Could not invert weight matrix" << endl;
 136     return;
 137   }
 138   cov[0]=(w[4]*w[8]-w[5]*w[7])/den;
 139   cov[1]=-(w[1]*w[8]-w[7]*w[2])/den;
 140   cov[2]=(w[1]*w[5]-w[4]*w[2])/den;
 141   cov[3]=(w[0]*w[8]-w[6]*w[2])/den;
 142   cov[4]=-(w[0]*w[5]-w[3]*w[2])/den;
 143   cov[5]=(w[0]*w[4]-w[1]*w[3])/den;
 144 }
 145
 146 //____________________________________________________________________________________________________________
 147 void AliITSUClusterLines::GetCovMatrix(Double_t cov[6]) {
 148   // Returns the covariance matrix (double precision)
 149   Double_t *w=fW;
 150   Double_t den=w[0]*(w[4]*w[8]-w[5]*w[7])-w[1]*(w[3]*w[8]-w[5]*w[6])+w[2]*(w[3]*w[7]-w[4]*w[6]);
 151   if(den==0) {
 152     cout << "Could not invert weight matrix" << endl;
 153     return;
 154   }
 155   cov[0]=(w[4]*w[8]-w[5]*w[7])/den;
 156   cov[1]=-(w[1]*w[8]-w[7]*w[2])/den;
 157   cov[2]=(w[1]*w[5]-w[4]*w[2])/den;
 158   cov[3]=(w[0]*w[8]-w[6]*w[2])/den;
 159   cov[4]=-(w[0]*w[5]-w[3]*w[2])/den;
 160   cov[5]=(w[0]*w[4]-w[1]*w[3])/den;
 161
 162 }
 163
 164 //____________________________________________________________________________________________________________
 165 Bool_t AliITSUClusterLines::IsEqual(const TObject* obj) const {
 166   // Comparison criteria between two clusters
 167   const AliITSUClusterLines *cl=(const AliITSUClusterLines*)obj;
 168   if(fLabels.size()==cl->GetSize()) return kTRUE;
 169   return kFALSE;
 170 }
 171