Reconstruction of multiplicity. Ipdated AlITSvertexerZ (Massimo)
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 28 Jun 2006 10:20:51 +0000 (10:20 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 28 Jun 2006 10:20:51 +0000 (10:20 +0000)
13 files changed:
ITS/AliITSVertexer.cxx
ITS/AliITSVertexer.h
ITS/AliITSVertexerZ.cxx
ITS/AliITSVertexerZ.h
STEER/AliESD.cxx
STEER/AliESD.h
STEER/AliMultiplicity.cxx [new file with mode: 0644]
STEER/AliMultiplicity.h [new file with mode: 0644]
STEER/AliReconstruction.cxx
STEER/AliVertexer.cxx
STEER/AliVertexer.h
STEER/ESDLinkDef.h
STEER/libESD.pkg

index 7e59cc7..d5e62ec 100644 (file)
@@ -2,6 +2,8 @@
 #include <AliITSVertexer.h>
 #include <AliRunLoader.h>
 #include <AliITSLoader.h>
+#include <AliMultiplicity.h>
+#include <AliITSMultReconstructor.h>
 
 ClassImp(AliITSVertexer)
 
@@ -57,6 +59,57 @@ AliITSVertexer& AliITSVertexer::operator=(const AliITSVertexer& /* vtxr */){
   return *this;
 }
 
+//______________________________________________________________________
+void AliITSVertexer::FindMultiplicity(Int_t evnumber){
+  // Invokes AliITSMultReconstructor to determine the
+  // charged multiplicity in the pixel layers
+  if(fMult){delete fMult; fMult = 0;}
+  Bool_t success=kTRUE;
+  if(!fCurrentVertex)success=kFALSE;
+  if(fCurrentVertex && fCurrentVertex->GetNContributors()<1)success=kFALSE;
+  if(!success){
+    AliWarning("Tracklets multiplicity not determined because the primary vertex was not found");
+    return;
+  }
+  AliITSMultReconstructor* multReco = new AliITSMultReconstructor();
+  AliRunLoader *rl =AliRunLoader::GetRunLoader();
+  AliITSLoader* itsLoader = (AliITSLoader*)rl->GetLoader("ITSLoader");
+  multReco->SetGeometry(itsLoader->GetITSgeom());
+  itsLoader->LoadRecPoints();
+  rl->GetEvent(evnumber);
+  TTree* itsClusterTree = itsLoader->TreeR();
+  if (!itsClusterTree) {
+    AliError(" Can't get the ITS cluster tree !\n");
+    return;
+  }
+  Double_t vtx[3];
+  fCurrentVertex->GetXYZ(vtx);
+  Float_t vtxf[3];
+  for(Int_t i=0;i<3;i++)vtxf[i]=vtx[i];
+  multReco->SetHistOn(kFALSE);
+  multReco->Reconstruct(itsClusterTree,vtxf,vtxf);
+  cout<<"======================================================="<<endl;
+  cout<<"Event number "<<evnumber<<"; tracklets= "<<multReco->GetNTracklets()<<endl;
+  Int_t notracks=multReco->GetNTracklets();
+  Float_t *trk = new Float_t [notracks];
+  Float_t *phi = new Float_t [notracks];
+  Float_t *dphi = new Float_t [notracks];
+  for(Int_t i=0;i<multReco->GetNTracklets();i++){
+    trk[i] = multReco->GetTracklet(i)[0];
+    phi[i] =  multReco->GetTracklet(i)[1];
+    dphi[i] = multReco->GetTracklet(i)[2];
+  }
+  fMult = new AliMultiplicity(notracks,trk,phi, dphi);
+  delete [] trk;
+  delete [] phi;
+  delete [] dphi;
+  for(Int_t i=0;i<multReco->GetNTracklets();i++){
+    cout<<i<<") theta= "<<fMult->GetTheta(i)<<", phi= "<<fMult->GetPhi(i)<<", DeltaPhi= "<<fMult->GetDeltaPhi(i)<<endl;
+  }
+  itsLoader->UnloadRecPoints();
+  delete multReco;
+  return;
+}
 
 //______________________________________________________________________
 void AliITSVertexer::WriteCurrentVertex(){
index 58a110b..01ac201 100644 (file)
@@ -21,6 +21,7 @@ class AliITSVertexer : public AliVertexer {
     // standard constructor     
     AliITSVertexer(TString filename); 
     virtual ~AliITSVertexer(){;}
+    virtual void FindMultiplicity(Int_t evnumber);
     virtual void WriteCurrentVertex();
  
 
index 0b9244a..8404c02 100644 (file)
@@ -14,9 +14,8 @@
  **************************************************************************/
 #include <AliITSVertexerZ.h>
 #include <TString.h>
-#include "AliITSLoader.h"
+#include <AliITSLoader.h>
 #include<TBranch.h>
-#include<TClonesArray.h>
 #include<TH1.h>
 #include<TTree.h>
 #include <AliITSgeom.h>
@@ -131,21 +130,20 @@ AliESDVertex* AliITSVertexerZ::FindVertexForCurrentEvent(Int_t evnumber){
       fDiffPhiMax=diffPhiMaxOrig;
     }
   }
-
+  FindMultiplicity(evnumber);
   return fCurrentVertex;
 }  
 
+
+
+
 //______________________________________________________________________
 void AliITSVertexerZ::VertexZFinder(Int_t evnumber){
   // Defines the AliESDVertex for the current event
   fCurrentVertex = 0;
   AliRunLoader *rl =AliRunLoader::GetRunLoader();
   AliITSLoader* itsLoader = (AliITSLoader*)rl->GetLoader("ITSLoader");
-  TDirectory * olddir = gDirectory;
-  rl->CdGAFile();
-  AliITSgeom* geom = (AliITSgeom*)gDirectory->Get("AliITSgeom");
-  olddir->cd();
-
+  AliITSgeom* geom = itsLoader->GetITSgeom();
   itsLoader->LoadRecPoints();
   rl->GetEvent(evnumber);
 
@@ -261,6 +259,8 @@ void AliITSVertexerZ::VertexZFinder(Int_t evnumber){
     }
     detTypeRec.ResetRecPoints();
   }
+  Int_t nolines=0;
   for(Int_t i=0;i<nrpL1;i++){ // loop on L1 RP
     Float_t r1=TMath::Sqrt(xc1[i]*xc1[i]+yc1[i]*yc1[i]); // radius L1 RP
     for(Int_t j=0;j<nrpL2;j++){ // loop on L2 RP
@@ -272,6 +272,15 @@ void AliITSVertexerZ::VertexZFinder(Int_t evnumber){
        fZCombf->Fill(zr0);
        fZCombc->Fill(zr0);
        fZCombv->Fill(zr0);
+       Double_t pA[3];
+       Double_t pB[3];
+       pA[0]=xc1[i];
+       pA[1]=yc1[i];
+       pA[2]=zc1[i];
+       pB[0]=xc2[j];
+       pB[1]=yc2[j];
+       pB[2]=zc2[j];
+       MakeTracklet(pA,pB,nolines);
       }
     }
   }
index 16ce67a..9d0945e 100644 (file)
@@ -45,6 +45,7 @@ class AliITSVertexerZ : public AliITSVertexer {
   Float_t GetBinWidthFine() const {return fStepFine;}
   void SetTolerance(Float_t tol = 20./10000.){fTolerance = tol;}
   Float_t GetTolerance() const {return fTolerance;}
+  virtual void MakeTracklet(Double_t * /* pA */, Double_t * /*pB */, Int_t & /* nolines */) {} // implemented in a derived class
 
  protected:
   AliITSVertexerZ(const AliITSVertexerZ& vtxr);
@@ -53,6 +54,7 @@ class AliITSVertexerZ : public AliITSVertexer {
   void VertexZFinder(Int_t evnumber);
   Float_t GetPhiMaxIter(Int_t i) const {return fPhiDiffIter[i];}
 
+
   Int_t fFirstL1;          // first module of the first pixel layer used
   Int_t fLastL1;           // last module of the first pixel layer used
   Int_t fFirstL2;          // first module of the second pixel layer used
index 9e09219..61fccb3 100644 (file)
@@ -44,6 +44,7 @@ AliESD::AliESD():
   fT0zVertex(0),
   fSPDVertex(),
   fPrimaryVertex(),
+  fSPDMult(),
   fT0timeStart(0),
   fTracks("AliESDtrack",15000),
   fHLTConfMapTracks("AliESDHLTtrack",25000),
@@ -127,6 +128,7 @@ void AliESD::Reset()
   fT0timeStart = 0;
   new (&fSPDVertex) AliESDVertex();
   new (&fPrimaryVertex) AliESDVertex();
+  new (&fSPDMult) AliMultiplicity();
   fTracks.Clear();
   fHLTConfMapTracks.Clear();
   fHLTHoughTracks.Clear();
@@ -159,6 +161,8 @@ void AliESD::Print(Option_t *) const
           fPrimaryVertex.GetXv(), fPrimaryVertex.GetXRes(),
           fPrimaryVertex.GetYv(), fPrimaryVertex.GetYRes(),
           fPrimaryVertex.GetZv(), fPrimaryVertex.GetZRes());
+    printf("SPD Multiplicity. Number of tracklets %d \n",
+           fSPDMult.GetNumberOfTracklets());
   printf("Event from reconstruction version %d \n",fRecoVersion);
   printf("Number of tracks: \n");
   printf("                 charged   %d\n", GetNumberOfTracks());
index e59e92b..3627373 100644 (file)
@@ -29,6 +29,7 @@
 #include "AliESDv0.h"
 #include "AliESDV0MI.h"
 #include "AliESDFMD.h"
+#include "AliMultiplicity.h"
 
 class AliESDfriend;
 
@@ -135,6 +136,11 @@ public:
   }
   const AliESDVertex *GetVertex() const {return &fSPDVertex;}
 
+  void SetMultiplicity(const AliMultiplicity *mul) {
+     new (&fSPDMult) AliMultiplicity(*mul);
+  }
+  const AliMultiplicity *GetMultiplicity() const {return &fSPDMult;}
+
   void SetPrimaryVertex(const AliESDVertex *vertex) {
      new (&fPrimaryVertex) AliESDVertex(*vertex);
   }
@@ -223,6 +229,7 @@ protected:
   Float_t      fT0zVertex;       // vertex z position estimated by the START
   AliESDVertex fSPDVertex;       // Primary vertex estimated by the SPD
   AliESDVertex fPrimaryVertex;   // Primary vertex estimated using ESD tracks
+  AliMultiplicity fSPDMult;      // SPD tracklet multiplicity
 
   Float_t      fT0timeStart;     // interaction time estimated by the START
   Float_t      fT0time[24];      // best TOF on each START PMT
@@ -247,7 +254,7 @@ protected:
  
   AliESDFMD *  fESDFMD; // FMD object containing rough multiplicity
 
-  ClassDef(AliESD,12)  //ESD class 
+  ClassDef(AliESD,13)  //ESD class 
 };
 #endif 
 
