]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
skipping injected signals in MC
authorjgrosseo <jgrosseo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 23 Jul 2012 16:06:06 +0000 (16:06 +0000)
committerjgrosseo <jgrosseo@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 23 Jul 2012 16:06:06 +0000 (16:06 +0000)
PWGCF/Correlations/Base/AliAnalyseLeadingTrackUE.cxx
PWGCF/Correlations/Base/AliAnalyseLeadingTrackUE.h
PWGCF/Correlations/DPhi/AliAnalysisTaskPhiCorrelations.cxx
PWGCF/Correlations/DPhi/AliAnalysisTaskPhiCorrelations.h

index c75b607990aaf7bc2d25ea5aec09ebba72ad292d..34a8d3654dbd7ca7c32eb8f3ac6aed8664e28485 100644 (file)
@@ -199,6 +199,88 @@ TObjArray*  AliAnalyseLeadingTrackUE::FindLeadingObjects(TObject *obj)
   }
 
 
+void AliAnalyseLeadingTrackUE::RemoveInjectedSignals(TObjArray* tracks, TObject* mcObj, Int_t maxLabel)
+{
+  // remove injected signals (primaries above <maxLabel>)
+  // <tracks> can be the following cases:
+  // a. tracks: in this case the label is taken and then case b.
+  // b. particles: the first stable mother is searched and checked if it is <= <maxLabel>
+  // <mcObj> can be AOD (TClonesArray) or ESD (AliMCEvent)
+  
+  TClonesArray* arrayMC = 0;
+  AliMCEvent* mcEvent = 0;
+  if (mcObj->InheritsFrom("AliMCEvent"))
+    mcEvent = static_cast<AliMCEvent*>(mcObj);
+  else if (mcObj->InheritsFrom("TClonesArray"))
+    arrayMC = static_cast<TClonesArray*>(mcObj);
+  else
+  {
+    arrayMC->Dump();
+    AliFatal("Invalid object passed");
+  }
+  
+  Int_t before = tracks->GetEntriesFast();
+
+  for (Int_t i=0; i<before; ++i) 
+  {
+    AliVParticle* part = (AliVParticle*) tracks->At(i);
+    
+    if (part->InheritsFrom("AliESDtrack") || part->InheritsFrom("AliAODTrack"))
+      part = ((mcEvent) ? mcEvent->GetTrack(TMath::Abs(part->GetLabel())) : (AliVParticle*)arrayMC->At(TMath::Abs(part->GetLabel())));
+      
+    AliVParticle* mother = part;
+    if (mcEvent)
+    {
+      while (!mcEvent->IsPhysicalPrimary(mother->GetLabel()))
+      {
+       if (((AliMCParticle*)mother)->GetMother() < 0)
+       {
+         mother = 0;
+         break;
+       }
+
+       mother = (AliMCParticle*) mcEvent->GetTrack(((AliMCParticle*)mother)->GetMother());
+       if (!mother)
+         break;
+      }
+    }
+    else
+    {
+      // find the primary mother
+      while (!((AliAODMCParticle*)mother)->IsPhysicalPrimary())
+      {
+       if (((AliAODMCParticle*)mother)->GetMother() < 0)
+       {
+         mother = 0;
+         break;
+       }
+         
+       mother = (AliVParticle*) arrayMC->At(((AliAODMCParticle*)mother)->GetMother());
+       if (!mother)
+         break;
+      }
+    }
+    
+    if (!mother)
+    {
+      Printf("WARNING: No mother found for particle %d:", part->GetLabel());
+      continue;
+    }
+   
+    if (mother->GetLabel() > maxLabel)
+    {
+//       Printf("Removing %d with label %d", i, part->GetLabel()); part->Dump();
+      TObject* object = tracks->RemoveAt(i);
+      if (tracks->IsOwner())
+       delete object;
+    }
+  }
+  tracks->Compress();
+  
+  AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast())); 
+}
+
 //-------------------------------------------------------------------
 TObjArray* AliAnalyseLeadingTrackUE::GetAcceptedParticles(TObject* obj, TObject* arrayMC, Bool_t onlyprimaries, Int_t particleSpecies, Bool_t useEtaPtCuts)
 {
index 6b54c4d1d2372b4eafc51112b93e423e54feceff..55de4ffbe568034ca02212b8e0ff184b7954c603 100644 (file)
@@ -57,7 +57,8 @@ class AliAnalyseLeadingTrackUE : public TObject {
   TObjArray*     GetAcceptedParticles(TObject* obj, TObject* arrayMC, Bool_t onlyprimaries = kTRUE, Int_t particleSpecies = -1, Bool_t useEtaPtCuts = kFALSE); 
   Bool_t         TriggerSelection(const TObject* obj);                                   // Select good triggers with AliPhysicsSelection class
   Bool_t         VertexSelection(const TObject* obj, Int_t ntracks, Double_t zed);       // Vertex selection: see implementation
-
+  void                  RemoveInjectedSignals(TObjArray* tracks, TObject* arrayMC, Int_t maxLabel);
+  
  private:
   Int_t          fDebug;             // debug flag
   Int_t          fFilterBit;         // track selection cuts
index aab3d332178e29a9c982e6ae5e95ce014f9d19d5..b508b396bb4106b65f949010bf7b1cd08a3efd87 100644 (file)
@@ -45,6 +45,8 @@
 #include "AliCentrality.h"
 #include "AliStack.h"
 #include "AliAODMCHeader.h"
+#include "AliGenCocktailEventHeader.h"
+#include "AliGenEventHeader.h"
 
 #include "AliEventPoolManager.h"
 
@@ -87,6 +89,7 @@ fTwoTrackEfficiencyStudy(kFALSE),
 fTwoTrackEfficiencyCut(0),
 fUseVtxAxis(kFALSE),
 fSkipTrigger(kFALSE),
+fInjectedSignals(kFALSE),
 // pointers to UE classes
 fAnalyseUE(0x0),
 fHistos(0x0),
@@ -326,6 +329,7 @@ void  AliAnalysisTaskPhiCorrelations::AddSettingsTree()
   settingsTree->Branch("fkTrackingEfficiency", "TH1D", &fkTrackingEfficiency);
   settingsTree->Branch("fMixingTracks", &fMixingTracks,"MixingTracks/I");
   settingsTree->Branch("fSkipTrigger", &fSkipTrigger,"SkipTrigger/O");
+  settingsTree->Branch("fInjectedSignals", &fInjectedSignals,"SkipTrigger/O");
   settingsTree->Branch("fRejectCentralityOutliers", &fRejectCentralityOutliers,"RejectCentralityOutliers/O");
   
   settingsTree->Fill();
@@ -398,12 +402,57 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
   Float_t weight = 1;
   if (fFillpT)
     weight = -1;
+  
+  // For productions with injected signals, figure out above which label to skip particles/tracks
+  Int_t skipParticlesAbove = 0;
+  if (fInjectedSignals)
+  {
+    AliGenEventHeader* eventHeader = 0;
+    Int_t headers = 0;
+
+    if (fMcEvent)
+    {
+      // ESD
+      AliHeader* header = (AliHeader*) fMcEvent->Header();
+      if (!header)
+       AliFatal("fInjectedSignals set but no MC header found");
+       
+      AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
+      if (!cocktailHeader)
+      {
+       header->Dump();
+       AliFatal("fInjectedSignals set but no MC cocktail header found");
+      }
+
+      headers = cocktailHeader->GetHeaders()->GetEntries();
+      eventHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
+    }
+    else
+    {
+      // AOD
+      AliAODMCHeader* header = (AliAODMCHeader*) fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
+      if (!header)
+       AliFatal("fInjectedSignals set but no MC header found");
+      
+      headers = header->GetNCocktailHeaders();
+      eventHeader = header->GetCocktailHeader(0);
+    }
+    
+    if (!eventHeader)
+      AliFatal("First event header not found");
     
+    skipParticlesAbove = eventHeader->NProduced();
+    AliInfo(Form("Injected signals in this event (%d headers). Keeping events of %s. Will skip particles/tracks above %d.", headers, eventHeader->ClassName(), skipParticlesAbove));
+  }
+  
   // Get MC primaries
   TObjArray* tmpList = fAnalyseUE->GetAcceptedParticles(mc, 0, kTRUE, -1, kTRUE);
+  if (fInjectedSignals)
+    fAnalyseUE->RemoveInjectedSignals(tmpList, mc, skipParticlesAbove);
   TObjArray* tracksMC = CloneAndReduceTrackList(tmpList);
   delete tmpList;
   
+  /*
   if (fAOD)
   {
     for (Int_t i=0; i<fArrayMC->GetEntriesFast(); i++)
@@ -414,6 +463,7 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
     for (Int_t i=0; i<fMcEvent->GetNumberOfTracks(); i++)
       ((TH1F*) fListOfHistos->FindObject("pids"))->Fill(fMcEvent->GetTrack(i)->PdgCode());
   }
+  */
   
   if (fFillOnlyStep0)
     zVtx = 0;
@@ -461,6 +511,13 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
         TObjArray* primMCParticles = fAnalyseUE->GetAcceptedParticles(mc, 0x0, kTRUE, particleSpecies, kTRUE);
         TObjArray* primRecoTracksMatched = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, particleSpecies, kTRUE);
         TObjArray* allRecoTracksMatched = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, particleSpecies, kTRUE);
