]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/CheckESD.C
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff
? search:
summary | shortlog | log | commit | commitdiff | tree
raw | inline | side by side
Initialization of some data members (Alberto)
[u/mrichter/AliRoot.git] / STEER / CheckESD.C
diff --git a/STEER/CheckESD.C b/STEER/CheckESD.C
index ef0b4a9a2040dd1a1d2ea40eb9f14e8ce8601c97..a918f27dcadffdb184f1d9ef64ba23cd93b4f28f 100644 (file)
--- a/STEER/CheckESD.C
+++ b/STEER/CheckESD.C
@@ -7,6 +7,7 @@
 #include <TCanvas.h>
 #include <TVector3.h>
 #include <TPDGCode.h>
+#include <TParticle.h>
 
 #include "AliRunLoader.h"
 #include "AliLoader.h"
@@ -15,6 +16,7 @@
 #include "AliStack.h"
 #include "AliHeader.h"
 #include "AliGenEventHeader.h"
+#include "AliPID.h"
 #else
 const Int_t kXiMinus = 3312;
 const Int_t kOmegaMinus = 3334;
@@ -155,6 +157,13 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
     Error("CheckESD", "opening ESD file %s failed", esdFileName);
     return kFALSE;
   }
+  AliESD* esd = new AliESD;
+  TTree* tree = (TTree*) esdFile->Get("esdTree");
+  if (!tree) {
+    Error("CheckESD", "no ESD tree found");
+    return kFALSE;
+  }
+  tree->SetBranchAddress("ESD", &esd);
 
   // efficienc and resolution histograms
   Int_t nBinsPt = 15;
@@ -176,15 +185,15 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
                                "#theta_{rec}-#theta_{sim} [mrad]", "N");
 
   // PID
-  Int_t partCode[AliESDtrack::kSPECIES] = 
+  Int_t partCode[AliPID::kSPECIES] = 
     {kElectron, kMuonMinus, kPiPlus, kKPlus, kProton};
-  const char* partName[AliESDtrack::kSPECIES+1] = 
+  const char* partName[AliPID::kSPECIES+1] = 
     {"electron", "muon", "pion", "kaon", "proton", "other"};
-  Double_t partFrac[AliESDtrack::kSPECIES] = 
+  Double_t partFrac[AliPID::kSPECIES] = 
     {0.01, 0.01, 0.85, 0.10, 0.05};
