Centrality Task (A. Toia)
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 29 Oct 2010 09:13:48 +0000 (09:13 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 29 Oct 2010 09:13:48 +0000 (09:13 +0000)
ANALYSIS/AliCentralitySelectionTask.cxx [new file with mode: 0644]

diff --git a/ANALYSIS/AliCentralitySelectionTask.cxx b/ANALYSIS/AliCentralitySelectionTask.cxx
new file mode 100644 (file)
index 0000000..3f4f8e0
--- /dev/null
@@ -0,0 +1,288 @@
+/**************************************************************************
+ * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+//*****************************************************
+//   Class AliCentralitySelectionTask
+//   Class to analyze determine centrality            
+//   author: Alberica Toia
+//*****************************************************
+
+#include <TTree.h>
+#include <TList.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TProfile.h>
+#include <TFile.h>
+#include <TString.h>
+#include <TCanvas.h>
+#include <iostream>
+
+#include "AliAnalysisManager.h"
+#include "AliVEvent.h"
+#include "AliESD.h"
+#include "AliESDEvent.h"
+#include "AliESDHeader.h"
+#include "AliESDInputHandler.h"
+#include "AliESDZDC.h"
+#include "AliESDFMD.h"
+#include "AliESDVZERO.h"
+#include "AliESDCentrality.h"
+#include "AliMultiplicity.h"
+#include "AliAODHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODVertex.h"
+#include "AliAODMCHeader.h"
+#include "AliMCEvent.h"
+#include "AliMCEventHandler.h"
+#include "AliMCParticle.h"
+#include "AliStack.h"
+#include "AliHeader.h"
+#include "AliAODMCParticle.h"
+#include "AliAnalysisTaskSE.h"
+#include "AliGenEventHeader.h"
+#include "AliGenHijingEventHeader.h"
+#include "AliPhysicsSelectionTask.h"
+#include "AliPhysicsSelection.h"
+#include "AliBackgroundSelection.h"
+#include "AliCentralitySelectionTask.h"
+
+ClassImp(AliCentralitySelectionTask)
+
+
+//________________________________________________________________________
+AliCentralitySelectionTask::AliCentralitySelectionTask():
+AliAnalysisTaskSE(),
+  fDebug(0),
+  fAnalysisInput("ESD"),
+  fIsMCInput(kFALSE),
+  fFile(0),
+  fCentfilename(0),
+  fMethod(0),
+  fCent(0),
+  fHtemp(0)
+{   
+  // Default constructor
+}   
+
+//________________________________________________________________________
+AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
+  AliAnalysisTaskSE(name),
+  fDebug(0),
+  fAnalysisInput("ESD"),
+  fIsMCInput(kFALSE),
+  fFile(0),
+  fCentfilename(0),
+  fMethod(0),
+  fCent(0),
+  fHtemp(0)
+{
+  // Default constructor
+}
+
+//________________________________________________________________________
+AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
+{
+  // Assignment operator
+  if (this!=&c) {
+    AliAnalysisTaskSE::operator=(c);
+  }
+  return *this;
+}
+
+//________________________________________________________________________
+AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
+  AliAnalysisTaskSE(ana),
+  fDebug(ana.fDebug),    
+  fAnalysisInput(ana.fDebug),
+  fIsMCInput(ana.fIsMCInput),
+  fFile(ana.fFile),
+  fCentfilename(ana.fCentfilename),
+  fMethod(ana.fMethod),
+  fCent(ana.fCent),
+  fHtemp(ana.fHtemp)
+{
+  // Copy Constructor  
+}
+//________________________________________________________________________
+ AliCentralitySelectionTask::~AliCentralitySelectionTask()
+ {
+   // Destructor
+ }  
+
+//________________________________________________________________________
+void AliCentralitySelectionTask::UserCreateOutputObjects()
+{  
+  // Create the output containers
+  if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
+}
+
+//________________________________________________________________________
+void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
+{ 
+  // Execute analysis for current event:
+  if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
+  
+  Float_t  zncEnergy;                 //  ZNC Energy
+  Float_t  zpcEnergy;                 //  ZPC Energy
+  Float_t  znaEnergy;                 //  ZNA Energy
+  Float_t  zpaEnergy;                 //  ZPA Energy
+  Float_t  zem1Energy;                //  ZEM1 Energy
+  Float_t  zem2Energy;                //  ZEM2 Energy
+  
+  Int_t    nTracks    = 0;            //  no. tracks
+  Int_t    nTracklets = 0;            //  no. tracklets
+  Int_t    nClusters[6];              //  no. clusters on 6 ITS layers
+  Int_t    nChips[2];                 //  no. chips on 2 SPD layers
+
+  Float_t  multV0A    = 0;            //  multiplicity from V0 reco side A
+  Float_t  multV0C    = 0;            //  multiplicity from V0 reco side C
+  Float_t  multFMDA   = 0;            //  multiplicity from FMD on detector A
+  Float_t  multFMDC   = 0;            //  multiplicity from FMD on detector C
+
+  AliESDCentrality *esdCent = 0;
+
+  if(fAnalysisInput.CompareTo("ESD")==0){
+
+    AliVEvent* event = InputEvent();
+    AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
+    
+    esdCent = esd->GetCentrality();
+
+    // ***** V0 info    
+    AliESDVZERO* esdV0 = esd->GetVZEROData();
+    multV0A=esdV0->GetMTotV0A();
+    multV0C=esdV0->GetMTotV0C();
+    
+    // ***** CB info (tracklets, clusters, chips)
+    nTracks    = event->GetNumberOfTracks();     
+
+    const AliMultiplicity *mult = esd->GetMultiplicity();
+
+    nTracklets = mult->GetNumberOfTracklets();
+
+    for(Int_t ilay=0; ilay<6; ilay++){
+      nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
+    }
+
+    for(Int_t ilay=0; ilay<2; ilay++){
+      nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
+    }
+    
+
+    // ***** FMD info
+    AliESDFMD *fmd = esd->GetFMDData();
+    Float_t totalMultA = 0;
+    Float_t totalMultC = 0;
+    const Float_t fFMDLowCut = 0.4;
+    
+    for(UShort_t det=1;det<=3;det++) {
+      Int_t nRings = (det==1 ? 1 : 2);
+      for (UShort_t ir = 0; ir < nRings; ir++) {         
+       Char_t   ring = (ir == 0 ? 'I' : 'O');
+       UShort_t nsec = (ir == 0 ? 20  : 40);
+       UShort_t nstr = (ir == 0 ? 512 : 256);
+       for(UShort_t sec =0; sec < nsec;  sec++)  {
+         for(UShort_t strip = 0; strip < nstr; strip++) {
+           
+           Float_t FMDmult = fmd->Multiplicity(det,ring,sec,strip);
+           if(FMDmult == 0 || FMDmult == AliESDFMD::kInvalidMult) continue;
+           
+           Float_t nParticles=0;
+           
+           if(FMDmult > fFMDLowCut) {
+             nParticles = 1.;
+           }
+           
+           if (det<3) totalMultA = totalMultA + nParticles;
+           else totalMultC = totalMultC + nParticles;
+           
+         }
+       }
+      }
+    }
+    multFMDA = totalMultA;
+    multFMDC = totalMultC;
+    
+    // ***** ZDC info
+    AliESDZDC *esdZDC = esd->GetESDZDC();
+    zncEnergy  = (Float_t) (esdZDC->GetZDCN1Energy());
+    zpcEnergy  = (Float_t) (esdZDC->GetZDCP1Energy());
+    znaEnergy  = (Float_t) (esdZDC->GetZDCN2Energy());
+    zpaEnergy  = (Float_t) (esdZDC->GetZDCP2Energy());
+    zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0));
+    zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1));
+    
+  }   
+  else if(fAnalysisInput.CompareTo("AOD")==0){
+    //AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
+    // to be implemented
+    printf("  AOD analysis not yet implemented!!!\n\n");
+    return;
+  }
+
+   // ***** Centrality Selection
+  if(fMethod.CompareTo("V0")==0){
+    fCent = fHtemp->GetBinContent(fHtemp->FindBin((multV0A+multV0C)));
+  }
+  if(fMethod.CompareTo("FMD")==0){
+    fCent = fHtemp->GetBinContent(fHtemp->FindBin((multFMDA+multFMDC)));
+  }
+  if(fMethod.CompareTo("TRACKS")==0) {
+    fCent = fHtemp->GetBinContent(fHtemp->FindBin(nTracks));
+  }
+  if(fMethod.CompareTo("TRACKLETS")==0){
+    fCent = fHtemp->GetBinContent(fHtemp->FindBin(nTracklets));
+  }
+  if(fMethod.CompareTo("CLUSTERS")==0) {
+    fCent = fHtemp->GetBinContent(fHtemp->FindBin(nClusters[0]));
+  }
+  printf(" **** centrality is %3.2f \n", fCent);
+  
+  esdCent->SetCentrality(fCent);
+}
+
+//________________________________________________________________________
+void AliCentralitySelectionTask::SetCentralityMethod(const char* x) 
+{
+  fMethod = x;
+
+  fFile  = new TFile(fCentfilename);
+  
+  if(fMethod.CompareTo("V0")==0)
+    fHtemp  = (TH1D*) (fFile->Get("hmultV0_percentile")); 
+  
+  if(fMethod.CompareTo("FMD")==0)
+    fHtemp  = (TH1D*) (fFile->Get("hmultFMD_percentile")); 
+  
+  if(fMethod.CompareTo("TRACKS")==0) 
+    fHtemp  = (TH1D*) (fFile->Get("hNtracks_percentile")); 
+  
+  if(fMethod.CompareTo("TRACKLETS")==0)
+    fHtemp  = (TH1D*) (fFile->Get("hNtracklets_percentile")); 
+  
+  if(fMethod.CompareTo("CLUSTERS")==0) 
+    fHtemp  = (TH1D*) (fFile->Get("hNclusters0_percentile")); 
+  
+}
+
+//________________________________________________________________________
+void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
+{
+  // Terminate analysis
+  fFile->Close();  
+}
+//________________________________________________________________________
+