diff --git a/STEER/AliMultiplicity.cxx b/STEER/AliMultiplicity.cxx
new file mode 100644 (file)
index 0000000..65b9c62
--- /dev/null
@@ -0,0 +1,85 @@
+#include "AliMultiplicity.h"
+
+ClassImp(AliMultiplicity)
+
+//______________________________________________________________________
+AliMultiplicity::AliMultiplicity():TObject() {
+  // Default Constructor
+  fNtracks = 0;
+  fTh = 0;
+  fPhi = 0;
+  fDeltPhi = 0;
+}
+
+//______________________________________________________________________
+AliMultiplicity::AliMultiplicity(Int_t ntr, Float_t *t,  Float_t *ph, Float_t *df):TObject() {
+// Standard constructor
+  fNtracks = ntr;
+  if(ntr>0){
+    fTh = new Float_t [ntr];
+    fPhi = new Float_t [ntr];
+    fDeltPhi = new Float_t [ntr];
+    for(Int_t i=0;i<fNtracks;i++){
+      fTh[i]=t[i];
+      fPhi[i]=ph[i];
+      fDeltPhi[i]=df[i];
+    }
+  }
+  else {
+    fTh = 0;
+    fPhi = 0;
+    fDeltPhi = 0;
+  }
+}
+
+//______________________________________________________________________
+AliMultiplicity::AliMultiplicity(const AliMultiplicity& m):TObject(m){
+  // copy constructor
+
+  Duplicate(m);
+
+}
+
+//______________________________________________________________________
+AliMultiplicity &AliMultiplicity::operator=(const AliMultiplicity& m){
+  // assignment operator
+  if(this == &m)return *this;
+  ((TObject *)this)->operator=(m);
+
+  if (fTh)delete [] fTh;
+  if(fPhi)delete [] fPhi;
+  if(fDeltPhi)delete [] fDeltPhi;
+
+  Duplicate(m);
+
+  return *this;
+}
+
+//______________________________________________________________________
+void AliMultiplicity::Duplicate(const AliMultiplicity& m){
+  // used by copy constructor and assignment operator
+  fNtracks = m.fNtracks;
+  if(fNtracks>0){
+    fTh = new Float_t [fNtracks];
+    fPhi = new Float_t [fNtracks];
+    fDeltPhi = new Float_t [fNtracks];
+  }
+  else {
+    fTh = 0;
+    fPhi = 0;
+    fDeltPhi = 0;
+  }
+  memcpy(fTh,m.fTh,fNtracks*sizeof(Float_t));
+  memcpy(fPhi,m.fPhi,fNtracks*sizeof(Float_t));
+  memcpy(fDeltPhi,m.fDeltPhi,fNtracks*sizeof(Float_t));
+}
+
+//______________________________________________________________________
+AliMultiplicity::~AliMultiplicity(){
+  // Destructor
+  if (fTh)delete [] fTh;
+  if(fPhi)delete [] fPhi;
+  if(fDeltPhi)delete [] fDeltPhi;
+}
+
+
diff --git a/STEER/AliMultiplicity.h b/STEER/AliMultiplicity.h
new file mode 100644 (file)
index 0000000..2e8e474
--- /dev/null
@@ -0,0 +1,39 @@
+#ifndef ALIMULTIPLICITY_H
+#define ALIMULTIPLICITY_H
+
+#include<TObject.h>
+
+////////////////////////////////////////////////////////
+////   Class containing multiplicity information      //
+////   to stored in the ESD                           //
+////////////////////////////////////////////////////////
+
+class AliMultiplicity : public TObject {
+
+ public:
+
+  AliMultiplicity();               // default constructor
+  // standard constructor
+  AliMultiplicity(Int_t ntr,Float_t *t=NULL, Float_t *ph=NULL, Float_t *df=NULL); 
+  AliMultiplicity(const AliMultiplicity& m);
+  AliMultiplicity& operator=(const AliMultiplicity& m);
+  virtual ~AliMultiplicity();
+  Int_t GetNumberOfTracklets() const {return fNtracks;}
+  Float_t GetTheta(Int_t i) const { if(i>=0 && i<fNtracks) {return fTh[i];}
+  else {Error("GetTheta","Invalid track number %d",i); return -9999.;}}
+  Float_t GetPhi(Int_t i) const { if(i>=0 && i<fNtracks) {return fPhi[i];}
+  else {Error("GetTheta","Invalid track number %d",i); return -9999.;}}
+  Float_t GetDeltaPhi(Int_t i) const {if(i>=0 && i<fNtracks) {return fDeltPhi[i];}
+  else {Error("GetDeltaPhi","Invalid track number %d",i); return -9999.;}}
+
+  protected:
+  void Duplicate(const AliMultiplicity &m);  // used by copy ctr.
+  Int_t fNtracks;         // Number of tracklets (=-1 when mult is not determined)
+  Float_t *fTh;           //[fNtracks] array with theta values
+  Float_t *fPhi;          //[fNtracks] array with phi values
+  Float_t *fDeltPhi;      //[fNtracks] array with delta phi values
+
+  ClassDef(AliMultiplicity,1);
+};
+
+#endif
index de43bf8..d26f1f4 100644 (file)
 #include "AliESD.h"
 #include "AliESDfriend.h"
 #include "AliESDVertex.h"
