Extend MC truth information with z-Vertex and pileup
authormfasel <mfasel@lbl.gov>
Mon, 29 Sep 2014 11:26:37 +0000 (13:26 +0200)
committermverweij <marta.verweij@cern.ch>
Mon, 29 Sep 2014 11:39:29 +0000 (13:39 +0200)
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskPtEMCalTrigger.cxx
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskPtEMCalTrigger.h

index d89db43..e2b3e7c 100644 (file)
@@ -180,19 +180,25 @@ namespace EMCalTriggerPtAnalysis {
                TAxis hclusteraxes[4];
                DefineAxis(hclusteraxes[0], "energy", "E (GeV)", ptbinning);
                DefineAxis(hclusteraxes[1], "zvertex", "z_{V} (cm)", zvertexBinning);
-               DefineAxis(hclusteraxes[2], "pileup", "Pileup rejection", 2, -0.5, 1.5);
+               DefineAxis(hclusteraxes[2], "pileup", "Pileup rejection", 2,1 -0.5, 1.5);
                DefineAxis(hclusteraxes[3], "mbtrigger", "Has MB trigger", 2, -0.5, 1.5);
                const TAxis *clusteraxes[4];
                for(int iaxis = 0; iaxis < 4; ++iaxis) clusteraxes[iaxis] = hclusteraxes + iaxis;
+               TAxis hpatchaxes[3];
+               DefineAxis(hpatchaxes[0], "energy", "Patch energy (GeV)", 100, 0., 100.);
+               DefineAxis(hpatchaxes[1], "eta", "#eta", etabinning);
+               DefineAxis(hpatchaxes[2], "phi", "#phi",  20, 0, 2 * TMath::Pi());
+               const TAxis *patchaxes[3];
+               for(int iaxis = 0; iaxis < 3; ++iaxis) patchaxes[iaxis] = hpatchaxes + iaxis;
                try{
                        std::string patchnames[] = {"Level0", "JetHigh", "JetLow", "GammaHigh", "GammaLow"};
                        for(std::string * triggerpatch = patchnames; triggerpatch < patchnames + sizeof(patchnames)/sizeof(std::string); ++triggerpatch){
-                               fHistos->CreateTH1(Form("Energy%s", triggerpatch->c_str()), Form("Patch energy for %s trigger patches", triggerpatch->c_str()), 100, 0., 100.);
-                               fHistos->CreateTH1(Form("EnergyMain%s", triggerpatch->c_str()), Form("Patch energy for main %s trigger patches", triggerpatch->c_str()), 100, 0., 100.);
+                               fHistos->CreateTHnSparse(Form("Energy%s", triggerpatch->c_str()), Form("Patch energy for %s trigger patches", triggerpatch->c_str()), 3, patchaxes);
+                               fHistos->CreateTHnSparse(Form("EnergyMain%s", triggerpatch->c_str()), Form("Patch energy for main %s trigger patches", triggerpatch->c_str()), 3, patchaxes);
                        }
 
                        // Create histogram for MC-truth
-                       fHistos->CreateTHnSparse("hMCtrueParticles", "Particle-based histogram for MC-true particles", 3, trackaxes);
+                       fHistos->CreateTHnSparse("hMCtrueParticles", "Particle-based histogram for MC-true particles", 5, trackaxes);
                        for(std::map<std::string,std::string>::iterator it = triggerCombinations.begin(); it != triggerCombinations.end(); ++it){
                                const std::string name = it->first, &title = it->second;
                                // Create event-based histogram
@@ -238,48 +244,40 @@ namespace EMCalTriggerPtAnalysis {
                // Common checks: Have SPD vertex and primary vertex from tracks, and both need to have at least one contributor
                AliDebug(1,Form("Number of calibrated clusters: %d", fCaloClusters->GetEntries()));
                AliDebug(1,Form("Number of matched tracks: %d", fTracks->GetEntries()));
-
                if(fMCEvent){
-                       for(int ipart = 0; ipart < fMCEvent->GetNumberOfTracks(); ipart++){
-                               // Select only physical primary particles
-                               AliVParticle *part = fMCEvent->GetTrack(ipart);
-                               if(!fEtaRange.IsInRange(part->Eta())) continue;
-                               if(!fPtRange.IsInRange(part->Pt())) continue;
-                               if(!fMCEvent->IsPhysicalPrimary(ipart)) continue;
-                               FillMCParticleHist(part);
-                       }
                        // Build always trigger strig from the trigger maker in case of MC
-               fUseTriggersFromTriggerMaker = kTRUE;
+                       fUseTriggersFromTriggerMaker = kTRUE;
                }
 
                // Loop over trigger patches, fill patch energy
                AliEmcalTriggerPatchInfo *triggerpatch(NULL);
                TIter patchIter(this->fTriggerPatchInfo);
                while((triggerpatch = dynamic_cast<AliEmcalTriggerPatchInfo *>(patchIter()))){
+                       double triggerpatchinfo[3] = {triggerpatch->GetPatchE(), triggerpatch->GetEtaCM(), triggerpatch->GetPhiCM()};
                        if(triggerpatch->IsJetHigh()){
-                               fHistos->FillTH1("EnergyJetHigh", triggerpatch->GetPatchE());
+                               fHistos->FillTHnSparse("EnergyJetHigh", triggerpatchinfo);
                                if(triggerpatch->IsMainTrigger())
-                                       fHistos->FillTH1("EnergyMainJetHigh", triggerpatch->GetPatchE());
+                                       fHistos->FillTHnSparse("EnergyMainJetHigh", triggerpatchinfo);
                        }
                        if(triggerpatch->IsJetLow()){
-                               fHistos->FillTH1("EnergyJetLow", triggerpatch->GetPatchE());
+                               fHistos->FillTHnSparse("EnergyJetLow", triggerpatchinfo);
                                if(triggerpatch->IsMainTrigger())
-                                       fHistos->FillTH1("EnergyMainJetLow", triggerpatch->GetPatchE());
+                                       fHistos->FillTHnSparse("EnergyMainJetLow", triggerpatchinfo);
                        }
                        if(triggerpatch->IsGammaHigh()){
-                               fHistos->FillTH1("EnergyGammaHigh", triggerpatch->GetPatchE());
+                               fHistos->FillTHnSparse("EnergyGammaHigh", triggerpatchinfo);
                                if(triggerpatch->IsMainTrigger())
-                                       fHistos->FillTH1("EnergyMainGammaHigh", triggerpatch->GetPatchE());
+                                       fHistos->FillTHnSparse("EnergyMainGammaHigh", triggerpatchinfo);
                        }
                        if(triggerpatch->IsGammaLow()){
-                               fHistos->FillTH1("EnergyGammaLow", triggerpatch->GetPatchE());
+                               fHistos->FillTHnSparse("EnergyGammaLow", triggerpatchinfo);
                                if(triggerpatch->IsMainTrigger())
-                                       fHistos->FillTH1("EnergyMainGammaLow", triggerpatch->GetPatchE());
+                                       fHistos->FillTHnSparse("EnergyMainGammaLow", triggerpatchinfo);
                        }
                        if(triggerpatch->IsLevel0()){
-                               fHistos->FillTH1("EnergyLevel0", triggerpatch->GetPatchE());
+                               fHistos->FillTHnSparse("EnergyLevel0", triggerpatchinfo);
                                if(triggerpatch->IsMainTrigger())
-                                       fHistos->FillTH1("EnergyMainLevel0", triggerpatch->GetPatchE());
+                                       fHistos->FillTHnSparse("EnergyMainLevel0", triggerpatchinfo);
                        }
                }
 
@@ -371,6 +369,18 @@ namespace EMCalTriggerPtAnalysis {
                                FillEventHist(it->c_str(), zv, isPileupEvent);
                }
 
+               // Fill MC truth
+               if(fMCEvent){
+                       for(int ipart = 0; ipart < fMCEvent->GetNumberOfTracks(); ipart++){
+                               // Select only physical primary particles
+                               AliVParticle *part = fMCEvent->GetTrack(ipart);
+                               if(!fEtaRange.IsInRange(part->Eta())) continue;
+                               if(!fPtRange.IsInRange(part->Pt())) continue;
+                               if(!fMCEvent->IsPhysicalPrimary(ipart)) continue;
+                               FillMCParticleHist(part, zv, isPileupEvent);
+                       }
+               }
+
                AliVTrack *track(NULL);
                // Loop over all tracks (No cuts applied)
                TIter allTrackIter(fTracks);
@@ -699,14 +709,18 @@ namespace EMCalTriggerPtAnalysis {
        }
 
        //______________________________________________________________________________
-       void AliAnalysisTaskPtEMCalTrigger::FillMCParticleHist(const AliVParticle * const track){
+       void AliAnalysisTaskPtEMCalTrigger::FillMCParticleHist(const AliVParticle * const track, double vz, bool isPileup){
                /*
                 * Fill histogram for MC-true particles with the information pt, eta and phi
                 *
                 * @param track: the Monte-Carlo track
                 */
-               double data[3] = {TMath::Abs(track->Pt()), track->Eta(), track->Phi()};
+               double data[5] = {TMath::Abs(track->Pt()), track->Eta(), track->Phi(), vz, 0.};
                fHistos->FillTHnSparse("hMCtrueParticles", data);
+               if(!isPileup){
+                       data[4] = 1.;
+                       fHistos->FillTHnSparse("hMCtrueParticles", data);
+               }
        }
 
        //______________________________________________________________________________
index 434d6cb..f5afc87 100644 (file)
@@ -53,7 +53,7 @@ private:
        void FillEventHist(const char *trigger, double vz, bool isPileup);
        void FillTrackHist(const char *trigger, const AliVTrack *track, double vz, bool isPileup, int cut, bool isMinBias);
        void FillClusterHist(const char *trigger, const AliVCluster *clust, bool isCalibrated, double vz, bool isPileup, bool isMinBias);
-       void FillMCParticleHist(const AliVParticle * const part);
+       void FillMCParticleHist(const AliVParticle * const part, double vz, bool isPileup);
        bool IsTrueTrack(const AliVTrack *const) const;
        TString BuildTriggerString();
        const AliVVertex *GetSPDVertex() const;