+       
+       if (fInjectedSignals)
+       {
+         fAnalyseUE->RemoveInjectedSignals(primMCParticles, mc, skipParticlesAbove);
+         fAnalyseUE->RemoveInjectedSignals(primRecoTracksMatched, mc, skipParticlesAbove);
+         fAnalyseUE->RemoveInjectedSignals(allRecoTracksMatched, mc, skipParticlesAbove);
+       }
       
         fHistos->FillTrackingEfficiency(primMCParticles, primRecoTracksMatched, allRecoTracksMatched, particleSpecies, centrality);
         
@@ -480,6 +537,8 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
       
       // Get MC primaries that match reconstructed track
       TObjArray* tracksRecoMatchedPrim = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, -1, kTRUE);
+      if (fInjectedSignals)
+       fAnalyseUE->RemoveInjectedSignals(tracksRecoMatchedPrim, mc, skipParticlesAbove);
       
       // (RECO-matched (quantities from MC particle) primary particles)
       // STEP 4
@@ -487,6 +546,8 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
       
       // Get MC primaries + secondaries that match reconstructed track
       TObjArray* tracksRecoMatchedAll = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, -1, kTRUE);
+      if (fInjectedSignals)
+       fAnalyseUE->RemoveInjectedSignals(tracksRecoMatchedAll, mc, skipParticlesAbove);
       
       // (RECO-matched (quantities from MC particle) all particles)
       // STEP 5