-  Int_t identified[AliESDtrack::kSPECIES+1][AliESDtrack::kSPECIES];
-  for (Int_t iGen = 0; iGen < AliESDtrack::kSPECIES+1; iGen++) {
-    for (Int_t iRec = 0; iRec < AliESDtrack::kSPECIES; iRec++) {
+  Int_t identified[AliPID::kSPECIES+1][AliPID::kSPECIES];
+  for (Int_t iGen = 0; iGen < AliPID::kSPECIES+1; iGen++) {
+    for (Int_t iRec = 0; iRec < AliPID::kSPECIES; iRec++) {
       identified[iGen][iRec] = 0;
     }
   }
@@ -212,7 +221,7 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
 
   // calorimeters
   TH1F* hEPHOS = CreateHisto("hEPHOS", "PHOS", 100, 0, 5, "E [GeV]", "N");
-  TH1F* hEEMCAL = CreateHisto("hEEMCAL", "EMCAL", 100, 0, 2, "E [GeV]", "N");
+  TH1F* hEEMCAL = CreateHisto("hEEMCAL", "EMCAL", 100, 0, 50, "E [GeV]", "N");
 
   // muons
   TH1F* hPtMUON = CreateHisto("hPtMUON", "MUON", 100, 0, 20, 
@@ -291,9 +300,7 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
     }
 
     // get the event summary data
-    char esdName[256]; 
-    sprintf(esdName, "ESD%d", iEvent);
-    AliESD* esd = (AliESD*) esdFile->Get(esdName);
+    tree->GetEvent(iEvent);
     if (!esd) {
       Error("CheckESD", "no ESD object found for event %d", iEvent);
       return kFALSE;
@@ -328,14 +335,14 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
       // PID
       if ((track->GetStatus() & AliESDtrack::kESDpid) == 0) continue;
       Int_t iGen = 5;
-      for (Int_t i = 0; i < AliESDtrack::kSPECIES; i++) {
+      for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
        if (TMath::Abs(particle->GetPdgCode()) == partCode[i]) iGen = i;
       }
-      Double_t probability[AliESDtrack::kSPECIES];
+      Double_t probability[AliPID::kSPECIES];
       track->GetESDpid(probability);
       Double_t pMax = 0;
       Int_t iRec = 0;
-      for (Int_t i = 0; i < AliESDtrack::kSPECIES; i++) {
+      for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
        probability[i] *= partFrac[i];
        if (probability[i] > pMax) {
          pMax = probability[i];
@@ -346,7 +353,7 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
       if (iGen == iRec) nIdentified++;
 
       // dE/dx and TOF
-      Double_t time[AliESDtrack::kSPECIES];
+      Double_t time[AliPID::kSPECIES];
       track->GetIntegratedTimes(time);
       if (iGen == iRec) {
        hDEdxRight->Fill(pTrack.Mag(), track->GetTPCsignal());
@@ -361,20 +368,6 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
       }
     }
 
-    // loop over calo tracks
-    for (Int_t iTrack = 0; iTrack < esd->GetNumberOfCaloTracks(); iTrack++) {
-      AliESDCaloTrack* caloTrack = esd->GetCaloTrack(iTrack);
-      TParticle* recParticle = caloTrack->GetRecParticle();
-      if (recParticle->InheritsFrom("AliPHOSRecParticle")) {
-       hEPHOS->Fill(recParticle->Energy());
-      } else if (recParticle->InheritsFrom("AliEMCALRecParticle")) {
-       hEEMCAL->Fill(recParticle->Energy());
-      } else {
-       Warning("CheckESD", "unknown calo particle");
-       recParticle->Dump();
-      }
-    }
-
     // loop over muon tracks
     {
     for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
@@ -389,18 +382,18 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
     // loop over V0s
     for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
       AliESDv0* v0 = esd->GetV0(iV0);
-      switch (v0->GetPdgCode()) {
-      case kK0Short   : hMassK0->Fill(v0->GetEffMass()); break;
-      case kLambda0   : hMassLambda->Fill(v0->GetEffMass()); break;
-      case kLambda0Bar: hMassLambdaBar->Fill(v0->GetEffMass()); break;
-      default   : break;
-      }
+      v0->ChangeMassHypothesis(kK0Short);
+      hMassK0->Fill(v0->GetEffMass());
+      v0->ChangeMassHypothesis(kLambda0);
+      hMassLambda->Fill(v0->GetEffMass());
+      v0->ChangeMassHypothesis(kLambda0Bar);
+      hMassLambdaBar->Fill(v0->GetEffMass());
 
       Int_t negLabel = TMath::Abs(esd->GetTrack(v0->GetNindex())->GetLabel());
       if (negLabel > stack->GetNtrack()) continue;     // background
       Int_t negMother = stack->Particle(negLabel)->GetMother(0);
       if (negMother < 0) continue;
-      Int_t posLabel = TMath::Abs(esd->GetTrack(v0->GetNindex())->GetLabel());
+      Int_t posLabel = TMath::Abs(esd->GetTrack(v0->GetPindex())->GetLabel());
       if (posLabel > stack->GetNtrack()) continue;     // background
       Int_t posMother = stack->Particle(posLabel)->GetMother(0);
       if (negMother != posMother) continue;
@@ -414,18 +407,18 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
     for (Int_t iCascade = 0; iCascade < esd->GetNumberOfCascades(); 
         iCascade++) {
       AliESDcascade* cascade = esd->GetCascade(iCascade);
-      switch (TMath::Abs(cascade->GetPdgCode())) {
-      case kXiMinus   : hMassXi->Fill(cascade->GetEffMass()); break;
-      case kOmegaMinus: hMassOmega->Fill(cascade->GetEffMass()); break;
-      default   : break;
-      }
+      Double_t v0q;
+      cascade->ChangeMassHypothesis(v0q,kXiMinus);
+      hMassXi->Fill(cascade->GetEffMass());
+      cascade->ChangeMassHypothesis(v0q,kOmegaMinus);
+      hMassOmega->Fill(cascade->GetEffMass());
 
       Int_t negLabel = TMath::Abs(esd->GetTrack(cascade->GetNindex())
                                  ->GetLabel());
       if (negLabel > stack->GetNtrack()) continue;     // background
       Int_t negMother = stack->Particle(negLabel)->GetMother(0);
       if (negMother < 0) continue;
-      Int_t posLabel = TMath::Abs(esd->GetTrack(cascade->GetNindex())
+      Int_t posLabel = TMath::Abs(esd->GetTrack(cascade->GetPindex())
                                  ->GetLabel());
       if (posLabel > stack->GetNtrack()) continue;     // background
       Int_t posMother = stack->Particle(posLabel)->GetMother(0);
@@ -442,6 +435,22 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
       selCascades.Remove(particle);
       nRecCascades++;
     }
+
+    // loop over the PHOS clusters
+    {
+    Int_t firstPHOSCluster = esd->GetFirstPHOSCluster();
+    Int_t lastPHOSCluster  = firstPHOSCluster + esd->GetNumberOfPHOSClusters();
+    for (Int_t iCluster=firstPHOSCluster; iCluster<lastPHOSCluster; iCluster++)
+      hEPHOS->Fill(esd->GetCaloCluster(iCluster)->GetClusterEnergy());
+    }
+
+    // loop over the EMCAL clusters
+    {
+    Int_t firstEMCALCluster = esd->GetFirstEMCALCluster();
+    Int_t lastEMCALCluster  = firstEMCALCluster + esd->GetNumberOfEMCALClusters();
+    for (Int_t iCluster=firstEMCALCluster; iCluster<lastEMCALCluster; iCluster++)
+      hEEMCAL->Fill(esd->GetCaloCluster(iCluster)->GetClusterEnergy());
+    }
   }
 
   // perform checks
@@ -512,13 +521,13 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
     }
 
     printf("%9s:", "gen\\rec");
-    for (Int_t iRec = 0; iRec < AliESDtrack::kSPECIES; iRec++) {
+    for (Int_t iRec = 0; iRec < AliPID::kSPECIES; iRec++) {
       printf("%9s", partName[iRec]);
     }
     printf("\n");
-    for (Int_t iGen = 0; iGen < AliESDtrack::kSPECIES+1; iGen++) {
+    for (Int_t iGen = 0; iGen < AliPID::kSPECIES+1; iGen++) {
       printf("%9s:", partName[iGen]);
-      for (Int_t iRec = 0; iRec < AliESDtrack::kSPECIES; iRec++) {
+      for (Int_t iRec = 0; iRec < AliPID::kSPECIES; iRec++) {
        printf("%9d", identified[iGen][iRec]);
       }
       printf("\n");
@@ -677,6 +686,7 @@ Bool_t CheckESD(const char* gAliceFileName = "galice.root",
   delete hMassXi;
   delete hMassOmega;
 
+  delete esd;
   esdFile->Close();
   delete esdFile;