Several bugs
[u/mrichter/AliRoot.git] / PWG4 / macros / ana.C
CommitLineData
d92b41ad 1/* $Id: $ */
d92b41ad 2//--------------------------------------------------
3// Example macro to do analysis with the
4// analysis classes in PWG4PartCorr
5// Can be executed with Root and AliRoot
6//
7// Pay attention to the options and definitions
8// set in the lines below
9//
10// Author : Gustavo Conesa Balbastre (INFN-LNF)
11//
12//-------------------------------------------------
2244659d 13enum anaModes {mLocal=0, mGRID=3};
14//mLocal = 0: Analyze locally files in your computer
15//mGRID = 3: Analyze files on GRID
d92b41ad 16
17//---------------------------------------------------------------------------
18//Settings to read locally several files, only for "mLocal" mode
19//The different values are default, they can be set with environmental
c8fe2783 20//variables: INDIR, PATTERN, NFILES, respectivelly
54769bc0 21
22char * kInDir = "/user/data/files/";
3e0577a2 23char * kPattern = ""; // Data are in files kInDir/kPattern+i
193828fd 24Int_t kFile = 1; // Number of files for local analysis
d92b41ad 25//---------------------------------------------------------------------------
26//Collection file for grid analysis
27char * kXML = "collection.xml";
28//---------------------------------------------------------------------------
29//Scale histograms from file. Change to kTRUE when xsection file exists
30//Put name of file containing xsection
31//Put number of events per ESD file
32//This is an specific case for normalization of Pythia files.
33const Bool_t kGetXSectionFromFileAndScale = kFALSE ;
34const char * kXSFileName = "pyxsec.root";
54769bc0 35
d92b41ad 36//---------------------------------------------------------------------------
37
193828fd 38//Set some default values, but used values are set in the code!
39
40Bool_t kMC = kFALSE; //With real data kMC = kFALSE
41TString kInputData = "ESD"; //ESD, AOD, MC, deltaAOD
42Int_t kYear = 2011;
43TString kCollision = "pp";
44Bool_t outAOD = kTRUE; //Some tasks doesnt need it.
2244659d 45TString kTreeName;
193828fd 46TString kPass = "";
47char kTrigger[1024];
48Int_t kRun = 0;
7175a03a 49
193828fd 50void ana(Int_t mode=mGRID)
d92b41ad 51{
52 // Main
cd54ca47 53
d92b41ad 54 //--------------------------------------------------------------------
55 // Load analysis libraries
56 // Look at the method below,
57 // change whatever you need for your analysis case
58 // ------------------------------------------------------------------
193828fd 59
2244659d 60 LoadLibraries() ;
193828fd 61 //gSystem->ListLibraries();
2244659d 62
d92b41ad 63 //-------------------------------------------------------------------------------------------------
64 //Create chain from ESD and from cross sections files, look below for options.
c8fe2783 65 //-------------------------------------------------------------------------------------------------
54769bc0 66
193828fd 67 // Set kInputData and kTreeName looking to the kINDIR
cd54ca47 68
193828fd 69 CheckInputData(mode);
cd54ca47 70
193828fd 71 // Check global analysis settings
cd54ca47 72
193828fd 73 CheckEnvironmentVariables();
74
75 printf("*********************************************\n");
76 printf("*** Input data < %s >, pass %s, tree < %s >, MC? < %d > ***\n",kInputData.Data(),kPass.Data(),kTreeName.Data(),kMC);
77 printf("*********************************************\n");
2244659d 78
193828fd 79 TChain * chain = new TChain(kTreeName) ;
d92b41ad 80 TChain * chainxs = new TChain("Xsection") ;
193828fd 81 CreateChain(mode, chain, chainxs);
54769bc0 82
193828fd 83 printf("*********************************************\n");
84 printf("number of entries # %lld, skipped %d\n", chain->GetEntries()) ;
85 printf("*********************************************\n");
86
87 if(!chain){
88 printf("STOP, no chain available\n");
89 return;
90 }
91
92 AliLog::SetGlobalLogLevel(AliLog::kError);//Minimum prints on screen
cd54ca47 93
193828fd 94 //--------------------------------------
95 // Make the analysis manager
96 //-------------------------------------
97 AliAnalysisManager *mgr = new AliAnalysisManager("Manager", "Manager");
98 //AliAnalysisManager::SetUseProgressBar(kTRUE);
99 //mgr->SetSkipTerminate(kTRUE);
100 //mgr->SetNSysInfo(1);
101
102 // MC handler
103 if((kMC || kInputData == "MC") && !kInputData.Contains("AOD")){
104 AliMCEventHandler* mcHandler = new AliMCEventHandler();
105 mcHandler->SetReadTR(kFALSE);//Do not search TrackRef file
106 mgr->SetMCtruthEventHandler(mcHandler);
107 if( kInputData == "MC") {
108 cout<<"MC INPUT EVENT HANDLER"<<endl;
109 mgr->SetInputEventHandler(NULL);
d92b41ad 110 }
193828fd 111 }
112
113 // AOD output handler
114 if(kInputData!="deltaAOD" && outAOD){
115 cout<<"Init output handler"<<endl;
116 AliAODHandler* aodoutHandler = new AliAODHandler();
117 aodoutHandler->SetOutputFileName("aod.root");
118 ////aodoutHandler->SetCreateNonStandardAOD();
119 mgr->SetOutputEventHandler(aodoutHandler);
120 }
121
122 //input
123
124 if(kInputData == "ESD"){
125 // ESD handler
126 AliESDInputHandler *esdHandler = new AliESDInputHandler();
127 esdHandler->SetReadFriends(kFALSE);
128 mgr->SetInputEventHandler(esdHandler);
129 cout<<"ESD handler "<<mgr->GetInputEventHandler()<<endl;
130 }
131 else if(kInputData.Contains("AOD")){
132 // AOD handler
133 AliAODInputHandler *aodHandler = new AliAODInputHandler();
134 mgr->SetInputEventHandler(aodHandler);
135 if(kInputData == "deltaAOD") aodHandler->AddFriend("deltaAODPartCorr.root");
136 cout<<"AOD handler "<<mgr->GetInputEventHandler()<<endl;
137 }
138 //mgr->RegisterExternalFile("deltaAODPartCorr.root");
139 //mgr->SetDebugLevel(1); // For debugging, do not uncomment if you want no messages.
140
141 TString outputFile = AliAnalysisManager::GetCommonFileName();
142
143 //-------------------------------------------------------------------------
144 //Define task, put here any other task that you want to use.
145 //-------------------------------------------------------------------------
146
147 // Physics selection
148 gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
149 AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(kMC);
150
151 // Centrality
152 // gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskCentrality.C");
153 // AliCentralitySelectionTask *taskCentrality = AddTaskCentrality();
154 // taskCentrality->SetPass(2); // remember to set the pass you are processing!!!
155
156 // Simple event counting tasks
157 AddTaskCounter(""); // All
158 AddTaskCounter("MB"); // Min Bias
159 if(!kMC){
160 AddTaskCounter("INT7"); // Min Bias
161 AddTaskCounter("EMC1"); // Trig Th > 1.5 GeV approx
162 AddTaskCounter("EMC7"); // Trig Th > 4-5 GeV
cd54ca47 163 AddTaskCounter("PHOS"); //
193828fd 164 }
165
193828fd 166 gROOT->LoadMacro("AddTaskPartCorr.C"); // $ALICE_ROOT/PWG4/macros
cd54ca47 167
193828fd 168 Bool_t kPrint = kFALSE;
169 Bool_t deltaAOD = kFALSE;
cd54ca47 170
171 // ------
172 // Tracks
173 // ------
193828fd 174
cd54ca47 175 // QA analysis to compare after clusterization and track isolation-correlation analysis
176 AliAnalysisTaskParticleCorrelation *ana = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
177 kYear,kRun,kCollision,"EMC7","");
193828fd 178
179 // -----
180 // EMCAL
cd54ca47 181 // -----
193828fd 182
183 gROOT->LoadMacro("AddTaskEMCALClusterize.C"); // $ALICE_ROOT/PWG4/CaloCalib/macros
184 Bool_t bTrackMatch = kTRUE;
cd54ca47 185 Int_t minEcell = 50; // 50 MeV (10 MeV used in reconstruction)
186 Int_t minEseed = 100; // 100 MeV
187 Int_t dTime = 40; // 40 ns difference in time of cells in cluster, open
188 Int_t wTime = 40; // 40 ns time window of cells in cluster, open
189 TString clTrigger = ""; // Do not select, do Min Bias and triggered
193828fd 190
cd54ca47 191 //Analysis with clusterizer V1
193828fd 192 AliAnalysisTaskEMCALClusterize * clv1 = AddTaskEMCALClusterize(kMC,"V1",clTrigger,kRun,kPass, bTrackMatch,
193 minEcell,minEseed,dTime,wTime);
194
195 TString arrayNameV1(Form("V1_Ecell%d_Eseed%d_DT%d_WT%d",minEcell,minEseed, dTime,wTime));
196 printf("Name of clusterizer array: %s\n",arrayNameV1.Data());
197
198 if(!kMC)
199 {
cd54ca47 200
201 AliAnalysisTaskParticleCorrelation *anav1tr = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
193828fd 202 kYear,kRun,kCollision,"EMC7",arrayNameV1);
54769bc0 203
cd54ca47 204 AliAnalysisTaskParticleCorrelation *anav1mb = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
193828fd 205 kYear,kRun,kCollision,"INT7",arrayNameV1);
206 }
207 else
208 {// No trigger (should be MB, but for single particle productions it does not work)
cd54ca47 209
193828fd 210 AliAnalysisTaskParticleCorrelation *anav1 = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
211 kYear,kRun,kCollision,"",arrayNameV1);
212 }
213
cd54ca47 214
215 //Analysis with clusterizer V2
216 AliAnalysisTaskEMCALClusterize * clv2 = AddTaskEMCALClusterize(kMC,"V2",clTrigger,kRun,kPass, bTrackMatch,
217 minEcell,minEseed,dTime,wTime);
218
219 TString arrayNameV2(Form("V2_Ecell%d_Eseed%d_DT%d_WT%d",minEcell,minEseed, dTime,wTime));
220 printf("Name of clusterizer array: %s\n",arrayNameV2.Data());
221
193828fd 222 if(!kMC)
223 {
cd54ca47 224
225 AliAnalysisTaskParticleCorrelation *anav2tr = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
226 kYear,kRun,kCollision,"EMC7",arrayNameV2);
54769bc0 227
cd54ca47 228 AliAnalysisTaskParticleCorrelation *anav2mb = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
229 kYear,kRun,kCollision,"INT7",arrayNameV2);
193828fd 230 }
231 else
cd54ca47 232 {// No trigger (should be MB, but for single particle productions it does not work)
233 AliAnalysisTaskParticleCorrelation *anav2 = AddTaskPartCorr(kInputData, "EMCAL", kPrint,kMC, deltaAOD, outputFile.Data(),
234 kYear,kRun,kCollision,"",arrayNameV2);
d92b41ad 235 }
cd54ca47 236
237 // -----
238 // PHOS
239 // -----
240 /*
241 //Add here PHOS tender or whatever is needed
242
243 if(!kMC)
244 {
245
246 AliAnalysisTaskParticleCorrelation *anav1tr = AddTaskPartCorr(kInputData, "PHOS", kPrint,kMC, deltaAOD, outputFile.Data(),
247 kYear,kRun,kCollision,"PHOS",""); // Change to PHOS trigger
248
249
250 AliAnalysisTaskParticleCorrelation *anav1mb = AddTaskPartCorr(kInputData, "PHOS", kPrint,kMC, deltaAOD, outputFile.Data(),
251 kYear,kRun,kCollision,"INT7","");
252
253 }
254 else
255 {// No trigger
256
257 AliAnalysisTaskParticleCorrelation *anav1mb = AddTaskPartCorr(kInputData, "PHOS", kPrint,kMC, deltaAOD, outputFile.Data(),
258 kYear,kRun,kCollision,"","");
259
260 }
261 */
193828fd 262
263
264 //-----------------------
265 // Run the analysis
266 //-----------------------
267 mgr->InitAnalysis();
268 mgr->PrintStatus();
269 mgr->StartAnalysis("local",chain);
270
271 cout <<" Analysis ended sucessfully "<< endl ;
272
d92b41ad 273
274}
275
193828fd 276//___________________
277void LoadLibraries()
278{
d92b41ad 279
280 //--------------------------------------
281 // Load the needed libraries most of them already loaded by aliroot
282 //--------------------------------------
283 gSystem->Load("libTree.so");
284 gSystem->Load("libGeom.so");
285 gSystem->Load("libVMC.so");
286 gSystem->Load("libXMLIO.so");
54769bc0 287 gSystem->Load("libMatrix.so");
288 gSystem->Load("libPhysics.so");
193828fd 289 gSystem->Load("libMinuit.so"); // Root + libraries to if reclusterization is done
290
291 gSystem->Load("libSTEERBase.so");
292 gSystem->Load("libGui.so"); // Root + libraries to if reclusterization is done
293 gSystem->Load("libCDB.so"); // Root + libraries to if reclusterization is done
294 gSystem->Load("libESD.so"); // Root + libraries to if reclusterization is done
295 gSystem->Load("libAOD.so");
296 gSystem->Load("libRAWDatabase.so"); // Root + libraries to if reclusterization is done
297 gSystem->Load("libProof.so");
298 gSystem->Load("libANALYSIS.so");
299 gSystem->Load("libSTEER.so"); // Root + libraries to if reclusterization is done
300
301 gSystem->Load("libRAWDatarec.so"); // Root + libraries to if reclusterization is done
302 gSystem->Load("libRAWDatasim.so"); // Root + libraries to if reclusterization is done
303 gSystem->Load("libVZERObase.so"); // Root + libraries to if reclusterization is done
304 gSystem->Load("libVZEROrec.so"); // Root + libraries to if reclusterization is done
305
306 gSystem->Load("libEMCALUtils");
307 //SetupPar("EMCALUtils");
308 gSystem->Load("libEMCALraw"); // Root + libraries to if reclusterization is done
309 gSystem->Load("libEMCALbase"); // Root + libraries to if reclusterization is done
310 gSystem->Load("libEMCALsim"); // Root + libraries to if reclusterization is done
311 gSystem->Load("libEMCALrec"); // Root + libraries to if reclusterization is done
312 //SetupPar("EMCALraw");
313 //SetupPar("EMCALbase");
314 //SetupPar("EMCALsim");
315 //SetupPar("EMCALrec");
316
317 gSystem->Load("libANALYSISalice.so");
318 //gSystem->Load("libTENDER.so");
319 //gSystem->Load("libTENDERSupplies.so");
320
321 gSystem->Load("libPHOSUtils");
322 gSystem->Load("libEMCALUtils");
323 gSystem->Load("libPWG4PartCorrBase");
324 gSystem->Load("libPWG4PartCorrDep");
325 gSystem->Load("libPWG4CaloCalib");
326 //SetupPar("PWG4PartCorrBase");
327 //SetupPar("PWG4PartCorrDep");
328 //SetupPar("PWG4CaloCalib");
d92b41ad 329}
330
331void SetupPar(char* pararchivename)
332{
333 //Load par files, create analysis libraries
334 //For testing, if par file already decompressed and modified
335 //classes then do not decompress.
193828fd 336
d92b41ad 337 TString cdir(Form("%s", gSystem->WorkingDirectory() )) ;
338 TString parpar(Form("%s.par", pararchivename)) ;
193828fd 339
d92b41ad 340 if ( gSystem->AccessPathName(pararchivename) ) {
341 TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
342 gROOT->ProcessLine(processline.Data());
343 }
344
345 TString ocwd = gSystem->WorkingDirectory();
346 gSystem->ChangeDirectory(pararchivename);
347
348 // check for BUILD.sh and execute
349 if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
350 printf("*******************************\n");
351 printf("*** Building PAR archive ***\n");
352 cout<<pararchivename<<endl;
353 printf("*******************************\n");
354
355 if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
356 Error("runProcess","Cannot Build the PAR Archive! - Abort!");
357 return -1;
358 }
359 }
360 // check for SETUP.C and execute
361 if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
362 printf("*******************************\n");
363 printf("*** Setup PAR archive ***\n");
364 cout<<pararchivename<<endl;
365 printf("*******************************\n");
366 gROOT->Macro("PROOF-INF/SETUP.C");
367 }
368
369 gSystem->ChangeDirectory(ocwd.Data());
370 printf("Current dir: %s\n", ocwd.Data());
371}
372
193828fd 373void CheckInputData(const anaModes mode){
374 //Sets input data and tree
375
376 TString ocwd = gSystem->WorkingDirectory();
377
378 //---------------------------------------
379 //Local files analysis
380 //---------------------------------------
381 if(mode == mLocal){
382 //If you want to add several ESD files sitting in a common directory INDIR
383 //Specify as environmental variables the directory (INDIR), the number of files
384 //to analyze (NFILES) and the pattern name of the directories with files (PATTERN)
cd54ca47 385
193828fd 386 if(gSystem->Getenv("INDIR"))
387 kInDir = gSystem->Getenv("INDIR") ;
388 else cout<<"INDIR not set, use default: "<<kInDir<<endl;
389
390 TString sindir(kInDir);
391 if (sindir.Contains("pass1")) kPass = "pass1";
392 else if(sindir.Contains("pass2")) kPass = "pass2";
393 else if(sindir.Contains("pass3")) kPass = "pass3";
394
395 if(gSystem->Getenv("PATTERN"))
396 kPattern = gSystem->Getenv("PATTERN") ;
397 else cout<<"PATTERN not set, use default: "<<kPattern<<endl;
398
399 cout<<"INDIR : "<<kInDir<<endl;
400 cout<<"NFILES : "<<kFile<<endl;
401
402 char fileE[120] ;
403 char fileA[120] ;
404 char fileG[120] ;
405 char fileEm[120] ;
406 for (Int_t event = 0 ; event < kFile ; event++) {
407 sprintf(fileE, "%s/%s%d/AliESDs.root", kInDir,kPattern,event) ;
408 sprintf(fileA, "%s/%s%d/AliAOD.root", kInDir,kPattern,event) ;
409 sprintf(fileG, "%s/%s%d/galice.root", kInDir,kPattern,event) ;
410 sprintf(fileEm, "%s/%s%d/embededAOD.root", kInDir,kPattern,event) ;
411
412 TFile * fAOD = 0 ;
413 //Check if file exists and add it, if not skip it
414 if ( TFile::Open(fileE)) {
415 kTreeName ="esdTree";
416 kInputData = "ESD";
417 if(TFile::Open(fileG)) kMC=kTRUE;
418 else kMC = kFALSE;
419 return;
420 }
421 else if(fAOD=TFile::Open(fileA)){
422 kTreeName ="aodTree";
423 kInputData = "AOD";
424 if(((TTree*) fAOD->Get("aodTree"))->GetBranch("mcparticles")) kMC=kTRUE;
425 else kMC = kFALSE;
426 return;
427 }
428 else if(fAOD=TFile::Open(fileEm)){
429 kTreeName ="aodTree";
430 kInputData = "AOD";
431 kMC=kTRUE;
432 return;
433 }
434 else if(TFile::Open(fileG)){
435 kTreeName ="TE";
436 kInputData = "MC";
437 kMC=kTRUE;
438 return;
439 }
440 }
441
442 }// local files analysis
443
444 //------------------------------
445 //GRID xml files
446 //-----------------------------
447 else if(mode == mGRID){
448 //Get colection file. It is specified by the environmental
449 //variable XML
450
451 if(gSystem->Getenv("XML") )
452 kXML = gSystem->Getenv("XML");
453 else
454 sprintf(kXML, "collection.xml") ;
455
456 if (!TFile::Open(kXML)) {
457 printf("No collection file with name -- %s -- was found\n",kXML);
458 return ;
459 }
460 else cout<<"XML file "<<kXML<<endl;
461
462 //Load necessary libraries and connect to the GRID
463 gSystem->Load("libNetx.so") ;
464 gSystem->Load("libRAliEn.so");
465 TGrid::Connect("alien://") ;
466
467 //Feed Grid with collection file
468 TGridCollection * collection = (TGridCollection*) TAlienCollection::Open(kXML);
469 if (! collection) {
470 AliError(Form("%s not found", kXML)) ;
471 return kFALSE ;
472 }
473 TGridResult* result = collection->GetGridResult("",0 ,0);
474
475 for (Int_t index = 0; index < result->GetEntries(); index++) {
476 TString alienURL = result->GetKey(index, "turl") ;
477 cout << "================== " << alienURL << endl ;
478
479 if (alienURL.Contains("pass1")) kPass = "pass1";
480 else if(alienURL.Contains("pass2")) kPass = "pass2";
481 else if(alienURL.Contains("pass3")) kPass = "pass3";
482
483 TFile * fAOD = 0 ;
484 //Check if file exists and add it, if not skip it
485 if (alienURL.Contains("AliESDs.root")) {
486 kTreeName ="esdTree";
487 kInputData = "ESD";
488 alienURL.ReplaceAll("AliESDs.root","galice.root");
489 if(TFile::Open(alienURL)) kMC=kTRUE;
490 else kMC = kFALSE;
491
492 kRun = AliAnalysisManager::GetRunFromAlienPath(alienURL.Data());
493 printf("Run number from alien path = %d\n",kRun);
494
495 return;
496 }
497 else if(alienURL.Contains("AliAOD.root")){
498 kTreeName ="aodTree";
499 kInputData = "AOD";
500 fAOD = TFile::Open(alienURL);
501 if(((TTree*) fAOD->Get("aodTree"))->GetBranch("mcparticles")) kMC=kTRUE;
502 else kMC = kFALSE;
503 return;
504 }
505 else if(alienURL.Contains("embededAOD.root")){
506 kTreeName ="aodTree";
507 kInputData = "AOD";
508 kMC=kTRUE;
509 return;
510 }
511 else if(alienURL.Contains("galice.root")){
512 kTreeName ="TE";
513 kInputData = "MC";
514 kMC=kTRUE;
515 return;
516 }
517 }
518 }// xml analysis
519
520 gSystem->ChangeDirectory(ocwd.Data());
521
522}
d92b41ad 523
193828fd 524//_____________________________________________________________________
525void CreateChain(const anaModes mode, TChain * chain, TChain * chainxs)
526{
d92b41ad 527 //Fills chain with data
528 TString ocwd = gSystem->WorkingDirectory();
529
d92b41ad 530 //---------------------------------------
531 //Local files analysis
532 //---------------------------------------
2244659d 533 if(mode == mLocal){
d92b41ad 534 //If you want to add several ESD files sitting in a common directory INDIR
535 //Specify as environmental variables the directory (INDIR), the number of files
c8fe2783 536 //to analyze (NFILES) and the pattern name of the directories with files (PATTERN)
193828fd 537
d92b41ad 538 if(gSystem->Getenv("INDIR"))
539 kInDir = gSystem->Getenv("INDIR") ;
540 else cout<<"INDIR not set, use default: "<<kInDir<<endl;
541
542 if(gSystem->Getenv("PATTERN"))
543 kPattern = gSystem->Getenv("PATTERN") ;
544 else cout<<"PATTERN not set, use default: "<<kPattern<<endl;
545
c8fe2783 546 if(gSystem->Getenv("NFILES"))
547 kFile = atoi(gSystem->Getenv("NFILES")) ;
548 else cout<<"NFILES not set, use default: "<<kFile<<endl;
d92b41ad 549
550 //Check if env variables are set and are correct
c8fe2783 551 if ( kInDir && kFile) {
552 printf("Get %d files from directory %s\n",kFile,kInDir);
d92b41ad 553 if ( ! gSystem->cd(kInDir) ) {//check if ESDs directory exist
193828fd 554 printf("%s does not exist\n", kInDir) ;
555 return ;
d92b41ad 556 }
193828fd 557
c8fe2783 558 //if(gSystem->Getenv("XSFILE"))
559 //kXSFileName = gSystem->Getenv("XSFILE") ;
560 //else cout<<" XS file name not set, use default: "<<kXSFileName<<endl;
561 char * kGener = gSystem->Getenv("GENER");
562 if(kGener) {
193828fd 563 cout<<"GENER "<<kGener<<endl;
564 if(!strcmp(kGener,"PYTHIA")) kXSFileName = "pyxsec.root";
565 else if(!strcmp(kGener,"HERWIG")) kXSFileName = "hexsec.root";
566 else cout<<" UNKNOWN GENER, use default: "<<kXSFileName<<endl;
c8fe2783 567 }
568 else cout<<" GENER not set, use default xs file name: "<<kXSFileName<<endl;
193828fd 569
c8fe2783 570 cout<<"INDIR : "<<kInDir<<endl;
571 cout<<"NFILES : "<<kFile<<endl;
572 cout<<"PATTERN : " <<kPattern<<endl;
573 cout<<"XSFILE : "<<kXSFileName<<endl;
193828fd 574
7175a03a 575 TString datafile="";
576 if(kInputData == "ESD") datafile = "AliESDs.root" ;
193828fd 577 else if(kInputData.Contains("AOD")) datafile = "AliAOD.root" ;//////////
7175a03a 578 else if(kInputData == "MC") datafile = "galice.root" ;
579
d92b41ad 580 //Loop on ESD files, add them to chain
581 Int_t event =0;
582 Int_t skipped=0 ;
583 char file[120] ;
584 char filexs[120] ;
585
c8fe2783 586 for (event = 0 ; event < kFile ; event++) {
193828fd 587 sprintf(file, "%s/%s%d/%s", kInDir,kPattern,event,datafile.Data()) ;
588 sprintf(filexs, "%s/%s%d/%s", kInDir,kPattern,event,kXSFileName) ;
589 TFile * fESD = 0 ;
590 //Check if file exists and add it, if not skip it
591 if ( fESD = TFile::Open(file)) {
592 if ( fESD->Get(kTreeName) ) {
593 printf("++++ Adding %s\n", file) ;
594 chain->AddFile(file);
595 chainxs->Add(filexs) ;
596 }
597 }
598 else {
599 printf("---- Skipping %s\n", file) ;
600 skipped++ ;
601 }
d92b41ad 602 }
d92b41ad 603 }
604 else {
605 TString input = "AliESDs.root" ;
606 cout<<">>>>>> No list added, take a single file <<<<<<<<< "<<input<<endl;
607 chain->AddFile(input);
608 }
609
610 }// local files analysis
611
612 //------------------------------
613 //GRID xml files
614 //-----------------------------
615 else if(mode == mGRID){
616 //Get colection file. It is specified by the environmental
617 //variable XML
d92b41ad 618
d92b41ad 619 //Feed Grid with collection file
d92b41ad 620 TGridCollection * collection = (TGridCollection*) TAlienCollection::Open(kXML);
621 if (! collection) {
622 AliError(Form("%s not found", kXML)) ;
623 return kFALSE ;
624 }
625 TGridResult* result = collection->GetGridResult("",0 ,0);
193828fd 626
d92b41ad 627 // Makes the ESD chain
628 printf("*** Getting the Chain ***\n");
629 for (Int_t index = 0; index < result->GetEntries(); index++) {
630 TString alienURL = result->GetKey(index, "turl") ;
631 cout << "================== " << alienURL << endl ;
632 chain->Add(alienURL) ;
633 alienURL.ReplaceAll("AliESDs.root",kXSFileName);
634 chainxs->Add(alienURL) ;
635 }
636 }// xml analysis
637
638 gSystem->ChangeDirectory(ocwd.Data());
193828fd 639
d92b41ad 640}
641
193828fd 642//_________________________________________________________________
d92b41ad 643void GetAverageXsection(TTree * tree, Double_t & xs, Float_t & ntr)
644{
645 // Read the PYTHIA statistics from the file pyxsec.root created by
646 // the function WriteXsection():
647 // integrated cross section (xsection) and
648 // the number of Pyevent() calls (ntrials)
649 // and calculate the weight per one event xsection/ntrials
650 // The spectrum calculated by a user should be
651 // multiplied by this weight, something like this:
652 // TH1F *userSpectrum ... // book and fill the spectrum
653 // userSpectrum->Scale(weight)
654 //
655 // Yuri Kharlov 19 June 2007
656 // Gustavo Conesa 15 April 2008
657 Double_t xsection = 0;
658 UInt_t ntrials = 0;
659 xs = 0;
660 ntr = 0;
661
662 Int_t nfiles = tree->GetEntries() ;
663 if (tree && nfiles > 0) {
664 tree->SetBranchAddress("xsection",&xsection);
665 tree->SetBranchAddress("ntrials",&ntrials);
666 for(Int_t i = 0; i < nfiles; i++){
667 tree->GetEntry(i);
668 xs += xsection ;
669 ntr += ntrials ;
670 cout << "xsection " <<xsection<<" ntrials "<<ntrials<<endl;
671 }
672
673 xs = xs / nfiles;
674 ntr = ntr / nfiles;
675 cout << "-----------------------------------------------------------------"<<endl;
676 cout << "Average of "<< nfiles<<" files: xsection " <<xs<<" ntrials "<<ntr<<endl;
677 cout << "-----------------------------------------------------------------"<<endl;
678 }
679 else cout << " >>>> Empty tree !!!! <<<<< "<<endl;
680
681}
682
683
193828fd 684//______________________________
685void CheckEnvironmentVariables()
686{
687
cd54ca47 688 sprintf(kTrigger,"");
193828fd 689
cd54ca47 690 Bool_t bRecalibrate = kFALSE;
691 Bool_t bBadChannel = kFALSE;
193828fd 692
cd54ca47 693 for (int i=0; i< gApplication->Argc();i++){
694
695#ifdef VERBOSEARGS
696
697 printf("Arg %d: %s\n",i,gApplication->Argv(i));
698
193828fd 699#endif
cd54ca47 700
701 if (!(strcmp(gApplication->Argv(i),"--trigger")))
702 sprintf(trigger,gApplication->Argv(i+1));
193828fd 703
704 if (!(strcmp(gApplication->Argv(i),"--recalibrate")))
705 bRecalibrate = atoi(gApplication->Argv(i+1));
cd54ca47 706
707 if (!(strcmp(gApplication->Argv(i),"--badchannel")))
708 bBadChannel = atoi(gApplication->Argv(i+1));
709
710 if (!(strcmp(gApplication->Argv(i),"--run"))){
711 TString sRun(gApplication->Argv(i+1));
712 if(sRun.Contains("LHC10")) {
713 kYear = 2010;
714 }
715 else {
716 if(kRun <=0){
717 kRun = atoi(gApplication->Argv(i+1));
718 }
719 else printf("** Run number already set to %d, do not set to %d\n",kRun,atoi(gApplication->Argv(i+1)));
720 }//numeric run
721 }//--run available
722
723 }// args loop
724
193828fd 725 if(!sRun.Contains("LHC10")){
726 if ( kRun < 140000) {
727 kYear = 2010;
728 if( kRun >= 136851 ) kCollision = "PbPb";
729 }
730 else{
731 kYear = 2011;
732 }
733 }
734
cd54ca47 735 if(kMC) sprintf(kTrigger,"");
736
737 printf("*********************************************\n");
738 //printf("*** Settings trigger %s, recalibrate %d, remove bad channels %d, year %d, collision %s, run %d ***\n",
739 // kTrigger,bRecalibrate,bBadChannel, kYear,kCollision.Data(), kRun);
740 printf("*** Settings year %d, collision %s, run %d ***\n",kYear,kCollision.Data(), kRun);
741 printf("*********************************************\n");
742
193828fd 743}
744
745//_____________________________________________________________
746void AddTaskCounter(const TString trigger = "MB")
747{
748
749 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
750
751 AliAnalysisTaskCounter * counter = new AliAnalysisTaskCounter(Form("Counter%s",trigger.Data()));
752 if(kRun > 140000 && kRun < 146900) counter ->RejectFastCluster();
753 if (kCollision=="pp" ) counter->SetZVertexCut(50.); //Open cut
754 else if(kCollision=="PbPb") counter->SetZVertexCut(10.); //Centrality defined in this range.
755
756 if(trigger=="EMC7"){
757 printf("counter trigger EMC7\n");
758 counter->SelectCollisionCandidates(AliVEvent::kEMC7);
759 }
760 else if (trigger=="INT7"){
761 printf("counter trigger INT7\n");
762 counter->SelectCollisionCandidates(AliVEvent::kINT7);
763 }
764 if(trigger=="EMC1"){
765 printf("counter trigger EMC1\n");
766 counter->SelectCollisionCandidates(AliVEvent::kEMC1);
767 }
768 else if(trigger=="MB"){
769 printf("counter trigger MB\n");
770 counter->SelectCollisionCandidates(AliVEvent::kMB);
771 }
cd54ca47 772 else if(trigger=="PHOS"){
773 printf("counter trigger PHOS\n");
774 counter->SelectCollisionCandidates(AliVEvent::kPHI7);
775 }
193828fd 776
777 TString outputFile = AliAnalysisManager::GetCommonFileName();
778 AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer();
779
780 AliAnalysisDataContainer *coutput =
781 mgr->CreateContainer(Form("Counter%s",trigger.Data()), TList::Class(), AliAnalysisManager::kOutputContainer, outputFile.Data());
782 mgr->AddTask(counter);
783 mgr->ConnectInput (counter, 0, cinput1);
784 mgr->ConnectOutput (counter, 1, coutput);
cd54ca47 785
193828fd 786}
787
788
789
790
d92b41ad 791