+#include "AliMultiplicity.h"
 #include "AliTracker.h"
 #include "AliVertexer.h"
 #include "AliVertexerTracks.h"
@@ -930,6 +931,10 @@ Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
     vertex = new AliESDVertex(vtxPos, vtxErr);
   }
   esd->SetVertex(vertex);
+  // if SPD multiplicity has been determined, it is stored in the ESD
+  AliMultiplicity *mult= fVertexer->GetMultiplicity();
+  if(mult)esd->SetMultiplicity(mult);
+
   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
     if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
   }  
index d8217f9..d3690d9 100644 (file)
@@ -31,7 +31,9 @@ ClassImp(AliVertexer)
 AliVertexer::AliVertexer() :
   fCurrentVertex(0),
   fFirstEvent(0),
-  fLastEvent(0)
+  fLastEvent(0),
+  fDebug(0),
+  fMult()
 {
   //
   // Default Constructor
@@ -62,10 +64,13 @@ AliVertexer& AliVertexer::operator=(const AliVertexer& /* vtxr */){
 //______________________________________________________________________
 AliVertexer::~AliVertexer() {
   // Default Destructor
+
+  if(fMult) delete fMult;
   // The objects pointed by the following pointers are not owned
   // by this class and are not deleted
 
     fCurrentVertex  = 0;
+
 }
 
 //______________________________________________________________________
index cb09afc..546acff 100644 (file)
@@ -2,6 +2,7 @@
 #define ALIVERTEXER_H
 
 #include<TObject.h>
+#include<AliMultiplicity.h>
 
 ///////////////////////////////////////////////////////////////////
 //                                                               //
@@ -25,8 +26,9 @@ class AliVertexer : public TObject {
     virtual ~AliVertexer(); 
     // computes the vertex for the current event
     virtual AliESDVertex* FindVertexForCurrentEvent(Int_t evnumb)=0; 
-    // computes the vetex for each event and stores it on file
+    // computes the vertex for each event and stores it on file
     virtual void FindVertices()= 0;
+    virtual AliMultiplicity* GetMultiplicity() const {return fMult;}
     virtual void PrintStatus() const = 0;
     virtual void SetDebug(Int_t debug = 0);
     virtual void SetFirstEvent(Int_t ev){fFirstEvent = ev;}
@@ -46,8 +48,9 @@ class AliVertexer : public TObject {
     Int_t fFirstEvent;          // First event to be processed by FindVertices
     Int_t fLastEvent;           // Last event to be processed by FindVertices 
     Int_t fDebug;               //! debug flag - verbose printing if >0
+    AliMultiplicity *fMult;     //! Multiplicity object
 
-  ClassDef(AliVertexer,1);
+  ClassDef(AliVertexer,2);
 };
 
 #endif
index fa1ae6a..d88d4b2 100644 (file)
@@ -54,6 +54,7 @@
 #pragma link C++ class  AliFMDFloatMap+;
 
 #pragma link C++ class  AliESDMultITS+;
+#pragma link C++ class  AliMultiplicity+;
 #pragma link C++ class  AliTracker+;
 #endif
 
index 1ed30b8..074ec1b 100644 (file)
@@ -11,7 +11,7 @@ SRCS = AliESD.cxx AliESDfriend.cxx\
        AliTrackPointArray.cxx AliCluster.cxx \
        AliESDFMD.cxx AliFMDMap.cxx AliFMDFloatMap.cxx \
        AliESDMultITS.cxx \
-       AliTracker.cxx
+       AliTracker.cxx AliMultiplicity.cxx
 
 
 HDRS:= $(SRCS:.cxx=.h)