extended task to read also AOD MC information, added check for kGridDataSet for runni...
authoresicking <esicking@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 25 Oct 2010 09:37:13 +0000 (09:37 +0000)
committeresicking <esicking@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 25 Oct 2010 09:37:13 +0000 (09:37 +0000)
PWG4/JetTasks/AliAnalysisTaskMinijet.cxx
PWG4/macros/AddTaskMinijet.C

index eb92716..be08db2 100644 (file)
@@ -18,6 +18,7 @@
 
 #include "AliAODEvent.h"
 #include "AliAODTrack.h"
+#include "AliAODMCParticle.h"
 
 #include "AliStack.h"
 #include "AliMCEvent.h"
@@ -273,8 +274,6 @@ void AliAnalysisTaskMinijet::UserExec(Option_t *)
     //AOD MCreading and analysis
     //------
     if (fUseMC){
-      Printf("Not yet implemented");
-      return;
       ntracks = LoopAODMC(&pt, &eta, &phi);//read tracks
       if(ntracks>0){
        Analyse(pt, eta, phi, ntracks, 0, 3);//analyse
@@ -523,43 +522,70 @@ Int_t AliAnalysisTaskMinijet::LoopAOD(Float_t **ptArray, Float_t ** etaArray,
 Int_t AliAnalysisTaskMinijet::LoopAODMC(Float_t **ptArray, Float_t ** etaArray, 
                                        Float_t ** phiArray)
 {
-  // Main loop
-  // Called for each event
+  // gives back the number of AOD MC particles and pointer to arrays with particle 
+  // properties (pt, eta, phi)
   
-  AliWarning("Not implemented yet !!!")
-  return 0;
-
-  AliMCEvent *mcEvent = (AliMCEvent*) MCEvent();
-  if (!mcEvent) {
-    Error("UserExec", "Could not retrieve MC event");
-    return 0;
+  //retreive MC particles from event
+  TClonesArray *mcArray = (TClonesArray*)fAODEvent->
+    FindListObject(AliAODMCParticle::StdBranchName());
+  if(!mcArray){
+    Printf("%s:%d No MC particle branch found",(char*)__FILE__,__LINE__);
+    return kFALSE;
   }
   
-  Int_t ntracks = mcEvent->GetNumberOfTracks();
+  Int_t ntracks = mcArray->GetEntriesFast();
   if(fDebug)Printf("MC particles: %d", ntracks);
 
-  *ptArray = new Float_t[ntracks]; 
-  *etaArray = new Float_t[ntracks]; 
-  *phiArray = new Float_t[ntracks]; 
-  
-  // Track loop
-  for (Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
-    AliMCParticle *track = dynamic_cast<AliMCParticle*>(mcEvent->GetTrack(iTracks));
+
+  // Track loop: chek how many particles will be accepted
+  Float_t vzMC=0.;
+  Int_t nAcceptedTracks=0;
+  for (Int_t it = 0; it < ntracks; it++) {
+    AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
     if (!track) {
-      Error("UserExec", "Could not receive track %d", iTracks);
+      Error("UserExec", "Could not receive track %d", it);
       continue;
     }
+    if(!track->IsPhysicalPrimary())continue;
+    if (track->Charge() == 0) continue;
+    if(TMath::Abs(track->Eta())>1. || track->Pt()<0.2 || track->Pt()>200) continue; //same cuts as in ESD filter
+    if(track->TestBit(16))nAcceptedTracks++;
+    if(nAcceptedTracks==1) vzMC= track->Zv(); // check only one time. (only one vertex per event allowed)
+  }
 
-    fHistPtMC->Fill(track->Pt());
+  //generate array with size of number of accepted tracks
+  *ptArray = new Float_t[nAcceptedTracks]; 
+  *etaArray = new Float_t[nAcceptedTracks]; 
+  *phiArray = new Float_t[nAcceptedTracks]; 
+  
+  if(nAcceptedTracks==0) return 0;
 
-    //fills arrays with track properties
-    (*ptArray)[iTracks]  = track->Pt();
-    (*etaArray)[iTracks] = track->Eta();
-    (*phiArray)[iTracks] = track->Phi();
+  //check vertex
+  if(TMath::Abs(vzMC)>fVertexZCut) return 0;
+  fVertexZ[3]->Fill(vzMC);
+  
 
-  } //track loop
+  // Track loop: fill arrays for accepted tracks 
+  Int_t iAcceptedTracks=0;
+  for (Int_t it = 0; it < ntracks; it++) {
+    AliAODMCParticle *track = (AliAODMCParticle*)mcArray->At(it);
+    if (!track) {
+      Error("UserExec", "Could not receive track %d", it);
+      continue;
+    }
+    if(!track->IsPhysicalPrimary())continue;
+    if (track->Charge() == 0) continue;
+    if(TMath::Abs(track->Eta())>0.9 || track->Pt()<0.2 || track->Pt()>200) continue;
+
+    if(track->TestBit(16)){
+      fHistPtMC->Fill(track->Pt());
+      (*ptArray)[iAcceptedTracks]  = track->Pt();
+      (*etaArray)[iAcceptedTracks] = track->Eta();
+      (*phiArray)[iAcceptedTracks++] = track->Phi();
+    }
+  }
   
-  return ntracks;
+  return nAcceptedTracks;
 
 } 
 
index 125324d..46170d0 100644 (file)
@@ -1,12 +1,24 @@
-AliAnalysisTaskMinijet* AddTaskMinijet(TString format="esd",Bool_t useMC = kFALSE, Bool_t IsHighMult=kTRUE)
+AliAnalysisTaskMinijet* AddTaskMinijet(TString format="esd",Bool_t useMC = kFALSE, TString kGridDataSet="LHC10e")
 {
+
+  //starting with periode LHC10e, there are also events triggered with High Mult trigger
+  //add stand alone task in case these events are availale
+  Bool_t IsHighMult= (!kGridDataSet.CompareTo("LHC10e") || !kGridDataSet.CompareTo("LHC10f"));
+  //extent with new periods (LHC10g,h,...)
+
   // Get the pointer to the existing analysis manager via the static access method.
   //==============================================================================
   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
   if (!mgr) mgr = new AliAnalysisManager("Analysis train");
-  
-  
-  // Configure analysis
+
+  // Set cuts (used for ESDs) 
+  //===========================================================================
+  AliESDtrackCuts* esdTrackCutsITSTPC=0x0;
+  if(!format.CompareTo("esd")){
+    esdTrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
+  }
+
+  // Configure tasks
   //===========================================================================
    
   //first task for min bias events
@@ -18,10 +30,15 @@ AliAnalysisTaskMinijet* AddTaskMinijet(TString format="esd",Bool_t useMC = kFALS
   task->SetEventAxis(1); //1=random track
   task->SetDebugLevel(0);
   task->SetMaxVertexZ(10.);
+  if(!format.CompareTo("esd")){
+    task->SetCuts(esdTrackCutsITSTPC);
+    task->SetMode(0);//0=reading ESDs
+    task->SelectCollisionCandidates(AliVEvent::kMB);//MB
+  }
 
   //second task for high multipliciy events
   AliAnalysisTaskMinijet *taskHM =0x0;
-  if(IsHighMult){
+  if(IsHighMult && !format.CompareTo("esd")){
     taskHM  = new AliAnalysisTaskMinijet("AliAnalysisTaskMinijet HighMult");
     taskHM->UseMC(useMC);
     taskHM->SetRadiusCut(0.7);
@@ -30,28 +47,13 @@ AliAnalysisTaskMinijet* AddTaskMinijet(TString format="esd",Bool_t useMC = kFALS
     taskHM->SetEventAxis(1); //1=random track
     taskHM->SetDebugLevel(0);
     taskHM->SetMaxVertexZ(10.);
+    taskHM->SetCuts(esdTrackCutsITSTPC);
+    taskHM->SetMode(0);//0=reading ESDs
+    taskHM->SelectCollisionCandidates(AliVEvent::kHighMult);//high mult triggered event
   }
 
-  //set cuts (used for ESDs) 
-  AliESDtrackCuts* esdTrackCutsITSTPC=0x0;
-
-  if(!format.CompareTo("esd")){
-    esdTrackCutsITSTPC = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
-    
-    //configure tasks
-    task->SetMode(0);//0=reading ESDs
-    task->SelectCollisionCandidates(AliVEvent::kMB);//MB
-    task->SetCuts(esdTrackCutsITSTPC);
-    if(IsHighMult){
-      taskHM->SetMode(0);//0=reading ESDs
-      taskHM->SelectCollisionCandidates(AliVEvent::kHighMult);//high mult triggered event
-      taskHM->SetCuts(esdTrackCutsITSTPC);
-    }
-  }
-
-  else if (!format.CompareTo("aod")){
+  if (!format.CompareTo("aod")){
     task->SetMode(1);// 1 = reading AODs
-    task->UseMC(kFALSE);//MC AODs are not implemented yet
   }
 
   //create output container and add task to mgr