]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/CreateAODfromKineTree.C
ITS cluster multiplicity and TPC standalone multiplicity in AODHeader
[u/mrichter/AliRoot.git] / STEER / CreateAODfromKineTree.C
index 62a1056816bd83b78ed0da201c7f5d073e260415..5271ffb4fbd52798d6fe267b9853baab005b41ff 100644 (file)
@@ -55,6 +55,15 @@ Int_t nK0 = 0;
 Int_t nKPlus = 0;
 Int_t nKMinus = 0;
 
+// global arrays and pointers
+Float_t p[3];
+Float_t x[3];
+Float_t *cov = NULL; // set to NULL because not provided
+Float_t *pid = NULL; // set to NULL because not provided
+AliAODVertex *primary = NULL; 
+AliAODVertex *secondary = NULL;
+AliAODTrack *currTrack = NULL;
+
 void CreateAODfromKineTree(const char *inFileName,
                           const char *outFileName) {
   printf("This macro works only correctly in comiled mode!\n");
@@ -68,7 +77,7 @@ void CreateAODfromKineTree(const char *inFileName,
   TFile *outFile = TFile::Open(outFileName, "RECREATE");
 
   // create the tree
-  TTree *aodTree = new TTree("AOD", "AliAOD tree");
+  TTree *aodTree = new TTree("aodTree", "AliAOD tree");
   aodTree->Branch(aod->GetList());
 
   AliRunLoader *runLoader;
@@ -112,28 +121,32 @@ void CreateAODfromKineTree(const char *inFileName,
     nKPlus = 0;
     nKMinus = 0;
 
-    // create the header
-    aod->AddHeader(new AliAODHeader(aliHeader->GetRun(),
-                                   0, // bunchX number
-                                   0, // orbit number
-                                   nTracks,
-                                   nPos,
-                                   nNeg,
-                                   -999, // mag. field
-                                   -999., // muon mag. field
-                                   -999., // centrality
-                                   -999, // ZDCN1Energy
-                                   -999, // ZDCP1Energy
-                                   -999, // ZDCN2Energy
-                                   -999, // ZDCP2Energy
-                                   -999, // ZDCEMEnergy
-                                   0, // TriggerMask
-                                   0, // TriggerCluster
-                                   0)); // EventType
-    
     // Access to the header
     AliAODHeader *header = aod->GetHeader();
 
+    Double_t emEnergy[2] = {-999., -999.};
+
+    // fill the header
+    *header = AliAODHeader(aliHeader->GetRun(),
+                          0, // bunchX number
+                          0, // orbit number
+                          0, // period number
+                          nTracks,
+                          nPos,
+                          nNeg,
+                          -999, // mag. field
+                          -999., // muon mag. field
+                          -999., // centrality
+                          -999., // ZDCN1Energy
+                          -999., // ZDCP1Energy
+                          -999., // ZDCN2Energy
+                          -999., // ZDCP2Energy
+                          emEnergy, // emEnergy
+                          0, // TriggerMask
+                          0, // TriggerCluster
+                          0, // EventType
+                          ""); // title
+  
     // Access to the AOD container of vertices
     TClonesArray &vertices = *(aod->GetVertices());
     jVertices=0;
@@ -148,12 +161,6 @@ void CreateAODfromKineTree(const char *inFileName,
  
     aod->ResetStd(nTracks, 1);
 
-    Float_t p[3];
-    Float_t x[3];
-    Float_t *covTr = NULL;
-    Float_t *pid = NULL;
-    AliAODVertex *primary = NULL; 
-    AliAODTrack* currTrack = NULL;
     // track loop
     for (Int_t iTrack = 0; iTrack < nPrims; ++iTrack) {
                                                                
@@ -177,11 +184,12 @@ void CreateAODfromKineTree(const char *inFileName,
                                                                kTRUE,
                                                                x,
                                                                kFALSE,
-                                                               covTr
+                                                               cov, 
                                                                (Short_t)-99,
                                                                0, // no ITSClusterMap
                                                                pid,
                                                                primary,
+                                                               kFALSE,  // no fit performed
                                                                kFALSE, // no fit preformed
                                                                AliAODTrack::kPrimary));
        currTrack = (AliAODTrack*)tracks.Last();
@@ -233,53 +241,50 @@ Int_t LoopOverSecondaries(TParticle *mother) {
     TClonesArray &vertices = *(aod->GetVertices());
     TClonesArray &tracks = *(aod->GetTracks());
 
-    Float_t pSec[3];
-    Float_t xSec[3];
-    Float_t *covSec = NULL;
-    Float_t *pidSec = NULL;
-    AliAODVertex *secondary = NULL;
-    AliAODTrack* currTrackSec = NULL;
-    
     for (Int_t iDaughter = mother->GetFirstDaughter(); iDaughter <= mother->GetLastDaughter(); iDaughter++) {
-      TParticle *partSec = stack->Particle(iDaughter);
-
-      pSec[0] = partSec->Px(); pSec[1] = partSec->Py(); pSec[2] = partSec->Pz();
-      xSec[0] = partSec->Vx(); xSec[1] = partSec->Vy(); xSec[2] = partSec->Vz();
+      TParticle *part = stack->Particle(iDaughter);
+      p[0] = part->Px(); 
+      p[1] = part->Py(); 
+      p[2] = part->Pz();
+      x[0] = part->Vx(); 
+      x[1] = part->Vy(); 
+      x[2] = part->Vz();
 
       if (iDaughter == mother->GetFirstDaughter()) {
        // add secondary vertex
        secondary = new(vertices[jVertices++])
-         AliAODVertex(xSec, NULL, -999., tracks.Last(), AliAODVertex::kUndef);
+         AliAODVertex(x, NULL, -999., tracks.Last(), AliAODVertex::kUndef);
        
-       SetVertexType(partSec, secondary);
+       SetVertexType(part, secondary);
       }
        
       // add secondary tracks
       secondary->AddDaughter(new(tracks[jTracks++]) AliAODTrack(0, // ID
                                                                0, // label
-                                                               pSec,
+                                                               p,
                                                                kTRUE,
-                                                               xSec,
+                                                               x,
                                                                kFALSE,
-                                                               covSec
+                                                               cov, 
                                                                (Short_t)-99,
                                                                0, // no cluster map available
-                                                               pidSec,
+                                                               pid,
                                                                secondary,
                                                                kFALSE, // no fit performed
+                                                               kFALSE, // no fit performed
                                                                AliAODTrack::kSecondary));
 
-      currTrackSec = (AliAODTrack*)tracks.Last();
-      SetChargeAndPID(partSec->GetPdgCode(), currTrackSec);
-      if (currTrackSec->Charge() != -99) {
-       if (currTrackSec->Charge() > 0) {
+      currTrack = (AliAODTrack*)tracks.Last();
+      SetChargeAndPID(part->GetPdgCode(), currTrack);
+      if (currTrack->Charge() != -99) {
+       if (currTrack->Charge() > 0) {
          nPos++;
-       } else if (currTrackSec->Charge() < 0) {
+       } else if (currTrack->Charge() < 0) {
          nNeg++;
        }           
       }
 
-      LoopOverSecondaries(partSec);
+      LoopOverSecondaries(part);
     }
     return 1;
   } else {
@@ -293,194 +298,249 @@ void SetChargeAndPID(Int_t pdgCode, AliAODTrack *track) {
   Float_t PID[10] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
 
   switch (pdgCode) {
+
   case 22: // gamma
     track->SetCharge(0);
-    PID[5] = 1.;
+    PID[AliAODTrack::kUnknown] = 1.;
     track->SetPID(PID);
     nGamma++;
     break;
 
-  case -11: // e- 
+  case 11: // e- 
     track->SetCharge(-1);
-    PID[0] = 1.;
+    PID[AliAODTrack::kElectron] = 1.;
     track->SetPID(PID);
     nElectron++;
     break;
     
-  case 11: // e+
+  case -11: // e+
     track->SetCharge(+1);
-    PID[0] = 1.;
+    PID[AliAODTrack::kElectron] = 1.;
     track->SetPID(PID);
     nPositron++;
     break;
     
-  case -13: // mu- 
+  case 13: // mu- 
     track->SetCharge(-1);
-    PID[1] = 1.;
+    PID[AliAODTrack::kMuon] = 1.;
     track->SetPID(PID);
     nMuon++;
     break;
     
-  case +13: // mu+
+  case -13: // mu+
     track->SetCharge(+1);
-    PID[1] = 1.;
+    PID[AliAODTrack::kMuon] = 1.;
     track->SetPID(PID);
     nAntiMuon++;
     break;
     
   case 111: // pi0
     track->SetCharge(0);
-    PID[6] = 1.;
+    PID[AliAODTrack::kUnknown] = 1.;
     track->SetPID(PID);
     nPi0++;
     break;
     
   case 211: // pi+
     track->SetCharge(+1);
-    PID[2] = 1.;
+    PID[AliAODTrack::kPion] = 1.;
     track->SetPID(PID);
     nPiPlus++;
     break;
     
   case -211: // pi-
     track->SetCharge(-1);
-    PID[2] = 1.;
+    PID[AliAODTrack::kPion] = 1.;
     track->SetPID(PID);
     nPiMinus++;
     break;
     
   case 130: // K0L
     track->SetCharge(0);
-    PID[8] = 1.;
+    PID[AliAODTrack::kUnknown] = 1.;
     track->SetPID(PID);
     nK0++;
     break;
     
   case 321: // K+
     track->SetCharge(+1);
-    PID[3] = 1.;
+    PID[AliAODTrack::kKaon] = 1.;
     track->SetPID(PID);
     nKPlus++;
     break;
     
   case -321: // K- 
     track->SetCharge(-1);
-    PID[3] = 1.;
+    PID[AliAODTrack::kKaon] = 1.;
     track->SetPID(PID);
     nKMinus++;
     break;
     
-  case 2122: // n
+  case 2112: // n
     track->SetCharge(0);
-    PID[7] = 1.;
+    PID[AliAODTrack::kUnknown] = 1.;
     track->SetPID(PID);
     nNeutron++;
     break;
     
   case 2212: // p
     track->SetCharge(+1);
-    PID[4] = 1.;
+    PID[AliAODTrack::kProton] = 1.;
     track->SetPID(PID);
     nProton++;
     break;
     
   case -2212: // anti-p
     track->SetCharge(-1);
-    PID[4] = 1.;
+    PID[AliAODTrack::kProton] = 1.;
     track->SetPID(PID);
     nAntiProton++;
     break;
 
   case 310: // K0S
     track->SetCharge(0);
-    PID[8] = 1.;
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
+    nK0++;
+    break;
+    
+  case 311: // K0
+    track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
+    nK0++;
+    break;
+    
+  case -311: // anti-K0
+    track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
     track->SetPID(PID);
     nK0++;
     break;
     
   case 221: // eta
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3122: // lambda
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3222: // Sigma+
     track->SetCharge(+1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3212: // Sigma0
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3112: // Sigma-
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3322: // Xi0
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 3312: // Xi-
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
-  case 3332: // Omega-
+  case 3334: // Omega-
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -2112: // n-bar
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3122: // anti-Lambda
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3222: // anti-Sigma-
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3212: // anti-Sigma0
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3112: // anti-Sigma+
     track->SetCharge(+1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3322: // anti-Xi0
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -3312: // anti-Xi+
     track->SetCharge(+1);
     break;
 
-  case -3332: // anti-Omega+
+  case -3334: // anti-Omega+
     track->SetCharge(+1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 411: // D+
     track->SetCharge(+1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -411: // D- 
     track->SetCharge(-1);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case 421: // D0
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   case -421: // anti-D0
     track->SetCharge(0);
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
     break;
 
   default : // unknown
     track->SetCharge(-99);
-  }
+    PID[AliAODTrack::kUnknown] = 1.;
+    track->SetPID(PID);
+ }
 
   return;
 }
@@ -499,29 +559,30 @@ void SetVertexType(TParticle *part, AliAODVertex *vertex) {
     Int_t lastPdgCode = stack->Particle(mother->GetLastDaughter())->GetPdgCode();
 
     if (!(pdgMother == 22 || pdgMother == 111 || pdgMother == 130 || 
-         TMath::Abs(pdgMother) == 2122 || pdgMother == 310 || pdgMother == 221 || 
+         TMath::Abs(pdgMother) == 2112 || pdgMother == 310 || pdgMother == 221 || 
          TMath::Abs(pdgMother) == 3122 || TMath::Abs(pdgMother) == 3322 || 
-         pdgMother == -3212 || TMath::Abs(pdgMother) == 421) // not neutral
+         pdgMother == -3212 || TMath::Abs(pdgMother) == 421 || 
+         TMath::Abs(pdgMother) == 311) // not neutral
        && (((firstPdgCode == 22 || firstPdgCode == 111 || firstPdgCode == 130 || 
-             TMath::Abs(firstPdgCode) == 2122 || firstPdgCode == 310 || 
+             TMath::Abs(firstPdgCode) == 2112 || firstPdgCode == 310 || 
              firstPdgCode == 221 || TMath::Abs(firstPdgCode) == 3122 || 
              TMath::Abs(firstPdgCode) == 3322 || firstPdgCode == -3212 || 
-             TMath::Abs(firstPdgCode) == 421) // neutral
+             TMath::Abs(firstPdgCode) == 421 || TMath::Abs(pdgMother) == 311) // neutral
             && !(lastPdgCode == 22 || lastPdgCode == 111 || lastPdgCode == 130 || 
-                 TMath::Abs(lastPdgCode) == 2122 || lastPdgCode == 310 || 
+                 TMath::Abs(lastPdgCode) == 2112 || lastPdgCode == 310 || 
                  lastPdgCode == 221 || TMath::Abs(lastPdgCode) == 3122 || 
                  TMath::Abs(lastPdgCode) == 3322 || lastPdgCode == -3212 || 
-                 TMath::Abs(lastPdgCode) == 421)) // not neutral
+                 TMath::Abs(lastPdgCode) == 421 || TMath::Abs(pdgMother) == 311)) // not neutral
            || !((firstPdgCode == 22 || firstPdgCode == 111 || firstPdgCode == 130 || 
-                 TMath::Abs(firstPdgCode) == 2122 || firstPdgCode == 310 || 
+                 TMath::Abs(firstPdgCode) == 2112 || firstPdgCode == 310 || 
                  firstPdgCode == 221 || TMath::Abs(firstPdgCode) == 3122 || 
                  TMath::Abs(firstPdgCode) == 3322 || firstPdgCode == -3212 || 
-                 TMath::Abs(firstPdgCode) == 421) // not neutral
+                 TMath::Abs(firstPdgCode) == 421 || TMath::Abs(pdgMother) == 311) // not neutral
                 && (lastPdgCode == 22 || lastPdgCode == 111 || lastPdgCode == 130 || 
-                    TMath::Abs(lastPdgCode) == 2122 || lastPdgCode == 310 || 
+                    TMath::Abs(lastPdgCode) == 2112 || lastPdgCode == 310 || 
                     lastPdgCode == 221 || TMath::Abs(lastPdgCode) == 3122 || 
                     TMath::Abs(lastPdgCode) == 3322 || lastPdgCode == -3212 || 
-                    TMath::Abs(lastPdgCode) == 421)))) { // neutral
+                    TMath::Abs(lastPdgCode) == 421 || TMath::Abs(pdgMother) == 311)))) { // neutral
       
       vertex->SetType(AliAODVertex::kKink);
       jKinks++;
@@ -534,20 +595,20 @@ void SetVertexType(TParticle *part, AliAODVertex *vertex) {
     Int_t lastPdgCode = stack->Particle(mother->GetLastDaughter())->GetPdgCode();
 
     if ((pdgMother == 22 || pdgMother == 111 || pdgMother == 130 || 
-        TMath::Abs(pdgMother) == 2122 || pdgMother == 310 || 
+        TMath::Abs(pdgMother) == 2112 || pdgMother == 310 || 
         pdgMother == 221 || TMath::Abs(pdgMother) == 3122 || 
         TMath::Abs(pdgMother) == 3322 || pdgMother == -3212 || 
-        TMath::Abs(pdgMother) == 421) // neutral
+        TMath::Abs(pdgMother) == 421 || TMath::Abs(pdgMother) == 311) // neutral
        && !(lastPdgCode == 22 || lastPdgCode == 111 || lastPdgCode == 130 || 
-            TMath::Abs(lastPdgCode) == 2122 || lastPdgCode == 310 || 
+            TMath::Abs(lastPdgCode) == 2112 || lastPdgCode == 310 || 
             lastPdgCode == 221 || TMath::Abs(lastPdgCode) == 3122 || 
             TMath::Abs(lastPdgCode) == 3322 || lastPdgCode == -3212 || 
-            TMath::Abs(lastPdgCode) == 421) // not neutral
+            TMath::Abs(lastPdgCode) == 421 || TMath::Abs(pdgMother) == 311) // not neutral
        && !(firstPdgCode == 22 || firstPdgCode == 111 || firstPdgCode == 130 || 
-            TMath::Abs(firstPdgCode) == 2122 || firstPdgCode == 310 || 
+            TMath::Abs(firstPdgCode) == 2112 || firstPdgCode == 310 || 
             firstPdgCode == 221 || TMath::Abs(firstPdgCode) == 3122 || 
             TMath::Abs(firstPdgCode) == 3322 || firstPdgCode == -3212 || 
-            TMath::Abs(firstPdgCode) == 421)) { // not neutral
+            TMath::Abs(firstPdgCode) == 421 || TMath::Abs(pdgMother) == 311)) { // not neutral
       
       vertex->SetType(AliAODVertex::kV0);
       jV0s++;
@@ -559,19 +620,19 @@ void SetVertexType(TParticle *part, AliAODVertex *vertex) {
     Int_t firstPdgCode = stack->Particle(mother->GetFirstDaughter())->GetPdgCode();
     Int_t lastPdgCode = stack->Particle(mother->GetLastDaughter())->GetPdgCode();
     
-    if ((TMath::Abs(pdgMother) == 3332 || TMath::Abs(pdgMother) == 3312 || TMath::Abs(pdgMother) == 3322) &&
+    if ((TMath::Abs(pdgMother) == 3334 || TMath::Abs(pdgMother) == 3312 || TMath::Abs(pdgMother) == 3322) &&
        (TMath::Abs(pdgPart) == 3122 || TMath::Abs(pdgPart) == 211 || TMath::Abs(pdgPart) == 321)
        && ((!(firstPdgCode == 22 || firstPdgCode == 111 || firstPdgCode == 130 || 
-              TMath::Abs(firstPdgCode) == 2122 || firstPdgCode == 310 || 
+              TMath::Abs(firstPdgCode) == 2112 || firstPdgCode == 310 || 
               firstPdgCode == 221 || TMath::Abs(firstPdgCode) == 3122 || 
               TMath::Abs(firstPdgCode) == 3322 || firstPdgCode == -3212 || 
-              TMath::Abs(firstPdgCode) == 421) // not neutral   
+              TMath::Abs(firstPdgCode) == 421 || TMath::Abs(pdgMother) == 311) // not neutral   
             && TMath::Abs(lastPdgCode) == 3122) // labmda or anti-lambda
            || ((!(lastPdgCode == 22 || lastPdgCode == 111 || lastPdgCode == 130 || 
-                  TMath::Abs(lastPdgCode) == 2122 || lastPdgCode == 310 || 
+                  TMath::Abs(lastPdgCode) == 2112 || lastPdgCode == 310 || 
                   lastPdgCode == 221 || TMath::Abs(lastPdgCode) == 3122 || 
                   TMath::Abs(lastPdgCode) == 3322 || lastPdgCode == -3212 || 
-                  TMath::Abs(lastPdgCode) == 421) // not neutral
+                  TMath::Abs(lastPdgCode) == 421 || TMath::Abs(pdgMother) == 311) // not neutral
                 && TMath::Abs(firstPdgCode) == 3122)))) { // lambda or anti-lambda
       vertex->SetType(AliAODVertex::kCascade);
       jCascades++;