]>
Commit | Line | Data |
---|---|---|
46f589c2 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | /* | |
17 | * Analysis task of the pt analysis on EMCal-triggered events | |
18 | * | |
19 | * Author: Markus Fasel | |
20 | */ | |
21 | ||
22 | #include <map> | |
23 | #include <cstring> | |
24 | #include <iostream> | |
25 | #include <memory> | |
26 | #include <vector> | |
27 | #include <string> | |
28 | #include <sstream> | |
29 | ||
bf1cb6ad | 30 | #include <TClonesArray.h> |
46f589c2 | 31 | #include <TDirectory.h> |
32 | #include <TH1.h> | |
33 | #include <THashList.h> | |
eed93770 | 34 | #include <TList.h> |
46f589c2 | 35 | #include <TKey.h> |
36 | #include <TMath.h> | |
37 | #include <TObjArray.h> | |
555bdd7f | 38 | #include <TObjString.h> |
46f589c2 | 39 | #include <TString.h> |
555bdd7f | 40 | #include <TVector2.h> |
46f589c2 | 41 | |
cdc26d91 | 42 | #include "AliAODEvent.h" |
46f589c2 | 43 | #include "AliESDEvent.h" |
cdc26d91 | 44 | #include "AliInputEventHandler.h" |
45 | #include "AliMCEvent.h" | |
46 | #include "AliParticleContainer.h" | |
47 | #include "AliVCluster.h" | |
48 | #include "AliVParticle.h" | |
49 | #include "AliVTrack.h" | |
50 | #include "AliVVertex.h" | |
46f589c2 | 51 | |
555bdd7f | 52 | #include "AliClusterContainer.h" |
53 | #include "AliEmcalJet.h" | |
64321abd | 54 | #include "AliEmcalPhysicsSelection.h" |
3fdf176c | 55 | #include "AliEmcalTriggerPatchInfo.h" |
46f589c2 | 56 | #include "AliEMCalHistoContainer.h" |
cdc26d91 | 57 | #include "AliEMCalPtTaskVTrackSelection.h" |
009a14cb | 58 | #include "AliEMCalPtTaskTrackSelectionAOD.h" |
cdc26d91 | 59 | #include "AliEMCalPtTaskTrackSelectionESD.h" |
555bdd7f | 60 | #include "AliJetContainer.h" |
61 | #include "AliParticleContainer.h" | |
7d1a1e9a | 62 | #include "AliPicoTrack.h" |
46f589c2 | 63 | #include "AliAnalysisTaskPtEMCalTrigger.h" |
64 | ||
65 | ClassImp(EMCalTriggerPtAnalysis::AliAnalysisTaskPtEMCalTrigger) | |
66 | ||
67 | namespace EMCalTriggerPtAnalysis { | |
68 | ||
555bdd7f | 69 | /* |
70 | * constants | |
71 | */ | |
72 | const Int_t AliAnalysisTaskPtEMCalTrigger::kNJetRadii = 4; | |
73 | const Double_t jetRadVals[AliAnalysisTaskPtEMCalTrigger::kNJetRadii] = {0.2, 0.3, 0.4, 0.5}; | |
74 | const Double_t *AliAnalysisTaskPtEMCalTrigger::kJetRadii = jetRadVals; | |
75 | ||
9f248f6f | 76 | //______________________________________________________________________________ |
77 | AliAnalysisTaskPtEMCalTrigger::AliAnalysisTaskPtEMCalTrigger(): | |
555bdd7f | 78 | AliAnalysisTaskEmcalJet(), |
9f248f6f | 79 | fHistos(NULL), |
80 | fListTrackCuts(NULL), | |
81 | fEtaRange(), | |
82 | fPtRange(), | |
306d380b | 83 | fEnergyRange(), |
84 | fVertexRange(), | |
555bdd7f | 85 | fJetContainersMC(), |
86 | fJetContainersData(), | |
9f248f6f | 87 | fSwapEta(kFALSE), |
88 | fUseTriggersFromTriggerMaker(kFALSE) | |
89 | { | |
90 | /* | |
91 | * Dummy constructor, initialising the values with default (NULL) values | |
92 | */ | |
93 | } | |
94 | ||
95 | //______________________________________________________________________________ | |
96 | AliAnalysisTaskPtEMCalTrigger::AliAnalysisTaskPtEMCalTrigger(const char *name): | |
555bdd7f | 97 | AliAnalysisTaskEmcalJet(name, kTRUE), |
9f248f6f | 98 | fHistos(NULL), |
99 | fListTrackCuts(NULL), | |
100 | fEtaRange(), | |
101 | fPtRange(), | |
306d380b | 102 | fEnergyRange(), |
103 | fVertexRange(), | |
555bdd7f | 104 | fJetContainersMC(), |
105 | fJetContainersData(), | |
9f248f6f | 106 | fSwapEta(kFALSE), |
107 | fUseTriggersFromTriggerMaker(kFALSE) | |
108 | { | |
109 | /* | |
110 | * Main constructor, setting default values for eta and zvertex cut | |
111 | */ | |
112 | ||
113 | fListTrackCuts = new TList; | |
114 | fListTrackCuts->SetOwner(false); | |
115 | ||
116 | // Set default cuts | |
117 | fEtaRange.SetLimits(-0.8, 0.8); | |
118 | fPtRange.SetLimits(0.15, 100.); | |
caff8cd1 | 119 | fEnergyRange.SetLimits(0., 1000.); |
2ae55f0a | 120 | fVertexRange.SetLimits(-10., 10.); |
9f248f6f | 121 | SetMakeGeneralHistograms(kTRUE); |
122 | } | |
123 | ||
124 | //______________________________________________________________________________ | |
125 | AliAnalysisTaskPtEMCalTrigger::~AliAnalysisTaskPtEMCalTrigger(){ | |
126 | /* | |
127 | * Destructor, deleting output | |
128 | */ | |
129 | //if(fTrackSelection) delete fTrackSelection; | |
130 | if(fHistos) delete fHistos; | |
131 | if(fListTrackCuts) delete fListTrackCuts; | |
132 | } | |
133 | ||
134 | //______________________________________________________________________________ | |
135 | void AliAnalysisTaskPtEMCalTrigger::UserCreateOutputObjects(){ | |
136 | /* | |
137 | * Create the list of output objects and define the histograms. | |
138 | * Also adding the track cuts to the list of histograms. | |
139 | */ | |
140 | AliAnalysisTaskEmcal::UserCreateOutputObjects(); | |
e99ef86a | 141 | SetCaloTriggerPatchInfoName("EmcalTriggers"); |
9f248f6f | 142 | fHistos = new AliEMCalHistoContainer("PtEMCalTriggerHistograms"); |
143 | fHistos->ReleaseOwner(); | |
144 | ||
e99ef86a | 145 | if(fJetContainersMC.GetEntries()){ |
146 | AliDebug(1,"Jet containers for MC truth:"); | |
147 | TObjString *contname(NULL); | |
148 | TIter contMCIter(&fJetContainersMC); | |
149 | while((contname = dynamic_cast<TObjString *>(contMCIter.Next()))) | |
150 | AliDebug(1, Form("Next container: %s", contname->String().Data())); | |
151 | } | |
152 | if(fJetContainersData.GetEntries()){ | |
153 | AliDebug(1, "Jet containers for Data:"); | |
154 | TObjString *contname(NULL); | |
155 | TIter contDataIter(&fJetContainersData); | |
156 | while((contname = dynamic_cast<TObjString *>(contDataIter.Next()))) | |
157 | AliDebug(1, Form("Next container: %s", contname->String().Data())); | |
158 | } | |
159 | if(fJetCollArray.GetEntries()){ | |
160 | AliDebug(1, "Jet containers attached to this task:"); | |
161 | AliJetContainer *cont(NULL); | |
162 | TIter contIter(&fJetCollArray); | |
163 | while((cont = dynamic_cast<AliJetContainer *>(contIter.Next()))) | |
164 | AliDebug(1, Form("Container: %s", cont->GetName())); | |
165 | } | |
166 | ||
9f248f6f | 167 | std::map<std::string, std::string> triggerCombinations; |
168 | const char *triggernames[12] = {"MinBias", "EMCJHigh", "EMCJLow", "EMCGHigh", "EMCGLow", "NoEMCal", "EMCHighBoth", "EMCHighGammaOnly", "EMCHighJetOnly", "EMCLowBoth", "EMCLowGammaOnly", "EMCLowJetOnly"}, | |
169 | *bitnames[4] = {"CINT7", "EMC7", "kEMCEGA", "kEMCEJE"}; | |
170 | // Define axes for the trigger correlation histogram | |
171 | const TAxis *triggeraxis[5]; memset(triggeraxis, 0, sizeof(const TAxis *) * 5); | |
172 | const TAxis *bitaxes[4]; memset(bitaxes, 0, sizeof(TAxis *) * 4); | |
173 | const char *binlabels[2] = {"OFF", "ON"}; | |
174 | TAxis mytrgaxis[5], mybitaxis[4]; | |
175 | for(int itrg = 0; itrg < 5; ++itrg){ | |
176 | DefineAxis(mytrgaxis[itrg], triggernames[itrg], triggernames[itrg], 2, -0.5, 1.5, binlabels); | |
177 | triggeraxis[itrg] = mytrgaxis+itrg; | |
178 | if(itrg < 4){ | |
179 | DefineAxis(mybitaxis[itrg], bitnames[itrg], bitnames[itrg], 2, -0.5, 1.5, binlabels); | |
180 | bitaxes[itrg] = mybitaxis+itrg; | |
181 | } | |
182 | } | |
183 | // Define names and titles for different triggers in the histogram container | |
184 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[0], "min. bias events")); | |
185 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[1], "jet-triggered events (high threshold)")); | |
186 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[2], "jet-triggered events (low threshold)")); | |
187 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[3], "gamma-triggered events (high threshold)")); | |
188 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[4], "gamma-triggered events (low threshold)")); | |
189 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[5], "non-EMCal-triggered events")); | |
190 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[6], "jet and gamma triggered events (high threshold)")); | |
191 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[7], "exclusively gamma-triggered events (high threshold)")); | |
192 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[8], "exclusively jet-triggered events (high threshold)")); | |
193 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[9], "jet and gamma triggered events (low threshold)")); | |
194 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[10], "exclusively gamma-triggered events (low threshold)")); | |
195 | triggerCombinations.insert(std::pair<std::string,std::string>(triggernames[11], "exclusively-triggered events (low threshold)")); | |
196 | // Define axes for the pt histogram | |
197 | // Dimensions: | |
198 | // 1. pt | |
199 | // 2. eta | |
200 | // 3. phi | |
201 | // 4. vertex | |
202 | // 5. pileup (0 = all events, 1 = after pileup rejection) | |
203 | // 6. track cuts (0 = no cuts; 1 = after std cuts) | |
204 | TArrayD ptbinning, zvertexBinning, etabinning, pileupaxis(3); | |
205 | pileupaxis[0] = -0.5; pileupaxis[1] = 0.5; pileupaxis[2] = 1.5; | |
206 | CreateDefaultPtBinning(ptbinning); | |
207 | CreateDefaultZVertexBinning(zvertexBinning); | |
208 | CreateDefaultEtaBinning(etabinning); | |
209 | TAxis htrackaxes[7]; | |
210 | DefineAxis(htrackaxes[0], "pt", "p_{t} (GeV/c)", ptbinning); | |
211 | DefineAxis(htrackaxes[1], "eta", "#eta", etabinning); | |
212 | DefineAxis(htrackaxes[2], "phi", "#phi", 20, 0, 2 * TMath::Pi()); | |
213 | DefineAxis(htrackaxes[3], "zvertex", "z_{V} (cm)", zvertexBinning); | |
214 | DefineAxis(htrackaxes[4], "pileup", "Pileup rejection", 2, -0.5, 1.5); | |
215 | DefineAxis(htrackaxes[5], "trackcuts", "Track Cuts", (fListTrackCuts ? fListTrackCuts->GetEntries() : 0) + 1, -0.5, (fListTrackCuts ? fListTrackCuts->GetEntries() : 0) + 0.5); | |
216 | DefineAxis(htrackaxes[6], "mbtrigger", "Has MB trigger", 2, -0.5, 1.5); | |
217 | const TAxis *trackaxes[7]; | |
218 | for(int iaxis = 0; iaxis < 7; ++iaxis) trackaxes[iaxis] = htrackaxes + iaxis; | |
219 | TAxis hclusteraxes[4]; | |
220 | DefineAxis(hclusteraxes[0], "energy", "E (GeV)", ptbinning); | |
221 | DefineAxis(hclusteraxes[1], "zvertex", "z_{V} (cm)", zvertexBinning); | |
222 | DefineAxis(hclusteraxes[2], "pileup", "Pileup rejection", 2,1 -0.5, 1.5); | |
223 | DefineAxis(hclusteraxes[3], "mbtrigger", "Has MB trigger", 2, -0.5, 1.5); | |
224 | const TAxis *clusteraxes[4]; | |
225 | for(int iaxis = 0; iaxis < 4; ++iaxis) clusteraxes[iaxis] = hclusteraxes + iaxis; | |
64321abd | 226 | TAxis hpatchenergyaxes[5]; |
99dfdc8c | 227 | DefineAxis(hpatchenergyaxes[0], "energy", "Patch energy (GeV)", 100, 0., 100); |
9f248f6f | 228 | DefineAxis(hpatchenergyaxes[1], "eta", "#eta", etabinning); |
229 | DefineAxis(hpatchenergyaxes[2], "phi", "#phi", 20, 0, 2 * TMath::Pi()); | |
230 | DefineAxis(hpatchenergyaxes[3], "isMain", "Main trigger", 2, -0.5, 1.5); | |
64321abd | 231 | DefineAxis(hpatchenergyaxes[4], "emcalgood", "EMCAL good event", 2, -0.5, 1.5); |
232 | const TAxis *patchenergyaxes[5]; | |
233 | for(int iaxis = 0; iaxis < 5; ++iaxis) patchenergyaxes[iaxis] = hpatchenergyaxes + iaxis; | |
234 | TAxis hpatchampaxes[5]; | |
fe20314d | 235 | DefineAxis(hpatchampaxes[0], "amplitude", "Patch energy (GeV)", 10000, 0., 10000.); |
8ce5fc44 | 236 | DefineAxis(hpatchampaxes[1], "eta", "#eta", etabinning); |
237 | DefineAxis(hpatchampaxes[2], "phi", "#phi", 20, 0, 2 * TMath::Pi()); | |
238 | DefineAxis(hpatchampaxes[3], "isMain", "Main trigger", 2, -0.5, 1.5); | |
64321abd | 239 | DefineAxis(hpatchampaxes[4], "emcalgood", "EMCAL good event", 2, -0.5, 1.5); |
240 | const TAxis *patchampaxes[5]; | |
241 | for(int iaxis = 0; iaxis < 5; ++iaxis) patchampaxes[iaxis] = hpatchampaxes + iaxis; | |
9f248f6f | 242 | try{ |
243 | std::string patchnames[] = {"Level0", "JetHigh", "JetLow", "GammaHigh", "GammaLow"}; | |
244 | for(std::string * triggerpatch = patchnames; triggerpatch < patchnames + sizeof(patchnames)/sizeof(std::string); ++triggerpatch){ | |
64321abd | 245 | fHistos->CreateTHnSparse(Form("Energy%s", triggerpatch->c_str()), Form("Patch energy for %s trigger patches", triggerpatch->c_str()), 5, patchenergyaxes, "s"); |
246 | fHistos->CreateTHnSparse(Form("EnergyRough%s", triggerpatch->c_str()), Form("Rough patch energy for %s trigger patches", triggerpatch->c_str()), 5, patchenergyaxes, "s"); | |
247 | fHistos->CreateTHnSparse(Form("Amplitude%s", triggerpatch->c_str()), Form("Patch amplitude for %s trigger patches", triggerpatch->c_str()), 5, patchampaxes, "s"); | |
9f248f6f | 248 | } |
249 | ||
250 | // Create histogram for MC-truth | |
fc0447f1 | 251 | fHistos->CreateTHnSparse("hMCtrueParticles", "Particle-based histogram for MC-true particles", 5, trackaxes, "s"); |
555bdd7f | 252 | if(fJetCollArray.GetEntries()){ |
253 | for(int irad = 0; irad < kNJetRadii; irad++){ | |
7d1a1e9a | 254 | fHistos->CreateTHnSparse(Form("hMCtrueParticlesRad%02d", int(kJetRadii[irad]*10)), |
555bdd7f | 255 | Form("Particle-based histogram for MC-true particles in Jets with radius %.1f", kJetRadii[irad]*10), 5, trackaxes, "s"); |
256 | } | |
257 | // histogram for isolated particles | |
258 | fHistos->CreateTHnSparse("hMCtrueParticlesIsolated", "Particle-based histogram for isolated MC-true particles", 5, trackaxes, "s"); | |
259 | } | |
9f248f6f | 260 | for(std::map<std::string,std::string>::iterator it = triggerCombinations.begin(); it != triggerCombinations.end(); ++it){ |
261 | const std::string name = it->first, &title = it->second; | |
262 | // Create event-based histogram | |
263 | fHistos->CreateTH2(Form("hEventHist%s", name.c_str()), Form("Event-based data for %s events; pileup rejection; z_{V} (cm)", title.c_str()), pileupaxis, zvertexBinning); | |
264 | // Create track-based histogram | |
fc0447f1 | 265 | fHistos->CreateTHnSparse(Form("hTrackHist%s", name.c_str()), Form("Track-based data for %s events", title.c_str()), 7, trackaxes, "s"); |
555bdd7f | 266 | fHistos->CreateTHnSparse(Form("hTrackInAcceptanceHist%s", name.c_str()), Form("Track-based data for %s events for tracks matched to EMCal clusters", title.c_str()), 7, trackaxes, "s"); |
fc0447f1 | 267 | fHistos->CreateTHnSparse(Form("hMCTrackHist%s", name.c_str()), Form("Track-based data for %s events with MC kinematics", title.c_str()), 7, trackaxes, "s"); |
555bdd7f | 268 | fHistos->CreateTHnSparse(Form("hMCTrackInAcceptanceHist%s", name.c_str()), Form("Track-based data for %s events with MC kinematics for tracks matched to EMCal clusters", title.c_str()), 7, trackaxes, "s"); |
9f248f6f | 269 | // Create cluster-based histogram (Uncalibrated and calibrated clusters) |
fc0447f1 | 270 | fHistos->CreateTHnSparse(Form("hClusterCalibHist%s", name.c_str()), Form("Calib. cluster-based histogram for %s events", title.c_str()), 4, clusteraxes, "s"); |
271 | fHistos->CreateTHnSparse(Form("hClusterUncalibHist%s", name.c_str()), Form("Uncalib. cluster-based histogram for %s events", title.c_str()), 4, clusteraxes, "s"); | |
555bdd7f | 272 | if(fJetCollArray.GetEntries()){ |
273 | // Create histograms for particles in different jetfinder radii, only track-based | |
274 | for(int irad = 0; irad < kNJetRadii; irad++){ | |
275 | fHistos->CreateTHnSparse(Form("hTrackHist%sRad%02d", name.c_str(), int(kJetRadii[irad]*10)), | |
276 | Form("Track-based data for %s events for tracks in jets with Radius %.1f", title.c_str(), kJetRadii[irad]), 7, trackaxes, "s"); | |
277 | fHistos->CreateTHnSparse(Form("hTrackInAcceptanceHist%sRad%02d", name.c_str(), int(kJetRadii[irad]*10)), | |
278 | Form("Track-based data for %s events for tracks matched to EMCal clusters in jets with Radius %.1f", title.c_str(), kJetRadii[irad]), | |
279 | 7, trackaxes, "s"); | |
280 | fHistos->CreateTHnSparse(Form("hMCTrackHist%sRad%02d", name.c_str(), int(kJetRadii[irad]*10)), | |
281 | Form("Track-based data for %s events with MC kinematics for tracks in jets with Radius %.1f", title.c_str(), kJetRadii[irad]), | |
282 | 7, trackaxes, "s"); | |
283 | fHistos->CreateTHnSparse(Form("hMCTrackInAcceptanceHist%sRad%02d", name.c_str(), int(kJetRadii[irad]*10)), | |
284 | Form("Track-based data for %s events with MC kinematics for tracks matched to EMCal clusters in jets with Radius %.1f", title.c_str(), kJetRadii[irad]), | |
285 | 7, trackaxes, "s"); | |
286 | fHistos->CreateTHnSparse(Form("hClusterCalibHist%sRad%02d", name.c_str(), int(kJetRadii[irad]*10)), | |
287 | Form("Calib. cluster-based histogram for %s events for clusters in Jets with radius %.1f", title.c_str(), kJetRadii[irad]), 4, clusteraxes, "s"); | |
288 | } | |
289 | } | |
290 | // Create also histograms for isolated particles | |
291 | fHistos->CreateTHnSparse(Form("hTrackHist%sIsolated", name.c_str()), Form("Track-based data for %s events for isolated tracks", title.c_str()), 7, trackaxes, "s"); | |
292 | fHistos->CreateTHnSparse(Form("hTrackInAcceptanceHist%sIsolated", name.c_str()), Form("Track-based data for %s events for isolated tracks matched to EMCal clusters", title.c_str()), 7, trackaxes, "s"); | |
293 | fHistos->CreateTHnSparse(Form("hMCTrackHist%sIsolated", name.c_str()), Form("Track-based data for %s events with MC kinematics for isolated tracks", title.c_str()), 7, trackaxes, "s"); | |
294 | fHistos->CreateTHnSparse(Form("hMCTrackInAcceptanceHist%sIsolated", name.c_str()), Form("Track-based data for %s events with MC kinematics for isolated tracks matched to EMCal clusters", title.c_str()), 7, trackaxes, "s"); | |
9f248f6f | 295 | } |
296 | fHistos->CreateTHnSparse("hEventTriggers", "Trigger type per event", 5, triggeraxis); | |
297 | fHistos->CreateTHnSparse("hEventsTriggerbit", "Trigger bits for the different events", 4, bitaxes); | |
64321abd | 298 | fHistos->CreateTH2("hTRUADC", "ADC value of the TRU", 10001, -0.5, 10000.5, 2, -0.5, 1.5); |
9f248f6f | 299 | } catch (HistoContainerContentException &e){ |
300 | std::stringstream errormessage; | |
301 | errormessage << "Creation of histogram failed: " << e.what(); | |
302 | AliError(errormessage.str().c_str()); | |
303 | } | |
304 | fOutput->Add(fHistos->GetListOfHistograms()); | |
305 | if(fListTrackCuts && fListTrackCuts->GetEntries()){ | |
306 | TIter cutIter(fListTrackCuts); | |
307 | AliEMCalPtTaskVTrackSelection *cutObject(NULL); | |
308 | while((cutObject = dynamic_cast<AliEMCalPtTaskVTrackSelection *>(cutIter()))){ | |
309 | AliESDtrackCuts *cuts = dynamic_cast<AliESDtrackCuts *>(cutObject->GetTrackCuts()); | |
310 | if(cuts){ | |
311 | cuts->DefineHistograms(); | |
312 | fOutput->Add(cuts); | |
313 | } | |
314 | } | |
315 | } | |
316 | PostData(1, fOutput); | |
317 | } | |
318 | ||
319 | //______________________________________________________________________________ | |
320 | Bool_t AliAnalysisTaskPtEMCalTrigger::Run(){ | |
321 | /* | |
322 | * Runs the event loop | |
323 | * | |
324 | * @param option: Additional options | |
325 | */ | |
326 | // Common checks: Have SPD vertex and primary vertex from tracks, and both need to have at least one contributor | |
327 | AliDebug(1,Form("Number of calibrated clusters: %d", fCaloClusters->GetEntries())); | |
328 | AliDebug(1,Form("Number of matched tracks: %d", fTracks->GetEntries())); | |
329 | if(fMCEvent){ | |
330 | // Build always trigger strig from the trigger maker in case of MC | |
331 | fUseTriggersFromTriggerMaker = kTRUE; | |
332 | } | |
333 | ||
64321abd | 334 | Bool_t emcalGood = fInputHandler->IsEventSelected() & AliEmcalPhysicsSelection::kEmcalOk; |
335 | ||
336 | // Fill TRU ADC | |
337 | // for debugging | |
338 | fCaloTriggers = fInputEvent->GetCaloTrigger(fIsEsd ? "EMCALTrigger" : "emcalTrigger"); | |
339 | if(fCaloTriggers){ | |
340 | fCaloTriggers->Reset(); | |
341 | Int_t adcAMP; | |
342 | while(fCaloTriggers->Next()){ | |
343 | fCaloTriggers->GetL1TimeSum(adcAMP); | |
344 | if(adcAMP == -1) adcAMP = 0; // Set ADC amps, which have a negative value, to 0 | |
345 | fHistos->FillTH2("hTRUADC", adcAMP, emcalGood ? 1. : 0.); | |
346 | } | |
347 | } | |
348 | ||
9f248f6f | 349 | // Loop over trigger patches, fill patch energy |
350 | AliEmcalTriggerPatchInfo *triggerpatch(NULL); | |
351 | TIter patchIter(this->fTriggerPatchInfo); | |
352 | while((triggerpatch = dynamic_cast<AliEmcalTriggerPatchInfo *>(patchIter()))){ | |
64321abd | 353 | double triggerpatchinfo[5] = {triggerpatch->GetPatchE(), triggerpatch->GetEtaGeo(), triggerpatch->GetPhiGeo(), triggerpatch->IsMainTrigger() ? 1. : 0., emcalGood ? 1. : 0.}; |
354 | double triggerpatchinfoamp[5] = {static_cast<double>(triggerpatch->GetADCAmp()), triggerpatch->GetEtaGeo(), triggerpatch->GetPhiGeo(), triggerpatch->IsMainTrigger() ? 1. : 0., emcalGood ? 1. : 0.}; | |
355 | double triggerpatchinfoer[5] = {triggerpatch->GetADCAmpGeVRough(), triggerpatch->GetEtaGeo(), triggerpatch->GetPhiGeo(), triggerpatch->IsMainTrigger() ? 1. : 0., emcalGood ? 1. : 0.}; | |
9f248f6f | 356 | if(triggerpatch->IsJetHigh()){ |
357 | fHistos->FillTHnSparse("EnergyJetHigh", triggerpatchinfo); | |
358 | fHistos->FillTHnSparse("AmplitudeJetHigh", triggerpatchinfoamp); | |
359 | fHistos->FillTHnSparse("EnergyRoughJetHigh", triggerpatchinfoer); | |
360 | } | |
361 | if(triggerpatch->IsJetLow()){ | |
362 | fHistos->FillTHnSparse("EnergyJetLow", triggerpatchinfo); | |
8ce5fc44 | 363 | fHistos->FillTHnSparse("AmplitudeJetLow", triggerpatchinfoamp); |
364 | fHistos->FillTHnSparse("EnergyRoughJetLow", triggerpatchinfoer); | |
9f248f6f | 365 | } |
366 | if(triggerpatch->IsGammaHigh()){ | |
367 | fHistos->FillTHnSparse("EnergyGammaHigh", triggerpatchinfo); | |
8ce5fc44 | 368 | fHistos->FillTHnSparse("AmplitudeGammaHigh", triggerpatchinfoamp); |
369 | fHistos->FillTHnSparse("EnergyRoughGammaHigh", triggerpatchinfoer); | |
9f248f6f | 370 | } |
371 | if(triggerpatch->IsGammaLow()){ | |
372 | fHistos->FillTHnSparse("EnergyGammaLow", triggerpatchinfo); | |
8ce5fc44 | 373 | fHistos->FillTHnSparse("AmplitudeGammaLow", triggerpatchinfoamp); |
374 | fHistos->FillTHnSparse("EnergyRoughGammaLow", triggerpatchinfoer); | |
9f248f6f | 375 | } |
376 | if(triggerpatch->IsLevel0()){ | |
377 | fHistos->FillTHnSparse("EnergyLevel0", triggerpatchinfo); | |
8ce5fc44 | 378 | fHistos->FillTHnSparse("AmplitudeLevel0", triggerpatchinfoamp); |
379 | fHistos->FillTHnSparse("EnergyRoughLevel0", triggerpatchinfoer); | |
9f248f6f | 380 | } |
381 | } | |
382 | ||
383 | const AliVVertex *vtxTracks = fInputEvent->GetPrimaryVertex(), | |
384 | *vtxSPD = GetSPDVertex(); | |
385 | if(!(vtxTracks && vtxSPD)) return false; | |
2ae55f0a | 386 | if(!fVertexRange.IsInRange(vtxTracks->GetZ())) return false; |
9f248f6f | 387 | if(vtxTracks->GetNContributors() < 1 || vtxSPD->GetNContributors() < 1) return false; |
388 | ||
389 | double triggers[5]; memset(triggers, 0, sizeof(double) * 5); | |
390 | double triggerbits[4]; memset(triggerbits, 0, sizeof(double) * 4); | |
391 | if(fInputHandler->IsEventSelected() & AliVEvent::kINT7){ | |
392 | triggers[0] = 1.; | |
393 | triggerbits[0] = 1.; | |
394 | } | |
395 | ||
396 | // check triggerbits | |
397 | if(fInputHandler->IsEventSelected() & AliVEvent::kEMC7){ | |
398 | triggerbits[1] = 1.; | |
399 | } | |
400 | if(fInputHandler->IsEventSelected() & AliVEvent::kEMCEGA){ | |
401 | triggerbits[2] = 1.; | |
402 | } | |
403 | if(fInputHandler->IsEventSelected() & AliVEvent::kEMCEJE){ | |
404 | triggerbits[3] = 1.; | |
405 | } | |
406 | try{ | |
407 | fHistos->FillTHnSparse("hEventsTriggerbit", triggerbits); | |
408 | } catch(HistoContainerContentException &e) { | |
409 | std::stringstream errormessage; | |
410 | errormessage << "Filling of histogram failed: " << e.what(); | |
411 | AliError(errormessage.str().c_str()); | |
412 | } | |
413 | ||
414 | std::vector<std::string> triggerstrings; | |
415 | // EMCal-triggered event, distinguish types | |
416 | TString trgstr(fUseTriggersFromTriggerMaker ? BuildTriggerString() : fInputEvent->GetFiredTriggerClasses()); | |
417 | AliDebug(1, Form("Triggerstring: %s\n", trgstr.Data())); | |
418 | if(trgstr.Contains("EJ1")){ | |
419 | triggerstrings.push_back("EMCJHigh"); | |
420 | triggers[1] = 1; | |
421 | if(trgstr.Contains("EG1")) | |
422 | triggerstrings.push_back("EMCHighBoth"); | |
423 | else | |
424 | triggerstrings.push_back("EMCHighJetOnly"); | |
425 | } | |
426 | if(trgstr.Contains("EJ2")){ | |
427 | triggerstrings.push_back("EMCJLow"); | |
428 | triggers[2] = 1; | |
429 | if(trgstr.Contains("EG2")) | |
430 | triggerstrings.push_back("EMCLowBoth"); | |
431 | else | |
432 | triggerstrings.push_back("EMCLowJetOnly"); | |
433 | } | |
434 | if(trgstr.Contains("EG1")){ | |
435 | triggerstrings.push_back("EMCGHigh"); | |
436 | triggers[3] = 1; | |
437 | if(!trgstr.Contains("EJ1")) | |
438 | triggerstrings.push_back("EMCHighGammaOnly"); | |
439 | } | |
440 | if(trgstr.Contains("EG2")){ | |
441 | triggerstrings.push_back("EMCGLow"); | |
442 | triggers[4] = 1; | |
443 | if(!trgstr.Contains("EJ2")) | |
444 | triggerstrings.push_back("EMCLowGammaOnly"); | |
445 | } | |
446 | ||
447 | try{ | |
448 | fHistos->FillTHnSparse("hEventTriggers", triggers); | |
449 | } catch (HistoContainerContentException &e){ | |
450 | std::stringstream errormessage; | |
451 | errormessage << "Filling of histogram failed: " << e.what(); | |
452 | AliError(errormessage.str().c_str()); | |
453 | } | |
454 | ||
455 | // apply event selection: Combine the Pileup cut from SPD with the other pA Vertex selection cuts. | |
456 | bool isPileupEvent = fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5.); | |
457 | isPileupEvent = isPileupEvent || (TMath::Abs(vtxTracks->GetZ() - vtxSPD->GetZ()) > 0.5); | |
458 | double covSPD[6]; vtxSPD->GetCovarianceMatrix(covSPD); | |
459 | isPileupEvent = isPileupEvent || (TString(vtxSPD->GetTitle()).Contains("vertexer:Z") && TMath::Sqrt(covSPD[5]) > 0.25); | |
460 | ||
461 | // Fill event-based histogram | |
462 | const double &zv = vtxTracks->GetZ(); | |
463 | if(triggers[0]) FillEventHist("MinBias", zv, isPileupEvent); | |
464 | if(!triggerstrings.size()) // Non-EMCal-triggered | |
465 | FillEventHist("NoEMCal", zv, isPileupEvent); | |
466 | else{ | |
467 | // EMCal-triggered events | |
468 | for(std::vector<std::string>::iterator it = triggerstrings.begin(); it != triggerstrings.end(); ++it) | |
469 | FillEventHist(it->c_str(), zv, isPileupEvent); | |
470 | } | |
471 | ||
555bdd7f | 472 | const AliEmcalJet *foundJet(NULL); |
473 | char histname[1024]; | |
474 | std::vector<double> coneradii; | |
475 | AliJetContainer *jetContMC(NULL), *jetContData(NULL); | |
476 | ||
9f248f6f | 477 | // Fill MC truth |
7d1a1e9a | 478 | AliVParticle *part(NULL); |
9f248f6f | 479 | if(fMCEvent){ |
7d1a1e9a | 480 | if(fJetCollArray.GetEntries() && |
481 | (jetContMC = dynamic_cast<AliJetContainer *>(fJetCollArray.FindObject((static_cast<TObjString *>(fJetContainersMC.At(0)))->String().Data())))){ | |
482 | // In case we have a jet array we loop over all MC selected particles in the particle container of the jet array | |
483 | TIter particles(jetContMC->GetParticleContainer()->GetArray()); | |
484 | while((part = dynamic_cast<AliVParticle *>(particles()))){ | |
485 | if(part->Charge() == 0) continue; | |
486 | if(!fEtaRange.IsInRange(part->Eta())) continue; | |
487 | if(!fPtRange.IsInRange(part->Pt())) continue; | |
488 | FillMCParticleHist("hMCtrueParticles", part, zv, isPileupEvent); | |
555bdd7f | 489 | |
555bdd7f | 490 | /* |
491 | * Jet part: Find track in jet container, | |
492 | * check according to number of particles in jet, and | |
493 | * check for different cone radii | |
494 | */ | |
495 | foundJet = FoundTrackInJet(part, jetContMC); | |
496 | if(foundJet && foundJet->GetNumberOfConstituents() > 1){ | |
497 | for(int irad = 0; irad < kNJetRadii; irad++){ | |
498 | if(IsInRadius(part, foundJet, kJetRadii[irad])){ | |
499 | sprintf(histname, "hMCtrueParticlesRad%02d", int(kJetRadii[irad]*10)); | |
500 | FillMCParticleHist(histname, part, zv, isPileupEvent); | |
501 | } | |
502 | } | |
503 | } else { | |
504 | // isolated track | |
505 | FillMCParticleHist("hMCtrueParticlesIsolated", part, zv, isPileupEvent); | |
506 | } | |
507 | } | |
7d1a1e9a | 508 | } else { |
509 | // Use MC Event | |
510 | for(int ipart = 0; ipart < fMCEvent->GetNumberOfTracks(); ipart++){ | |
511 | // Select only physical primary particles | |
512 | part = fMCEvent->GetTrack(ipart); | |
513 | if(part->Charge() == 0) continue; | |
514 | if(!fEtaRange.IsInRange(part->Eta())) continue; | |
515 | if(!fPtRange.IsInRange(part->Pt())) continue; | |
516 | if(!fMCEvent->IsPhysicalPrimary(ipart)) continue; | |
517 | FillMCParticleHist("hMCtrueParticles", part, zv, isPileupEvent); | |
518 | } | |
9f248f6f | 519 | } |
520 | } | |
521 | ||
522 | AliVTrack *track(NULL); | |
7d1a1e9a | 523 | AliPicoTrack *picoTrack(NULL); |
524 | TObject *containerObject(NULL); | |
9f248f6f | 525 | // Loop over all tracks (No cuts applied) |
526 | TIter allTrackIter(fTracks); | |
7d1a1e9a | 527 | while((containerObject = dynamic_cast<TObject *>(allTrackIter()))){ |
528 | if((picoTrack = dynamic_cast<AliPicoTrack *>(containerObject))){ | |
529 | track = picoTrack->GetTrack(); | |
530 | } else | |
531 | track = dynamic_cast<AliVTrack *>(containerObject); | |
9f248f6f | 532 | if(!IsTrueTrack(track)) continue; |
533 | if(!fEtaRange.IsInRange(track->Eta())) continue; | |
534 | if(!fPtRange.IsInRange(track->Pt())) continue; | |
535 | if(triggers[0]) FillTrackHist("MinBias", track, zv, isPileupEvent, 0, triggers[0]); | |
536 | if(!triggerstrings.size()) // Non-EMCal-triggered | |
537 | FillTrackHist("NoEMCal", track, zv, isPileupEvent, 0, triggers[0]); | |
538 | else { | |
539 | // EMCal-triggered events | |
540 | for(std::vector<std::string>::iterator it = triggerstrings.begin(); it != triggerstrings.end(); ++it) | |
541 | FillTrackHist(it->c_str(), track, zv, isPileupEvent, 0, triggers[0]); | |
542 | } | |
543 | } | |
544 | ||
545 | // Now apply track selection cuts | |
546 | // allow for several track selections to be tested at the same time | |
547 | // each track selection gets a different cut ID starting from 1 | |
548 | // cut ID 0 is reserved for the case of no cuts | |
555bdd7f | 549 | // In case we have a jet container attached we check whether the track is |
550 | // found in a jet, and check for different cone radii around a jet | |
9f248f6f | 551 | if(fListTrackCuts && fListTrackCuts->GetEntries()){ |
552 | for(int icut = 0; icut < fListTrackCuts->GetEntries(); icut++){ | |
553 | AliEMCalPtTaskVTrackSelection *trackSelection = static_cast<AliEMCalPtTaskVTrackSelection *>(fListTrackCuts->At(icut)); | |
554 | TIter trackIter(trackSelection->GetAcceptedTracks(fTracks)); | |
555 | while((track = dynamic_cast<AliVTrack *>(trackIter()))){ | |
556 | //if(!IsTrueTrack(track)) continue; | |
557 | if(!fEtaRange.IsInRange(track->Eta())) continue; | |
558 | if(!fPtRange.IsInRange(track->Pt())) continue; | |
555bdd7f | 559 | coneradii.clear(); |
560 | if(this->fJetCollArray.GetEntries() && | |
561 | (jetContData = dynamic_cast<AliJetContainer *>(this->fJetCollArray.FindObject((static_cast<TObjString *>(this->fJetContainersData.At(0)))->String().Data())))){ | |
562 | foundJet = FoundTrackInJet(track, jetContData); | |
563 | if(foundJet){ | |
564 | for(int irad = 0; irad < kNJetRadii; irad++){ | |
565 | if(IsInRadius(track, foundJet, kJetRadii[irad])) coneradii.push_back(kJetRadii[irad]); | |
566 | } | |
567 | } | |
568 | } | |
569 | if(triggers[0]){ | |
570 | FillTrackHist("MinBias", track, zv, isPileupEvent, icut + 1, triggers[0]); | |
571 | if(coneradii.size()){ | |
572 | for(std::vector<double>::iterator radIter = coneradii.begin(); radIter != coneradii.end(); radIter++) | |
573 | FillTrackHist("MinBias", track, zv, isPileupEvent, icut + 1, triggers[0], *radIter); | |
574 | } | |
575 | } | |
576 | if(!triggerstrings.size()){ // Non-EMCal-triggered | |
9f248f6f | 577 | FillTrackHist("NoEMCal", track, zv, isPileupEvent, icut + 1, triggers[0]); |
555bdd7f | 578 | for(std::vector<double>::iterator radIter = coneradii.begin(); radIter != coneradii.end(); radIter++) |
579 | FillTrackHist("NoEMCal", track, zv, isPileupEvent, icut + 1, triggers[0], *radIter); | |
580 | } else { | |
9f248f6f | 581 | // EMCal-triggered events |
555bdd7f | 582 | for(std::vector<std::string>::iterator it = triggerstrings.begin(); it != triggerstrings.end(); ++it){ |
9f248f6f | 583 | FillTrackHist(it->c_str(), track, zv, isPileupEvent, icut + 1, triggers[0]); |
555bdd7f | 584 | for(std::vector<double>::iterator radIter = coneradii.begin(); radIter != coneradii.end(); radIter++) |
585 | FillTrackHist(it->c_str(), track, zv, isPileupEvent, icut + 1, triggers[0], *radIter); | |
586 | } | |
9f248f6f | 587 | } |
588 | } | |
589 | } | |
590 | } | |
591 | ||
592 | // Next step we loop over the (uncalibrated) emcal clusters and fill histograms with the cluster energy | |
593 | const AliVCluster *clust(NULL); | |
594 | for(int icl = 0; icl < fInputEvent->GetNumberOfCaloClusters(); icl++){ | |
595 | clust = fInputEvent->GetCaloCluster(icl); | |
596 | if(!clust->IsEMCAL()) continue; | |
caff8cd1 | 597 | if(!fEnergyRange.IsInRange(clust->E())) continue; |
555bdd7f | 598 | if(triggers[0]) FillClusterHist("hClusterUncalibHistMinBias", clust, zv, isPileupEvent, triggers[0]); |
599 | if(!triggerstrings.size()){ // Non-EMCal-triggered | |
600 | FillClusterHist("hClusterUncalibHistNoEMCal", clust, zv, isPileupEvent, triggers[0]); | |
601 | } else{ | |
9f248f6f | 602 | for(std::vector<std::string>::iterator it = triggerstrings.begin(); it != triggerstrings.end(); ++it){ |
555bdd7f | 603 | sprintf(histname, "hClusterUncalibHist%s", it->c_str()); |
604 | FillClusterHist(histname, clust, zv, isPileupEvent, triggers[0]); | |
9f248f6f | 605 | } |
606 | } | |
607 | } | |
608 | ||
609 | if(fCaloClusters){ | |
610 | TIter clustIter(fCaloClusters); | |
611 | while((clust = dynamic_cast<const AliVCluster *>(clustIter()))){ | |
612 | if(!clust->IsEMCAL()) continue; | |
caff8cd1 | 613 | if(!fEnergyRange.IsInRange(clust->E())) continue; |
555bdd7f | 614 | if(triggers[0]) FillClusterHist("hClusterCalibHistMinBias", clust, zv, isPileupEvent, triggers[0]); |
9f248f6f | 615 | if(!triggerstrings.size()) // Non-EMCal-triggered |
555bdd7f | 616 | FillClusterHist("hClusterCalibHistNoEMCal", clust, zv, isPileupEvent, triggers[0]); |
9f248f6f | 617 | else{ |
618 | for(std::vector<std::string>::iterator it = triggerstrings.begin(); it != triggerstrings.end(); ++it){ | |
555bdd7f | 619 | sprintf(histname, "hClusterCalibHist%s", it->c_str()); |
620 | FillClusterHist(histname, clust, zv, isPileupEvent, triggers[0]); | |
9f248f6f | 621 | } |
622 | } | |
623 | } | |
624 | } | |
625 | ||
626 | PostData(1, fOutput); | |
627 | return true; | |
628 | } | |
629 | ||
630 | //______________________________________________________________________________ | |
631 | void AliAnalysisTaskPtEMCalTrigger::CreateDefaultPtBinning(TArrayD &binning) const{ | |
632 | /* | |
633 | * Creating the default pt binning. | |
634 | * | |
635 | * @param binning: Array where to store the results. | |
636 | */ | |
637 | std::vector<double> mybinning; | |
638 | std::map<double,double> definitions; | |
639 | definitions.insert(std::pair<double,double>(2.5, 0.1)); | |
640 | definitions.insert(std::pair<double,double>(7., 0.25)); | |
641 | definitions.insert(std::pair<double,double>(15., 0.5)); | |
642 | definitions.insert(std::pair<double,double>(25., 1.)); | |
643 | definitions.insert(std::pair<double,double>(40., 2.5)); | |
2ae55f0a | 644 | definitions.insert(std::pair<double,double>(50., 5.)); |
645 | definitions.insert(std::pair<double,double>(100., 10.)); | |
9f248f6f | 646 | double currentval = 0; |
647 | for(std::map<double,double>::iterator id = definitions.begin(); id != definitions.end(); ++id){ | |
648 | double limit = id->first, binwidth = id->second; | |
2ae55f0a | 649 | while(currentval < limit){ |
9f248f6f | 650 | currentval += binwidth; |
651 | mybinning.push_back(currentval); | |
652 | } | |
653 | } | |
654 | binning.Set(mybinning.size()); | |
655 | int ib = 0; | |
656 | for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it) | |
657 | binning[ib++] = *it; | |
658 | } | |
659 | ||
660 | //______________________________________________________________________________ | |
661 | void AliAnalysisTaskPtEMCalTrigger::CreateDefaultZVertexBinning(TArrayD &binning) const { | |
662 | /* | |
663 | * Creating default z-Vertex binning. | |
664 | * | |
665 | * @param binning: Array where to store the results. | |
666 | */ | |
667 | std::vector<double> mybinning; | |
2ae55f0a | 668 | double currentval = -10; |
9f248f6f | 669 | mybinning.push_back(currentval); |
2ae55f0a | 670 | while(currentval <= 10.){ |
671 | currentval += 5.; | |
9f248f6f | 672 | mybinning.push_back(currentval); |
673 | } | |
674 | binning.Set(mybinning.size()); | |
675 | int ib = 0; | |
676 | for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it) | |
677 | binning[ib++] = *it; | |
678 | } | |
679 | ||
680 | //______________________________________________________________________________ | |
681 | void AliAnalysisTaskPtEMCalTrigger::CreateDefaultEtaBinning(TArrayD& binning) const { | |
682 | /* | |
683 | * Creating default z-Vertex binning. | |
684 | * | |
685 | * @param binning: Array where to store the results. | |
686 | */ | |
687 | std::vector<double> mybinning; | |
688 | double currentval = -0.8; | |
689 | mybinning.push_back(currentval); | |
690 | while(currentval <= 0.8){ | |
691 | currentval += 0.1; | |
692 | mybinning.push_back(currentval); | |
693 | } | |
694 | binning.Set(mybinning.size()); | |
695 | int ib = 0; | |
696 | for(std::vector<double>::iterator it = mybinning.begin(); it != mybinning.end(); ++it) | |
697 | binning[ib++] = *it; | |
698 | } | |
699 | ||
700 | //______________________________________________________________________________ | |
701 | void AliAnalysisTaskPtEMCalTrigger::DefineAxis(TAxis& axis, const char* name, | |
702 | const char* title, const TArrayD& binning, const char** labels) { | |
703 | /* | |
704 | * Define an axis with a given binning | |
705 | * | |
706 | * @param axis: Axis to be defined | |
707 | * @param name: Name of the axis | |
708 | * @param title: Title of the axis | |
709 | * @param binning: axis binning | |
710 | * @param labels (@optional): array of bin labels | |
711 | */ | |
712 | axis.Set(binning.GetSize()-1, binning.GetArray()); | |
713 | axis.SetName(name); | |
714 | axis.SetTitle(title); | |
715 | if(labels){ | |
716 | for(int ib = 1; ib <= axis.GetNbins(); ++ib) | |
717 | axis.SetBinLabel(ib, labels[ib-1]); | |
718 | } | |
719 | } | |
720 | ||
721 | //______________________________________________________________________________ | |
722 | void AliAnalysisTaskPtEMCalTrigger::DefineAxis(TAxis& axis, const char* name, | |
723 | const char* title, int nbins, double min, double max, | |
724 | const char** labels) { | |
725 | /* | |
726 | * Define an axis with number of bins from min to max | |
727 | * | |
728 | * @param axis: Axis to be defined | |
729 | * @param name: Name of the axis | |
730 | * @param title: Title of the axis | |
731 | * @param nbins: Number of bins | |
732 | * @param min: lower limit of the axis | |
733 | * @param max: upper limit of the axis | |
734 | * @param labels (@optional): array of bin labels | |
735 | */ | |
736 | axis.Set(nbins, min, max); | |
737 | axis.SetName(name); | |
738 | axis.SetTitle(title); | |
739 | if(labels){ | |
740 | for(int ib = 1; ib <= axis.GetNbins(); ++ib) | |
741 | axis.SetBinLabel(ib, labels[ib-1]); | |
742 | } | |
743 | } | |
744 | ||
745 | ||
746 | //______________________________________________________________________________ | |
747 | void AliAnalysisTaskPtEMCalTrigger::FillEventHist(const char* trigger, | |
748 | double vz, bool isPileup) { | |
749 | /* | |
750 | * Fill event-based histogram | |
751 | * | |
752 | * @param trigger: name of the trigger configuration to be processed | |
753 | * @param vz: z-position of the vertex | |
754 | * @param isPileup: signalises if the event is flagged as pileup event | |
755 | */ | |
756 | char histname[1024]; | |
757 | sprintf(histname, "hEventHist%s", trigger); | |
758 | try{ | |
759 | fHistos->FillTH2(histname, 0., vz); | |
760 | } catch (HistoContainerContentException &e){ | |
761 | std::stringstream errormessage; | |
762 | errormessage << "Filling of histogram failed: " << e.what(); | |
763 | AliError(errormessage.str().c_str()); | |
764 | } | |
765 | if(!isPileup){ | |
766 | try{ | |
767 | fHistos->FillTH2(histname, 1., vz); | |
768 | } catch(HistoContainerContentException &e){ | |
769 | std::stringstream errormessage; | |
770 | errormessage << "Filling of histogram failed: " << e.what(); | |
771 | AliError(errormessage.str().c_str()); | |
772 | } | |
773 | } | |
774 | } | |
775 | ||
776 | //______________________________________________________________________________ | |
777 | void AliAnalysisTaskPtEMCalTrigger::FillTrackHist(const char* trigger, | |
555bdd7f | 778 | const AliVTrack* track, double vz, bool isPileup, int cut, bool isMinBias, double jetradius) { |
9f248f6f | 779 | /* |
780 | * Fill track-based histogram with corresponding information | |
781 | * | |
782 | * @param trigger: name of the trigger | |
783 | * @param track: ESD track selected | |
784 | * @param vz: z-position of the vertex | |
785 | * @param isPileup: flag event as pileup event | |
786 | * @param cut: id of the cut (0 = no cut) | |
787 | */ | |
788 | double etasign = fSwapEta ? -1. : 1.; | |
789 | double data[7] = {TMath::Abs(track->Pt()), etasign * track->Eta(), track->Phi(), vz, 0, static_cast<double>(cut), isMinBias ? 1. : 0.}; | |
790 | double dataMC[7] = {0., 0., 0., vz, 0, static_cast<double>(cut), isMinBias ? 1. : 0.}; | |
791 | AliVParticle *assocMC(NULL); | |
792 | if(fMCEvent && (assocMC = fMCEvent->GetTrack(TMath::Abs(track->GetLabel())))){ | |
793 | dataMC[0] = TMath::Abs(assocMC->Pt()); | |
794 | dataMC[1] = etasign * assocMC->Eta(); | |
795 | dataMC[2] = assocMC->Phi(); | |
796 | } | |
797 | char histname[1024], histnameAcc[1024], histnameMC[1024], histnameMCAcc[1024]; | |
798 | sprintf(histname, "hTrackHist%s", trigger); | |
799 | sprintf(histnameAcc, "hTrackInAcceptanceHist%s", trigger); | |
800 | sprintf(histnameMC, "hMCTrackHist%s", trigger); | |
801 | sprintf(histnameMCAcc, "hMCTrackInAcceptanceHist%s", trigger); | |
555bdd7f | 802 | if(jetradius > 0.){ |
803 | char *hnames[] = {histname, histnameAcc, histnameMC, histnameMCAcc}; | |
a745bb5b | 804 | for(unsigned int iname = 0; iname < sizeof(hnames)/sizeof(char *); iname++){ |
555bdd7f | 805 | char *myhname = hnames[iname]; |
806 | sprintf(myhname, "%sRad%02d", myhname, int(jetradius * 10.)); | |
807 | } | |
808 | } | |
9f248f6f | 809 | Bool_t isEMCAL = kFALSE; |
810 | if(track->IsEMCAL()){ | |
811 | // Check if the cluster is matched to only one track | |
812 | AliVCluster *emcclust(NULL); | |
813 | AliDebug(2, Form("cluster id: %d\n", track->GetEMCALcluster())); | |
814 | if(fCaloClusters) { | |
815 | AliDebug(2, "Using calibrated clusters"); | |
816 | emcclust = dynamic_cast<AliVCluster *>(fCaloClusters->At(track->GetEMCALcluster())); | |
817 | } else { | |
818 | AliDebug(2, "Using uncalibrated clusters"); | |
819 | emcclust = fInputEvent->GetCaloCluster(track->GetEMCALcluster()); | |
820 | } | |
821 | if(!emcclust) AliError("Null pointer to EMCal cluster"); | |
822 | if(emcclust && emcclust->GetNTracksMatched() <= 1){ | |
823 | isEMCAL = kTRUE; | |
824 | } | |
825 | } | |
826 | try{ | |
827 | fHistos->FillTHnSparse(histname, data); | |
828 | if(fMCEvent) fHistos->FillTHnSparse(histnameMC, dataMC); | |
829 | if(isEMCAL){ | |
830 | fHistos->FillTHnSparse(histnameAcc, data); | |
831 | if(fMCEvent) fHistos->FillTHnSparse(histnameMCAcc, dataMC); | |
832 | } | |
833 | } catch (HistoContainerContentException &e){ | |
834 | std::stringstream errormessage; | |
835 | errormessage << "Filling of histogram failed: " << e.what(); | |
836 | AliError(errormessage.str().c_str()); | |
837 | } | |
838 | if(!isPileup){ | |
839 | data[4] = 1; | |
840 | dataMC[4] = 1; | |
841 | try{ | |
842 | fHistos->FillTHnSparse(histname, data); | |
843 | if(fMCEvent) fHistos->FillTHnSparse(histnameMC, dataMC); | |
844 | if(isEMCAL){ | |
845 | fHistos->FillTHnSparse(histnameAcc, data); | |
846 | if(fMCEvent) fHistos->FillTHnSparse(histnameMCAcc, dataMC); | |
8738a2f4 | 847 | } |
9f248f6f | 848 | } catch (HistoContainerContentException &e){ |
849 | std::stringstream errormessage; | |
850 | errormessage << "Filling of histogram failed: " << e.what(); | |
851 | AliError(errormessage.str().c_str()); | |
852 | } | |
853 | } | |
854 | } | |
855 | ||
856 | //______________________________________________________________________________ | |
555bdd7f | 857 | void AliAnalysisTaskPtEMCalTrigger::FillClusterHist(const char* histname, |
858 | const AliVCluster* clust, double vz, bool isPileup, bool isMinBias) { | |
9f248f6f | 859 | /* |
860 | * Fill cluster-based histogram with corresponding information | |
861 | * | |
862 | * @param trigger: name of the trigger | |
863 | * @param cluster: the EMCal cluster information | |
864 | * @param vz: z-position of the vertex | |
865 | * @param isPileup: flag event as pileup event | |
866 | */ | |
867 | double data[4] = {clust->E(), vz, 0, isMinBias ? 1. : 0.}; | |
9f248f6f | 868 | try{ |
869 | fHistos->FillTHnSparse(histname, data); | |
870 | } catch (HistoContainerContentException &e){ | |
871 | std::stringstream errormessage; | |
872 | errormessage << "Filling of histogram failed: " << e.what(); | |
873 | AliError(errormessage.str().c_str()); | |
874 | } | |
875 | if(!isPileup){ | |
876 | data[2] = 1.; | |
877 | try{ | |
878 | fHistos->FillTHnSparse(histname, data); | |
879 | } catch (HistoContainerContentException &e){ | |
880 | std::stringstream errormessage; | |
881 | errormessage << "Filling of histogram failed: " << e.what(); | |
882 | AliError(errormessage.str().c_str()); | |
883 | } | |
884 | } | |
885 | } | |
886 | ||
887 | //______________________________________________________________________________ | |
555bdd7f | 888 | void AliAnalysisTaskPtEMCalTrigger::FillMCParticleHist(const char *histname, const AliVParticle * const track, double vz, bool isPileup){ |
9f248f6f | 889 | /* |
890 | * Fill histogram for MC-true particles with the information pt, eta and phi | |
891 | * | |
892 | * @param track: the Monte-Carlo track | |
893 | */ | |
894 | double data[5] = {TMath::Abs(track->Pt()), track->Eta(), track->Phi(), vz, 0.}; | |
555bdd7f | 895 | fHistos->FillTHnSparse(histname, data); |
9f248f6f | 896 | if(!isPileup){ |
897 | data[4] = 1.; | |
555bdd7f | 898 | fHistos->FillTHnSparse(histname, data); |
9f248f6f | 899 | } |
900 | } | |
901 | ||
902 | //______________________________________________________________________________ | |
903 | bool AliAnalysisTaskPtEMCalTrigger::IsTrueTrack(const AliVTrack *const track) const{ | |
904 | /* | |
905 | * Check if the track has an associated MC particle, and that the particle is a physical primary | |
906 | * In case of data we do not do the selection at that step (always return true) | |
907 | * | |
908 | * @param track: Track to check | |
909 | * @result: true primary track (true or false) | |
910 | */ | |
911 | if(!fMCEvent) return true; | |
912 | AliVParticle *mcassociate = fMCEvent->GetTrack(TMath::Abs(track->GetLabel())); | |
913 | if(!mcassociate) return false; | |
914 | return fMCEvent->IsPhysicalPrimary(TMath::Abs(track->GetLabel())); | |
915 | } | |
916 | ||
917 | //______________________________________________________________________________ | |
918 | void AliAnalysisTaskPtEMCalTrigger::AddESDTrackCuts(AliESDtrackCuts* trackCuts) { | |
919 | /* | |
920 | * Add new track cuts to the task | |
921 | * | |
922 | * @param trackCuts: Object of type AliESDtrackCuts | |
923 | */ | |
924 | fListTrackCuts->AddLast(new AliEMCalPtTaskTrackSelectionESD(trackCuts)); | |
925 | } | |
926 | ||
927 | //______________________________________________________________________________ | |
928 | void AliAnalysisTaskPtEMCalTrigger::AddCutsForAOD(AliESDtrackCuts* trackCuts, UInt_t filterbits) { | |
929 | /* | |
930 | * Add new track cuts to the task | |
931 | * | |
932 | * @param trackCuts: Object of type AliESDtrackCuts | |
933 | */ | |
934 | fListTrackCuts->AddLast(new AliEMCalPtTaskTrackSelectionAOD(trackCuts, filterbits)); | |
935 | } | |
936 | ||
937 | ||
938 | //______________________________________________________________________________ | |
939 | TString AliAnalysisTaskPtEMCalTrigger::BuildTriggerString() { | |
940 | /* | |
941 | * Build trigger string from the trigger maker | |
942 | * | |
943 | * @return: blank-separated string of fired trigger classes | |
944 | */ | |
945 | AliDebug(1, "trigger checking"); | |
946 | TString result = ""; | |
947 | if(HasTriggerType(kJ1)) result += "EJ1 "; | |
948 | if(HasTriggerType(kJ2)) result += "EJ2 "; | |
949 | if(HasTriggerType(kG1)) result += "EG1 "; | |
950 | if(HasTriggerType(kG2)) result += "EG2 "; | |
951 | return result; | |
952 | } | |
953 | ||
954 | //______________________________________________________________________________ | |
955 | const AliVVertex* AliAnalysisTaskPtEMCalTrigger::GetSPDVertex() const { | |
956 | /* | |
957 | * Accessor for the SPD vertex, creating transparency for ESDs and AODs | |
958 | * | |
959 | * @return: the spd vertex | |
960 | */ | |
961 | AliESDEvent *esd = dynamic_cast<AliESDEvent *>(fInputEvent); | |
962 | if(esd){ | |
963 | return esd->GetPrimaryVertexSPD(); | |
964 | } else { | |
965 | AliAODEvent *aod = dynamic_cast<AliAODEvent *>(fInputEvent); | |
966 | if(aod){ | |
967 | return aod->GetPrimaryVertexSPD(); | |
968 | } | |
969 | } | |
970 | return NULL; | |
971 | } | |
cdc26d91 | 972 | |
555bdd7f | 973 | //______________________________________________________________________________ |
974 | const AliEmcalJet * AliAnalysisTaskPtEMCalTrigger::FoundTrackInJet( | |
a745bb5b | 975 | const AliVParticle * const track, AliJetContainer *const jets) const |
555bdd7f | 976 | { |
977 | /* | |
978 | * Correlate track to reconstructed jet | |
979 | * | |
980 | * @param track: particle to be checked | |
981 | * @param jets: container of recontructed jets | |
982 | * @return: The matched jet (NULL if not found) | |
983 | */ | |
984 | const AliEmcalJet *result = NULL; | |
a745bb5b | 985 | const AliEmcalJet *tmp = jets->GetNextAcceptJet(0); |
986 | while(tmp){ | |
555bdd7f | 987 | if(TrackInJet(track, tmp, jets->GetParticleContainer())){ |
988 | result = tmp; | |
989 | break; | |
990 | } | |
a745bb5b | 991 | tmp = jets->GetNextAcceptJet(); |
555bdd7f | 992 | } |
993 | return result; | |
994 | } | |
995 | ||
996 | //______________________________________________________________________________ | |
997 | bool AliAnalysisTaskPtEMCalTrigger::IsInRadius(const AliVParticle *const track, const AliEmcalJet *reconstructedJet, Double_t radius) const { | |
998 | /* | |
999 | * Check if track is in radius around a given jet | |
1000 | * | |
1001 | * @param track: Track to check | |
1002 | * @param reconstructed jet: jet to probe | |
1003 | * @param radius: cone radius | |
1004 | * @return: result of the test (true if track is inside the cone radius, false otherwise) | |
1005 | */ | |
1006 | return reconstructedJet->DeltaR(track) < radius; | |
1007 | } | |
1008 | ||
1009 | //______________________________________________________________________________ | |
1010 | bool AliAnalysisTaskPtEMCalTrigger::IsInRadius(const AliVCluster *const clust, const AliEmcalJet *reconstructedJet, Double_t radius) const { | |
1011 | /* | |
1012 | * Check if track is in radius around a given jet | |
1013 | * | |
1014 | * @param track: Track to check | |
1015 | * @param reconstructed jet: jet to probe | |
1016 | * @param radius: cone radius | |
1017 | * @return: result of the test (true if track is inside the cone radius, false otherwise) | |
1018 | */ | |
1019 | double vertexPos[3]; | |
1020 | fInputEvent->GetPrimaryVertex()->GetXYZ(vertexPos); | |
1021 | TLorentzVector clustVect; | |
1022 | clust->GetMomentum(clustVect, vertexPos); | |
1023 | ||
1024 | Double_t dPhi = reconstructedJet->Phi() - clustVect.Phi(); | |
1025 | Double_t dEta = reconstructedJet->Eta() - clustVect.Eta(); | |
1026 | dPhi = TVector2::Phi_mpi_pi ( dPhi ); | |
1027 | ||
1028 | double deltaR = TMath::Sqrt ( dPhi * dPhi + dEta * dEta ); | |
1029 | return deltaR < radius; | |
1030 | } | |
1031 | ||
1032 | //______________________________________________________________________________ | |
1033 | bool AliAnalysisTaskPtEMCalTrigger::TrackInJet(const AliVParticle* const track, | |
1034 | const AliEmcalJet* reconstructedJet, const AliParticleContainer * const particles) const { | |
1035 | /* | |
1036 | * Check whether track is among the jet constituents | |
1037 | * | |
1038 | * @param track: track to check | |
1039 | * @param reconstructedJet: reconstructed jet to check | |
1040 | * @param tracks: container with tracks used for jetfinding | |
1041 | * @return: true if found, false otherwise | |
1042 | */ | |
1043 | bool found = false; | |
7d1a1e9a | 1044 | const AliPicoTrack *picotmp(NULL); |
555bdd7f | 1045 | const AliVParticle *tmp(NULL); |
1046 | for(int ipart = 0; ipart < reconstructedJet->GetNumberOfTracks(); ipart++){ | |
1047 | tmp = dynamic_cast<const AliVParticle *>(reconstructedJet->TrackAt(ipart, particles->GetArray())); | |
7d1a1e9a | 1048 | if((picotmp = dynamic_cast<const AliPicoTrack *>(tmp))) // handle pico tracks |
1049 | tmp = picotmp->GetTrack(); | |
555bdd7f | 1050 | if(!tmp->Compare(track)){ |
1051 | found = true; | |
1052 | break; | |
1053 | } | |
1054 | } | |
1055 | return found; | |
1056 | } | |
1057 | ||
1058 | //______________________________________________________________________________ | |
a745bb5b | 1059 | const AliEmcalJet* AliAnalysisTaskPtEMCalTrigger::FoundClusterInJet(const AliVCluster* const clust, AliJetContainer* const jets) const { |
555bdd7f | 1060 | /* |
1061 | * Check whether a cluster is in a radius around a given jet | |
1062 | * | |
1063 | * @param clust: The cluster to check | |
1064 | * @param reconstructedJet: reconstructed jet to check | |
1065 | * @return: the jet containing the cluster (null otherwise) | |
1066 | */ | |
1067 | const AliEmcalJet *result = NULL; | |
a745bb5b | 1068 | const AliEmcalJet *tmp = jets->GetNextAcceptJet(0); |
1069 | while(tmp){ | |
555bdd7f | 1070 | if(ClusterInJet(clust, tmp, jets->GetClusterContainer())){ |
1071 | result = tmp; | |
1072 | break; | |
1073 | } | |
a745bb5b | 1074 | tmp =jets->GetNextAcceptJet(); |
555bdd7f | 1075 | } |
1076 | return result; | |
1077 | } | |
1078 | ||
1079 | //______________________________________________________________________________ | |
1080 | bool AliAnalysisTaskPtEMCalTrigger::ClusterInJet(const AliVCluster* const clust, | |
1081 | const AliEmcalJet* reconstructedJet, const AliClusterContainer* const clusters) const { | |
1082 | /* | |
1083 | * Check whether cluster is among the jet constituents | |
1084 | * | |
1085 | * @param track: track to check | |
1086 | * @param reconstructedJet: reconstructed jet to check | |
1087 | * @param clusters: the cluster container | |
1088 | * @return: true if found, false otherwise | |
1089 | */ | |
1090 | bool found = false; | |
1091 | const AliVCluster *tmp(NULL); | |
1092 | for(int ipart = 0; ipart < reconstructedJet->GetNumberOfTracks(); ipart++){ | |
1093 | tmp = dynamic_cast<const AliVCluster *>(reconstructedJet->ClusterAt(ipart, clusters->GetArray())); | |
1094 | if(!tmp->Compare(clust)){ | |
1095 | found = true; | |
1096 | break; | |
1097 | } | |
1098 | } | |
1099 | return found; | |
1100 | } | |
1101 | ||
1102 | ||
1103 | //______________________________________________________________________________ | |
86c0a4aa | 1104 | void EMCalTriggerPtAnalysis::AliAnalysisTaskPtEMCalTrigger::AddJetContainerName(const Char_t* contname, Bool_t isMC) { |
555bdd7f | 1105 | /* |
1106 | * Add new Jet input container to the analysis task | |
1107 | * | |
1108 | * @param contname: Name of the container | |
1109 | * @param isMC: Defines whether the container is for MC truth or not | |
1110 | */ | |
1111 | TList &mycontainer = isMC ? fJetContainersMC : fJetContainersData; | |
1112 | mycontainer.Add(new TObjString(contname)); | |
1113 | } | |
bf1cb6ad | 1114 | } |
5bdebc68 | 1115 |