@@ -494,7 +555,9 @@ void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
       
       // Get RECO tracks
       TObjArray* tracks = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, -1, kTRUE);
-      
+      if (fInjectedSignals)
+       fAnalyseUE->RemoveInjectedSignals(tracks, mc, skipParticlesAbove);
+     
       // (RECO all tracks)
       // STEP 6
       if (!fSkipStep6)
index 6d36aa1234f33b121f5a23886d9d4bd9ca509b14..a03b78a41c1883af624524e13b9d1f1003a7c452 100644 (file)
@@ -68,6 +68,7 @@ class  AliAnalysisTaskPhiCorrelations : public AliAnalysisTask
     virtual    void    SetTwoTrackEfficiencyCut(Float_t value = 0.02) { fTwoTrackEfficiencyCut = value; }
     virtual    void    SetUseVtxAxis(Bool_t flag) { fUseVtxAxis = flag; }
     virtual     void    SetSkipTrigger(Bool_t flag) { fSkipTrigger = flag; }
+    virtual     void    SetInjectedSignals(Bool_t flag) { fInjectedSignals = flag; }
     
     // histogram settings
     void SetTrackingEfficiency( const TH1D* hist) { fkTrackingEfficiency = hist; }
@@ -116,6 +117,7 @@ class  AliAnalysisTaskPhiCorrelations : public AliAnalysisTask
     Float_t            fTwoTrackEfficiencyCut;   // enable two-track efficiency cut
     Bool_t             fUseVtxAxis;              // use z vtx as axis (needs 7 times more memory!)
     Bool_t             fSkipTrigger;             // skip trigger selection
+    Bool_t             fInjectedSignals;         // check header to skip injected signals in MC
     
     // Pointers to external UE classes
     AliAnalyseLeadingTrackUE*     fAnalyseUE;      //! points to class containing common analysis algorithms
@@ -161,7 +163,7 @@ class  AliAnalysisTaskPhiCorrelations : public AliAnalysisTask
     
     Bool_t fFillpT;                // fill sum pT instead of number density
     
-    ClassDef( AliAnalysisTaskPhiCorrelations, 12); // Analysis task for delta phi correlations
+    ClassDef( AliAnalysisTaskPhiCorrelations, 13); // Analysis task for delta phi correlations
   };
 
 class AliDPhiBasicParticle : public AliVParticle