ITS/AliITSVertexer.cxx

   1 #include "AliLog.h"
   2 #include "AliMultiplicity.h"
   3 #include "AliITSgeomTGeo.h"
   4 #include "AliITSDetTypeRec.h"
   5 #include "AliITSVertexer.h"
   6 #include "AliITSLoader.h"
   7 #include "AliITSMultReconstructor.h"
   8 #include "AliITSRecPointContainer.h"
   9 #include "AliRunLoader.h"
  10
  11 const Float_t AliITSVertexer::fgkPipeRadius = 3.0;
  12
  13 ClassImp(AliITSVertexer)
  14
  15 //////////////////////////////////////////////////////////////////////
  16 // Base class for primary vertex reconstruction                     //
  17 // AliESDVertexer is a class for full 3D primary vertex finding     //
  18 // derived classes: AliITSvertexer3D, AliITSVertexerZ.              //
  19 //                  AliITSVertexerCosmics                           //
  20 //////////////////////////////////////////////////////////////////////
  21
  22 /* $Id$ */
  23
  24 //______________________________________________________________________
  25 AliITSVertexer::AliITSVertexer():AliVertexer(),
  26 fLadders(NULL),
  27 fLadOnLay2(0),
  28 fComputeMultiplicity(kFALSE),
  29 fDetTypeRec(NULL),
  30 fMinTrackletsForPilup(0),
  31 fIsPileup(0),
  32 fNTrpuv(-2),
  33 fZpuv(-9999999.),
  34 fNoVertices(0),
  35 fVertArray(NULL),
  36 fFirstEvent(0),
  37 fLastEvent(-1)
  38 {
  39   // Default Constructor
  40   SetLaddersOnLayer2();
  41   SetMinTrackletsForPilup();
  42   for(Int_t i=0; i<kNSPDMod;i++) fUseModule[i]=kTRUE;
  43 }
  44
  45 //______________________________________________________________________
  46 AliITSVertexer::~AliITSVertexer() {
  47   // Destructor
  48   delete [] fLadders;
  49   delete []fVertArray;
  50 }
  51
  52 //______________________________________________________________________
  53 void AliITSVertexer::ResetVertex(){
  54   // Resets vertex related data members
  55   if(fNoVertices > 0){
  56     if(fVertArray) delete []fVertArray;
  57     fVertArray = NULL;
  58     fNoVertices = 0;
  59   }
  60   fIsPileup=kFALSE;
  61   fNTrpuv=-2;
  62   fZpuv=-99999.;
  63
  64 }
  65 //______________________________________________________________________
  66 void AliITSVertexer::FindMultiplicity(TTree *itsClusterTree){
  67   // Invokes AliITSMultReconstructor to determine the
  68   // charged multiplicity in the pixel layers
  69   if(fMult){delete fMult; fMult = 0;}
  70
  71   Bool_t success=kTRUE;
  72   Bool_t cosmics=kFALSE;
  73   if(!fCurrentVertex)success=kFALSE;
  74   if(fCurrentVertex && fCurrentVertex->GetNContributors()<1)success=kFALSE;
  75   if(fCurrentVertex && strstr(fCurrentVertex->GetTitle(),"cosmics")) {
  76     success=kFALSE;
  77     cosmics=kTRUE;
  78   }
  79
  80   // get the FastOr bit mask
  81   TBits fastOrFiredMap = fDetTypeRec->GetFastOrFiredMap();
  82   TBits firedChipMap = fDetTypeRec->GetFiredChipMap(itsClusterTree);
  83
  84   AliITSMultReconstructor multReco;
  85
  86   if(!success){
  87     if(!cosmics) {
  88       AliDebug(1,"Tracklets multiplicity not determined because the primary vertex was not found");
  89       AliDebug(1,"Just counting the number of cluster-fired chips on the SPD layers");
  90     }
  91     if (!itsClusterTree) {
  92       AliError(" Invalid ITS cluster tree !\n");
  93       return;
  94     }
  95     multReco.LoadClusterFiredChips(itsClusterTree);
  96     Short_t nfcL1 = multReco.GetNFiredChips(0);
  97     Short_t nfcL2 = multReco.GetNFiredChips(1);
  98     fMult = new AliMultiplicity(0,0,0,0,0,0,0,0,0,0,0,nfcL1,nfcL2,fastOrFiredMap);
  99     fMult->SetFiredChipMap(firedChipMap);
 100     AliITSRecPointContainer* rcont = AliITSRecPointContainer::Instance();
 101     fMult->SetITSClusters(0,rcont->GetNClustersInLayer(1,itsClusterTree));
 102     for(Int_t kk=2;kk<=6;kk++){
 103       fMult->SetITSClusters(kk-1,rcont->GetNClustersInLayerFast(kk));
 104     }
 105     return;
 106   }
 107
 108   if (!itsClusterTree) {
 109     AliError(" Invalid ITS cluster tree !\n");
 110     return;
 111   }
 112   Double_t vtx[3];
 113   fCurrentVertex->GetXYZ(vtx);
 114   Float_t vtxf[3];
 115   for(Int_t i=0;i<3;i++)vtxf[i]=vtx[i];
 116   multReco.SetHistOn(kFALSE);
 117   multReco.Reconstruct(itsClusterTree,vtxf,vtxf);
 118   Int_t notracks=multReco.GetNTracklets();
 119   Float_t *tht = new Float_t [notracks];
 120   Float_t *phi = new Float_t [notracks];
 121   Float_t *dtht = new Float_t [notracks];
 122   Float_t *dphi = new Float_t [notracks];
 123   Int_t *labels = new Int_t[notracks];
 124   Int_t *labelsL2 = new Int_t[notracks];
 125   for(Int_t i=0;i<multReco.GetNTracklets();i++){
 126     tht[i] = multReco.GetTracklet(i)[0];
 127     phi[i] =  multReco.GetTracklet(i)[1];
 128     dtht[i] = multReco.GetTracklet(i)[3];
 129     dphi[i] = multReco.GetTracklet(i)[2];
 130     labels[i] = static_cast<Int_t>(multReco.GetTracklet(i)[4]);
 131     labelsL2[i] = static_cast<Int_t>(multReco.GetTracklet(i)[5]);
 132   }
 133   Int_t nosingleclus=multReco.GetNSingleClusters();
 134   Float_t *ths = new Float_t [nosingleclus];
 135   Float_t *phs = new Float_t [nosingleclus];
 136   Int_t *labelss = new Int_t [nosingleclus];
 137   for(Int_t i=0;i<nosingleclus;i++){
 138     ths[i] = multReco.GetCluster(i)[0];
 139     phs[i] = multReco.GetCluster(i)[1];
 140     labelss[i] = (Int_t)multReco.GetCluster(i)[2];
 141   }
 142   Short_t nfcL1 = multReco.GetNFiredChips(0);
 143   Short_t nfcL2 = multReco.GetNFiredChips(1);
 144   fMult = new AliMultiplicity(notracks,tht,phi,dtht,dphi,labels,labelsL2,nosingleclus,ths,phs,labelss,nfcL1,nfcL2,fastOrFiredMap);
 145   fMult->SetFiredChipMap(firedChipMap);
 146   AliITSRecPointContainer* rcont = AliITSRecPointContainer::Instance();
 147   fMult->SetITSClusters(0,rcont->GetNClustersInLayer(1,itsClusterTree));
 148   for(Int_t kk=2;kk<=6;kk++){
 149     fMult->SetITSClusters(kk-1,rcont->GetNClustersInLayerFast(kk));
 150   }
 151   delete [] tht;
 152   delete [] phi;
 153   delete [] dtht;
 154   delete [] dphi;
 155   delete [] ths;
 156   delete [] phs;
 157   delete [] labels;
 158   delete [] labelsL2;
 159   delete [] labelss;
 160
 161   return;
 162 }
 163
 164 //______________________________________________________________________
 165 void AliITSVertexer::SetLaddersOnLayer2(Int_t ladwid){
 166   // Calculates the array of ladders on layer 2 to be used with a
 167   // given ladder on layer 1
 168   if(ladwid == fLadOnLay2 && fLadders)return;
 169   fLadOnLay2=ladwid;
 170   Int_t ladtot1=AliITSgeomTGeo::GetNLadders(1);
 171   if(fLadders) delete [] fLadders;
 172   fLadders=new UShort_t[ladtot1];
 173
 174
 175   Double_t pos1[3],pos2[3];
 176   Int_t mod1=AliITSgeomTGeo::GetModuleIndex(2,1,1);
 177   AliITSgeomTGeo::GetTranslation(mod1,pos1);  // position of the module in the MRS
 178   Double_t phi0=TMath::ATan2(pos1[1],pos1[0]);
 179   if(phi0<0) phi0+=2*TMath::Pi();
 180   Int_t mod2=AliITSgeomTGeo::GetModuleIndex(2,2,1);
 181   AliITSgeomTGeo::GetTranslation(mod2,pos2);
 182   Double_t phi2=TMath::ATan2(pos2[1],pos2[0]);
 183   if(phi2<0) phi2+=2*TMath::Pi();
 184   Double_t deltaPhi= phi0-phi2; // phi width of a layer2 module
 185
 186   for(Int_t i= 0; i<ladtot1;i++){
 187     Int_t modlad= AliITSgeomTGeo::GetModuleIndex(1,i+1,1);
 188     Double_t posmod[3];
 189     AliITSgeomTGeo::GetTranslation(modlad,posmod);
 190     Double_t phimod=TMath::ATan2(posmod[1],posmod[0]);
 191     if(phimod<0) phimod+=2*TMath::Pi();
 192     Double_t phi1= phimod+deltaPhi*double(fLadOnLay2);
 193     if(phi1<0) phi1+=2*TMath::Pi();
 194     if(phi1>2*TMath::Pi()) phi1-=2*TMath::Pi();
 195     Double_t philad1=phi0-phi1;
 196     UShort_t lad1;
 197     Double_t ladder1=(philad1)/(deltaPhi) +1.;
 198     if(ladder1<1){ladder1=40+ladder1;}
 199     lad1=int(ladder1+0.5);
 200     fLadders[i]=lad1;
 201   }
 202 }
 203
 204
 205 //______________________________________________________________________
 206 void AliITSVertexer::Init(TString filename){
 207   // Initialize the vertexer in case of
 208   // analysis of an entire file
 209   AliRunLoader *rl = AliRunLoader::Instance();
 210   if(!rl){
 211     AliFatal("Run Loader not found");
 212     return;
 213   }
 214   if (fLastEvent < 0) SetLastEvent(rl->GetNumberOfEvents()-1);
 215
 216   AliITSLoader* itsloader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
 217   if(!filename.Contains("default"))itsloader->SetVerticesFileName(filename);
 218   if(!filename.Contains("null"))itsloader->LoadVertices("recreate");
 219 }
 220
 221 //______________________________________________________________________
 222 void AliITSVertexer::WriteCurrentVertex(){
 223   // Write the current AliVertex object to file fOutFile
 224   AliRunLoader *rl = AliRunLoader::Instance();
 225   AliITSLoader* itsLoader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
 226   fCurrentVertex->SetName("Vertex");
 227   //  const char * name = fCurrentVertex->GetName();
 228   //  itsLoader->SetVerticesContName(name);
 229   itsLoader->PostVertex(fCurrentVertex);
 230   itsLoader->WriteVertices();
 231 }
 232
 233 //______________________________________________________________________
 234 void AliITSVertexer::FindVertices(){
 235   // computes the vertices of the events in the range FirstEvent - LastEvent
 236
 237   AliRunLoader *rl = AliRunLoader::Instance();
 238   AliITSLoader* itsloader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
 239   itsloader->LoadRecPoints("read");
 240   for(Int_t i=fFirstEvent;i<=fLastEvent;i++){
 241     rl->GetEvent(i);
 242     TTree* cltree = itsloader->TreeR();
 243     FindVertexForCurrentEvent(cltree);
 244     if(fCurrentVertex){
 245       WriteCurrentVertex();
 246     }
 247     else {
 248       AliDebug(1,Form("Vertex not found for event %d",i));
 249     }
 250   }
 251 }