limit printouts
[u/mrichter/AliRoot.git] / PWG4 / macros / AddTaskJets.C
CommitLineData
0651dd18 1AliJetReader *CreateJetReader(Char_t *jr); // Common config\r
76f9015f 2AliJetFinder *CreateJetFinder(Char_t *jf,Float_t radius = -1);\r
0651dd18 3\r
76f9015f 4AliAnalysisTaskJets *AddTaskJets(Char_t *jr, Char_t *jf,Float_t radius = -1); // for the new AF\r
0651dd18 5\r
9902ce1e 6AliAnalysisTaskJets *AddTaskJets(){\r
7 // fills the standard "jets" branch in the AOD\r
8 // need the ESDFilter to run before, to access the AODtracks\r
9 // Tracks selected by the first Filter (1<<0)\r
10 // needs to be adapted for different cuts\r
11 \r
12 // UA1 as standard chosen, since it is the most robust and simple JF\r
13 // R = 0.4 suffficient to provide accurate jet axis for correlation studies\r
14 // energy resolution suffers a little\r
15 // Acceptance of jets not limited by the Jet Finder but should be done\r
16 // by user to abs(eta) < 0.5 \r
17\r
18 return AddTaskJets("AOD","UA1",0.4);\r
19\r
20}\r
21\r
22\r
23\r
b486b6c5 24Int_t AddTaskJetsDelta(char *nonStdFile = ""){\r
25\r
26 // Adds a whole set of jet finders all to be written\r
27 // to a delta AOD\r
28 \r
29 // this can in principle be done also with on the fly \r
30 // if we have an ESD input jet fidner task does automatically fetch the ouput aod\r
31\r
32 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
33 if (!mgr) {\r
34 ::Error("AddTaskJetsDelta", "No analysis manager to connect to.");\r
35 return 0;\r
36 } \r
37 \r
38 // Check the analysis type using the event handlers connected to the analysis manager.\r
39 //==============================================================================\r
40 if (!mgr->GetInputEventHandler()) {\r
41 ::Error("AddTaskJetsDelta", "This task requires an input event handler");\r
42 return 0;\r
43 }\r
44\r
45 AliAODHandler *aodh = (AliAODHandler*)mgr->GetOutputEventHandler();\r
46 if (!aodh) {\r
47 ::Error("AddTaskJetsDelta", "This task needs an output event handler");\r
48 return 0;\r
49 } \r
50\r
f6d6265b 51 TString outFile(nonStdFile);\r
b486b6c5 52\r
53 TString type = mgr->GetInputEventHandler()->GetDataType();\r
54\r
55 AliAnalysisTaskJets *jetana = 0;\r
56 Int_t iCount = 0;\r
57\r
58\r
8a4430e1 59 const char *cJF[7] = {"UA1","UA1","UA1","CDF","DA","SISCONE","FASTJET"};\r
60 const Float_t radius[7] = { 0.4, 0.7, 1.0, 0.7, 0.7, 0.4, 0.4};\r
61 const UInt_t flag[7] = { 6, 7, 7, 7, 7, 7, 7};\r
b486b6c5 62 // flag first bit AOD, second bit AODMC2 third bit AODMC2\r
63 // i.e. 7 all, 6 only MC2 and MC\r
64 // this stay at three\r
b486b6c5 65 const char *cReader[3] = {"AOD","AODMC","AODMC2"}; \r
66\r
8a4430e1 67 for(int i = 0; i< 7;i++){\r
b486b6c5 68 for(int ib = 0;ib<3;ib++){\r
69 if(flag[i]&(1<<ib)){\r
70 jetana = AddTaskJets(cReader[ib],cJF[i],radius[i]);\r
71 if(jetana){\r
72 char *cRadius = "";\r
8a4430e1 73 if(radius[i]>0)cRadius = Form("%02d",(int)((radius[i]+0.01)*10.)); // add an offset beacuse of precision\r
74 Printf(cRadius);\r
b486b6c5 75 jetana->SetNonStdBranch(Form("jets%s_%s%s",cReader[ib],cJF[i],cRadius));\r
f6d6265b 76 if(outFile.Length()>0)jetana->SetNonStdOutputFile(outFile.Data());\r
b486b6c5 77 iCount++;\r
78 }\r
79 }\r
80 }\r
81 }\r
82 \r
83 Printf("Added %d JetFinders",iCount);\r
84}\r
85\r
86\r
87\r
88\r
89\r
90\r
76f9015f 91AliAnalysisTaskJets *AddTaskJets(Char_t *jr, Char_t *jf, Float_t radius)\r
672f1183 92{\r
0651dd18 93 // Creates a jet finder task, configures it and adds it to the analysis manager.\r
672f1183 94\r
95 // Get the pointer to the existing analysis manager via the static access method.\r
96 //==============================================================================\r
97 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
98 if (!mgr) {\r
99 ::Error("AddTaskJets", "No analysis manager to connect to.");\r
100 return NULL;\r
101 } \r
102 \r
103 // Check the analysis type using the event handlers connected to the analysis manager.\r
104 //==============================================================================\r
105 if (!mgr->GetInputEventHandler()) {\r
106 ::Error("AddTaskJets", "This task requires an input event handler");\r
107 return NULL;\r
108 }\r
109\r
0511798f 110 AliAODHandler *aodh = (AliAODHandler*)mgr->GetOutputEventHandler();\r
111 if (!aodh) {\r
112 ::Error("AddTaskJets", "This task needs an output event handler");\r
113 return NULL;\r
114 } \r
115\r
116\r
672f1183 117 // Create the task and configure it.\r
118 //===========================================================================\r
119 AliAnalysisTaskJets *jetana;\r
0651dd18 120 AliJetReader *er = CreateJetReader(jr);\r
672f1183 121 // Define jet header and jet finder\r
76f9015f 122 AliJetFinder *jetFinder = CreateJetFinder(jf,radius);\r
0651dd18 123\r
672f1183 124 if (jetFinder){\r
672f1183 125 if (er) jetFinder->SetJetReader(er);\r
126 }\r
127\r
76f9015f 128 char *cRadius = "";\r
8a4430e1 129 if(radius>0)cRadius = Form("%02d",(int)((radius+0.01)*10.));\r
76f9015f 130\r
b486b6c5 131 jetana = new AliAnalysisTaskJets(Form("JetAnalysis%s_%s%s",jr,jf,cRadius));\r
06847475 132 TString type = mgr->GetInputEventHandler()->GetDataType();\r
133 if (type == "AOD") jetana->SetNonStdBranch(Form("jets%s",jf));\r
06847475 134 \r
be60d38b 135 TString c_jr(jr);\r
136 c_jr.ToLower();\r
137 TString c_jf(jf);\r
138 c_jf.ToLower();\r
139\r
b486b6c5 140 AliAnalysisDataContainer *cout_jet = mgr->CreateContainer(Form("jethist_%s_%s%s",c_jr.Data(),c_jf.Data(),cRadius), TList::Class(),\r
561a15fa 141 AliAnalysisManager::kOutputContainer, Form("%s:PWG4_jethist_%s_%s%s",AliAnalysisManager::GetCommonFileName(),\r
142 c_jr.Data(),c_jf.Data(),cRadius));\r
672f1183 143 // Connect jet finder to task.\r
144 jetana->SetJetFinder(jetFinder);\r
145 jetana->SetConfigFile("");\r
19027d3e 146 jetana->SetDebugLevel(0);\r
672f1183 147 mgr->AddTask(jetana);\r
148\r
149 // Create ONLY the output containers for the data produced by the task.\r
150 // Get and connect other common input/output containers via the manager as below\r
151 //==============================================================================\r
152 mgr->ConnectInput (jetana, 0, mgr->GetCommonInputContainer());\r
153// AOD output slot will be used in a different way in future\r
154 mgr->ConnectOutput (jetana, 0, mgr->GetCommonOutputContainer());\r
155 mgr->ConnectOutput (jetana, 1, cout_jet);\r
156 \r
157 return jetana;\r
158}\r
0651dd18 159\r
76f9015f 160AliJetFinder *CreateJetFinder(Char_t *jf,Float_t radius){\r
70144f07 161 AliJetFinder *jetFinder = 0;\r
0651dd18 162\r
163 switch (jf) {\r
164 case "CDF":\r
165 AliCdfJetHeader *jh = new AliCdfJetHeader();\r
166 jh->SetRadius(0.7);\r
8a4430e1 167 if(radius>0)jh->SetRadius(radius); \r
0651dd18 168 jetFinder = new AliCdfJetFinder();\r
0651dd18 169 if (jh) jetFinder->SetJetHeader(jh);\r
170 break;\r
171\r
172 case "DA":\r
173 AliDAJetHeader *jh=new AliDAJetHeader();\r
174 jh->SetComment("DA jet code with default parameters");\r
175 jh->SelectJets(kTRUE);\r
4be01860 176 jh->SetRadius(0.7);\r
8a4430e1 177 if(radius>0)jh->SetRadius(radius);\r
0651dd18 178 jetFinder = new AliDAJetFinder();\r
179 if (jh) jetFinder->SetJetHeader(jh);\r
180 break;\r
181\r
70144f07 182 case "FASTJET":\r
b486b6c5 183 // DEFAULT is ANTI KT\r
70144f07 184 AliFastJetHeaderV1 *jh = new AliFastJetHeaderV1();\r
185 jh->SetRparam(0.4); // setup parameters \r
186 if(radius>0)jh->SetRparam(radius);\r
187 jh->SetAlgorithm(2); // antikt from fastjet/JetDefinition.hh\r
0651dd18 188 jetFinder = new AliFastJetFinder();\r
0651dd18 189 if (jh) jetFinder->SetJetHeader(jh);\r
190 break;\r
191\r
8a4430e1 192 case "FASTKT":\r
193 AliFastJetHeaderV1 *jh = new AliFastJetHeaderV1();\r
194 jh->SetRparam(0.4); // setup parameters \r
195 if(radius>0)jh->SetRparam(radius);\r
196 jh->SetAlgorithm(0); // kt from fastjet/JetDefinition.hh\r
197 jetFinder = new AliFastJetFinder();\r
198 if (jh) jetFinder->SetJetHeader(jh);\r
199 break;\r
200\r
201 case "SISCONE":\r
202 AliSISConeJetHeader * jh = new AliSISConeJetHeader();\r
203\r
204 jh->SetJetEtaMax(1.5);\r
205 jh->SetJetEtaMin(-1.5);\r
206\r
207 //siscone parameters\r
208 jh->SetConeRadius(0.4); // default cone radius\r
209 if(radius>0)jh->SetConeRadius(radius); // cone radius\r
210\r
211 jh->SetOverlapThreshold(0.75); // overlap parameter, between 0 and 1 excluded!! 0.75 value is advised\r
212 jh->SetPtProtojetMin(0); // pt min of protojets\r
213 jh->SetMinJetPt(10); // Ptmin of jets (GeV)\r
214\r
215 //do you want to subtract BG (0 = no, 1 = yes)\r
216 jh->SetBGMode(0);\r
217\r
218 //for background\r
219 jh->SetRapRange( -0.9, 0.9); // rapidity range for subtracting background must be < ghostmaxrap-0.95*R\r
220 jh->SetPhiRange(0 , 2*TMath::Pi()); // phi range for subtracting background\r
221 \r
222 //to determine jets area\r
223 jh->SetBGAlgorithm(1); // algorithm for rho calculation : 0 = kT, 1 = Cambridge\r
224 jh->SetGhostEtaMax(4); // eta max where ghosts are generated \r
225 jh->SetGhostArea(0.05); // area of a ghost \r
226 jh->SetMeanGhostKt(1e-100); // average transverse momentum of the ghosts.\r
227 jh->SetAreaTypeNumber(4); // from 1 to 5 : 1 = active_area, 2 = active_area_explicit_ghosts, 3 = one_ghost_passive_area, 4 = passive_area, 5 = voronoi_area\r
228 jetFinder = new AliSISConeJetFinder();\r
229 if (jh) jetFinder->SetJetHeader(jh);\r
230 break;\r
231\r
0651dd18 232 case "UA1":\r
233 AliUA1JetHeaderV1 *jh=new AliUA1JetHeaderV1();\r
234 jh->SetComment("UA1 jet code with default parameters");\r
ba146e16 235 jh->BackgMode(0);\r
236 jh->SetRadius(0.4);\r
76f9015f 237 if(radius>0)jh->SetRadius(radius);\r
0651dd18 238 jh->SetEtSeed(4.);\r
c40ab9fa 239 jh->SetEtSeed(4.);\r
240 jh->SetNAcceptJets(6);\r
0651dd18 241 jh->SetLegoNbinPhi(432);\r
242 jh->SetLegoNbinEta(274);\r
243 jh->SetLegoEtaMin(-2);\r
244 jh->SetLegoEtaMax(+2);\r
9902ce1e 245 jh->SetMinJetEt(5.);\r
ba146e16 246 jh->SetJetEtaMax(1.5);\r
247 jh->SetJetEtaMin(-1.5);\r
0651dd18 248\r
249 jetFinder = new AliUA1JetFinderV1();\r
250 if (jh) jetFinder->SetJetHeader(jh);\r
251 break;\r
252\r
ba146e16 253 case "UA1MC":\r
254 AliUA1JetHeaderV1 *jh=new AliUA1JetHeaderV1();\r
255 jh->SetComment("UA1 jet code with default MC parameters");\r
256 jh->BackgMode(0);\r
257 jh->SetRadius(1.0);\r
76f9015f 258 if(radius>0)jh->SetRadius(radius);\r
ba146e16 259 jh->SetEtSeed(4.);\r
c40ab9fa 260 jh->SetNAcceptJets(6);\r
ba146e16 261 jh->SetLegoNbinPhi(432);\r
262 jh->SetLegoNbinEta(274);\r
263 jh->SetLegoEtaMin(-2);\r
264 jh->SetLegoEtaMax(+2);\r
9902ce1e 265 jh->SetMinJetEt(5.);\r
ba146e16 266 jh->SetJetEtaMax(1.5);\r
267 jh->SetJetEtaMin(-1.5);\r
ba146e16 268 jetFinder = new AliUA1JetFinderV1();\r
269 if (jh) jetFinder->SetJetHeader(jh);\r
270 break;\r
70144f07 271 default:\r
272 Printf("\n >>>>>>> AddTaskJets Error Wrong jet finder selected\n");\r
273 break;\r
0651dd18 274 }\r
275\r
276 return jetFinder;\r
277\r
278}\r
279\r
280AliJetReader *CreateJetReader(Char_t *jr){\r
281 AliJetReader *er = 0;\r
282\r
283 switch (jr) {\r
284 case "MC":\r
285 AliJetKineReaderHeader *jrh = new AliJetKineReaderHeader();\r
286 jrh->SetComment("MC full Kinematics");\r
287 jrh->SetFastSimTPC(kFALSE);\r
288 jrh->SetFastSimEMCAL(kFALSE);\r
289 jrh->SetPtCut(0.);\r
290 jrh->SetFiducialEta(-2.1,2.1); // to take all MC particles default is 0 .9 \r
291 // Define reader and set its header \r
292 er = new AliJetKineReader();\r
293 er->SetReaderHeader(jrh);\r
294 break;\r
ba146e16 295 case "MC2":\r
296 AliJetKineReaderHeader *jrh = new AliJetKineReaderHeader();\r
5dbee319 297 jrh->SetComment("MC full Kinematics spearate config charged only");\r
ba146e16 298 jrh->SetFastSimTPC(kFALSE);\r
299 jrh->SetFastSimEMCAL(kFALSE);\r
5dbee319 300 jrh->SetChargedOnly(kTRUE);\r
ba146e16 301 jrh->SetPtCut(0.);\r
302 jrh->SetFiducialEta(-2.1,2.1); // to take all MC particles default is 0 .9 \r
303 // Define reader and set its header \r
304 er = new AliJetKineReader();\r
305 er->SetReaderHeader(jrh);\r
306 break;\r
0651dd18 307 case "ESD":\r
308 AliJetESDReaderHeader *jrh = new AliJetESDReaderHeader();\r
309 jrh->SetComment("Testing");\r
310 jrh->SetFirstEvent(0);\r
311 jrh->SetLastEvent(1000);\r
312 jrh->SetPtCut(0.);\r
313 jrh->SetReadSignalOnly(kFALSE);\r
314 // Define reader and set its header \r
315 er = new AliJetESDReader();\r
316 er->SetReaderHeader(jrh);\r
317 break;\r
318\r
319 case "AOD":\r
320 AliJetAODReaderHeader *jrh = new AliJetAODReaderHeader();\r
321 jrh->SetComment("AOD Reader");\r
322 jrh->SetPtCut(0.);\r
8a4430e1 323 jrh->SetTestFilterMask(16); // Change this one for a different set of cuts\r
0651dd18 324 // Define reader and set its header\r
325 er = new AliJetAODReader();\r
326 er->SetReaderHeader(jrh);\r
327 break;\r
0c43aaa9 328 case "AODMC":\r
329 AliJetAODReaderHeader *jrh = new AliJetAODReaderHeader();\r
330 jrh->SetComment("AOD MC Reader");\r
331 jrh->SetPtCut(0.);\r
332 jrh->SetFiducialEta(-2.1,2.1); // to take all MC particles default is 0.9\r
333 jrh->SetReadAODMC(1);// 1 all primary MC , 2 all primary charged\r
334 // Define reader and set its header\r
335 er = new AliJetAODReader();\r
336 er->SetReaderHeader(jrh);\r
337 break;\r
5dbee319 338 case "AODMC2":\r
339 AliJetAODReaderHeader *jrh = new AliJetAODReaderHeader();\r
340 jrh->SetComment("AOD MC Reader");\r
341 jrh->SetPtCut(0.);\r
342 jrh->SetFiducialEta(-2.1,2.1); // to take all MC particles default is 0.9\r
343 jrh->SetReadAODMC(2);// 1 all primary MC , 2 all primary charged\r
344 // Define reader and set its header\r
345 er = new AliJetAODReader();\r
346 er->SetReaderHeader(jrh);\r
347 break;\r
348\r
0651dd18 349 default:\r
350 ::Error("AddTaskJets", "Wrong jet reader selected\n");\r
351 return 0;\r
352 }\r
353\r
354 return er;\r
355\r
356}\r