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