Compatibility with the Root trunk
[u/mrichter/AliRoot.git] / PWGJE / UserTasks / AliAnalysisTaskCheckSingleTrackJetRejection.cxx
CommitLineData
635693ee 1//#include <string.h>
2//#include <TStyle.h>
3#include <list>
4#include <string>
5
6#include "TTree.h"
635693ee 7#include "TCanvas.h"
8#include "AliAnalysisTask.h"
9#include "AliInputEventHandler.h"
10#include "AliESDtrack.h"
11#include "AliAODVertex.h"
12#include "AliAODCluster.h"
13
14#include <TROOT.h>
15#include <TRandom.h>
16#include <TSystem.h>
17#include <TInterpreter.h>
18#include <TChain.h>
19#include <TFile.h>
20#include <TKey.h>
21#include <TH1F.h>
22#include <TH2F.h>
23#include <TH3F.h>
24#include <TProfile.h>
25#include <TList.h>
26#include <TLorentzVector.h>
27#include <TClonesArray.h>
28#include <TRefArray.h>
29
30#include "TDatabasePDG.h"
31#include "AliAnalysisManager.h"
32#include "AliJetFinder.h"
33#include "AliJetHeader.h"
34#include "AliJetReader.h"
35#include "AliJetReaderHeader.h"
36#include "AliUA1JetHeaderV1.h"
37#include "AliSISConeJetHeader.h"
38#include "AliESDEvent.h"
39#include "AliAODEvent.h"
40#include "AliAODHandler.h"
1d478650 41#include "AliAODInputHandler.h"
635693ee 42#include "AliAODTrack.h"
43#include "AliVParticle.h"
44#include "AliAODJet.h"
45#include "AliAODJetEventBackground.h"
46#include "AliMCParticle.h"
47#include "AliAODMCParticle.h"
48#include "AliMCEventHandler.h"
49#include "AliMCEvent.h"
50#include "AliStack.h"
51
52#include "AliAODHeader.h"
53#include "AliAODMCHeader.h"
54#include "AliJetKineReaderHeader.h"
55#include "AliGenCocktailEventHeader.h"
56#include "AliInputEventHandler.h"
57#include "AliGenEventHeader.h"
58#include "AliGenDPMjetEventHeader.h"
59
60
61
62#include "AliAnalysisTaskCheckSingleTrackJetRejection.h"
63#include "AliAnalysisTaskPhiCorrelations.h"
64#include "AliAnalysisHelperJetTasks.h"
65#include "AliPWG4HighPtQAMC.h"
66
67
68
69
c64cb1f6 70using std::cout;
71using std::endl;
72
635693ee 73ClassImp(AliAnalysisTaskCheckSingleTrackJetRejection)
74
75 //________________________________________________________________________
76 AliAnalysisTaskCheckSingleTrackJetRejection::AliAnalysisTaskCheckSingleTrackJetRejection():
77 AliAnalysisTaskSE(),
78 fUseAODInput(kFALSE),
79 fFillAOD(kFALSE),
80 fNonStdFile(""),
81 fJetBranch("jets"),
82 fAODIn(0x0),
83 fAODOut(0x0),
84 fAODExtension(0x0),
1ae43f37 85 JFAlg("ANTIKT"),
86 Radius(0.4),
87 Filtermask(256),
88 BackM(0),
89 TrackPtcut(0.15),
90 SkipCone(0),
91 IsMC(kTRUE),
635693ee 92 fHistList(0x0),
93 fxsec(0.),
94 ftrial(1.),
1ae43f37 95 fJetRecEtaWindow(0.5),
96 fMinJetPt(0.15),
635693ee 97 fH1Xsec(0x0),
1ae43f37 98 fH1Trials(0x0),
99 fH1Events(0x0),
100
101 fH2jetMCAKT04_Jetpt_maxpt(0x0),
102 fH2jetAKT04_Jetpt_maxpt (0x0)
635693ee 103
1ae43f37 104{
635693ee 105 for(int j=0;j<6;j++){
106 fH1jetMCAKT04_pt [j]=0;
107 fH1jetAKT04_pt [j]=0;
108
109 fH2jetMCAKT04_Eratio [j]=0;
110 fH1jetMCAKT04_match [j]=0;
111 }
112
113 // Default constructor
114}
115
116//________________________________________________________________________
117AliAnalysisTaskCheckSingleTrackJetRejection::AliAnalysisTaskCheckSingleTrackJetRejection(const char *name):
118 AliAnalysisTaskSE(name),
119 fUseAODInput(kFALSE),
120 fFillAOD(kFALSE),
121 fNonStdFile(""),
122 fJetBranch("jets"),
123 fAODIn(0x0),
124 fAODOut(0x0),
125 fAODExtension(0x0),
1ae43f37 126 JFAlg("ANTIKT"),
127 Radius(0.4),
128 Filtermask(256),
129 BackM(0),
130 TrackPtcut(0.15),
131 SkipCone(0),
132 IsMC(kTRUE),
635693ee 133 fHistList(0x0),
134 fxsec(0.),
135 ftrial(1.),
1ae43f37 136 fJetRecEtaWindow(0.5),
137 fMinJetPt(0.15),
635693ee 138 fH1Xsec(0x0),
1ae43f37 139 fH1Trials(0x0),
140 fH1Events(0x0),
141 fH2jetMCAKT04_Jetpt_maxpt(0x0),
142 fH2jetAKT04_Jetpt_maxpt (0x0)
635693ee 143
144{
145
635693ee 146 for(int j=0;j<6;j++){
147 fH1jetMCAKT04_pt [j]=0;
148 fH1jetAKT04_pt [j]=0;
149
150 fH2jetMCAKT04_Eratio [j]=0;
151 fH1jetMCAKT04_match [j]=0;
152 }
153
154 // Constructor
155 // Define input and output slots here
156 // Input slot #0 works with a TChain
157 DefineInput(0, TChain::Class());
158 // Output slot #0 id reserved by the base class for AOD
159 // Output slot #1 writes into a TH1 container
160 DefineOutput(1, TList::Class());
161
162}
163
164//________________________________________________________________________
165void AliAnalysisTaskCheckSingleTrackJetRejection::UserCreateOutputObjects()
166{
167 // Create histograms
168 // Called once
169
170 if (!fHistList){ fHistList = new TList();fHistList->SetOwner(kTRUE); cout<<"TList is created for output "<<endl;}
171
172 Bool_t oldStatus = TH1::AddDirectoryStatus();
173 TH1::AddDirectory(kFALSE);
174
175 char *histname;
176 if(IsMC){
177 fH1Xsec = new TProfile("Xsec","Xsec",1,0,1);
178 fH1Trials = new TH1F ("Trials","Trials",1,0,1);
179 histname = Form("fH2jetMCAKT04_Jetpt_maxpt");
180 fH2jetMCAKT04_Jetpt_maxpt = new TH2F(histname,histname,400,0,400,400,0,400);
181 histname = Form("fH2jetAKT04_Jetpt_maxpt");
182 fH2jetAKT04_Jetpt_maxpt = new TH2F(histname,histname,400,0,400,400,0,400);
183 for(int j=0;j<6;j++){
184 histname = Form("fH1jetMCAKT04_pt%d",j);
185 fH1jetMCAKT04_pt[j] = new TH1F(histname,histname,200,0,200);
186 histname = Form("fH2jetAKT04_pt%d",j);
187 fH1jetAKT04_pt[j] = new TH1F(histname,histname,200,0,200);
188
189 histname = Form("fH2jetMCAKT04_Eratio%d",j);
190 fH2jetMCAKT04_Eratio[j] = new TH2F(histname,histname,200,0,200,100,0,30);
191 histname = Form("fH1jetMCAKT04_match%d",j);
192 fH1jetMCAKT04_match[j] = new TH1F(histname,histname,200,0,200);
193 }
194 fHistList->Add(fH1Xsec);
195 fHistList->Add(fH1Trials);
196 fHistList->Add(fH2jetMCAKT04_Jetpt_maxpt);
197 fHistList->Add(fH2jetAKT04_Jetpt_maxpt);
198 for(int j=0;j<6;j++){
199 fHistList->Add(fH1jetAKT04_pt[j]);
200 fHistList->Add(fH1jetMCAKT04_pt[j]);
201 fHistList->Add(fH2jetMCAKT04_Eratio[j]);
202 fHistList->Add(fH1jetMCAKT04_match[j]);
203 }
204 }
205 else {
206 fH1Events = new TH1F("Events","Number of Events",1,0,1);
207 histname = Form("fH2jetAKT04_Jetpt_maxpt");
208 fH2jetAKT04_Jetpt_maxpt = new TH2F(histname,histname,400,0,400,400,0,400);
209 for(int j=0;j<6;j++){
210 histname = Form("fH2jetAKT04_pt%d",j);
211 fH1jetAKT04_pt[j] = new TH1F(histname,histname,200,0,200);
212 }
213 fHistList->Add(fH1Events);
214 fHistList->Add(fH2jetAKT04_Jetpt_maxpt);
215 for(int j=0;j<6;j++){
216 fHistList->Add(fH1jetAKT04_pt[j]);
217 }
218 }
219
220
221 // =========== Switch on Sumw2 for all histos ===========
222 for (Int_t i=0; i<fHistList->GetEntries(); ++i)
223 {
224 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
225 if (h1)
226 {
227 h1->Sumw2();
228 continue;
229 }
230 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
231 if(hn)hn->Sumw2();
232 }
233 TH1::AddDirectory(oldStatus);
234
235
236 PostData(1,fHistList);
237
238}
239
240//----------------------------------------------------------------------
241void AliAnalysisTaskCheckSingleTrackJetRejection::Init()
242{
243 // Initialization
244 if (fDebug) printf("AnalysisTaskCheckSingleTrackJetRejection::Init() \n");
245
246}
247
248Bool_t AliAnalysisTaskCheckSingleTrackJetRejection::Notify()
249{
250 fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
251 fAODOut = AODEvent();
252 if(fNonStdFile.Length()!=0){
253 AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
254 fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
255 if(fAODExtension){
256 if(fDebug>1)Printf("AODExtension found for %s ",fNonStdFile.Data());
257 }
258 }
259
260 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
261 fxsec=0.;
262 ftrial=1.;
263
264 if(tree){
265 TFile *curfile = tree->GetCurrentFile();
266 if(!curfile){
267 Error("Notify","No current file");
268 return kFALSE;
269 }
270
271 if(IsMC){
272 PythiaInfoFromFile(curfile->GetName(),fxsec,ftrial);
273 //cout<<" Xsec "<<fxsec<<" trial "<<ftrial<<endl;
274 fH1Xsec ->Fill(0.,fxsec);
275 fH1Trials->Fill(0.,ftrial);
276 }
277 else{
f4b5b842 278 //Float_t totalEvent;
279 //totalEvent = GetTotalEvents(curfile->GetName());
280 //fH1Events->Fill(0.,totalEvent);
635693ee 281 }
282
283 }
284
285 printf("Reading File %s ",fInputHandler->GetTree()->GetCurrentFile()->GetName());
286 return kTRUE;
287}
288
289void AliAnalysisTaskCheckSingleTrackJetRejection::FinishTaskOutput()
290{
291}
292
293
294
295//________________________________________________________________________
296void AliAnalysisTaskCheckSingleTrackJetRejection::UserExec(Option_t *)
297{
298
299
300 // Main loop (called each event)
301 // Execute analysis for current event
1d478650 302
303 AliAODEvent *fAOD;
304 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
305 if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
306 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
307 if(fUseAODInput) fAODIn = fAOD;
308 if (fDebug > 1) Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
309 }
310 else {
311 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
312 if( handler && handler->InheritsFrom("AliAODHandler") ) {
313 fAOD = ((AliAODHandler*)handler)->GetAOD();
314 fAODIn = fAOD;
315 if (fDebug > 1) Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
316 }
317 }
318
319 if(!fAODIn && !fUseAODInput){ // case we have AOD in input & output and want jets from output
320 TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
321 if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
322 fAODIn = ((AliAODHandler*)outHandler)->GetAOD();
323 if (fDebug > 1) Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
324 }
325 }
326
28c36b48 327 //if(fNonStdFile.Length()!=0){
328 // // case we have an AOD extension - fetch the jets from the extended output
329 // AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
330 // fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
331 // if(!fAODExtension){
332 // if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
333 // }
334 //}
1d478650 335 //fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
336 if (!fAODIn) {
337 Printf("ERROR %s : fAODIn not available",(char*)__FILE__);
338 return;
339 }
340
f4b5b842 341 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
342 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
343 if(!(inputHandler->IsEventSelected() & AliVEvent::kMB)){
344 if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
345 return;
346 }
347 if(!IsMC)fH1Events->Fill(0);
348
635693ee 349
1d478650 350
635693ee 351 // start jet analysis
352
353 Double_t Jet_n [20];
354 Double_t Jet_pt [20][10000];
355 Double_t Jet_eta[20][10000];
356 Double_t Jet_phi[20][10000];
357
358 for(int i=0;i<20;i++){
359 Jet_n[i]=0;
360 for(int j=0;j<10000;j++){
361 Jet_pt[i][j]=0.;
362 Jet_phi[i][j]=999.;
363 Jet_eta[i][j]=999.;
364 }
365 }
366
1d478650 367 // get AOD event from input/ouput
635693ee 368 TString cAdd = "";
369 cAdd += Form("%02d_",(int)((Radius+0.01)*10.));
370 cAdd += Form("B%d",(int)BackM);
371 cAdd += Form("_Filter%05d",Filtermask);
372 cAdd += Form("_Cut%05d",(int)(1000.*TrackPtcut));
373 cAdd += Form("_Skip%02d",SkipCone);
374 TString Branchname_gen,Branchname_rec;
652845bc 375 //Branchname_gen = Form("clustersMCKINE2_%s%s",JFAlg.Data(),cAdd.Data());
376 Branchname_gen = Form("clustersAODMC2_%s%s",JFAlg.Data(),cAdd.Data());
635693ee 377 Branchname_rec = Form("clustersAOD_%s%s",JFAlg.Data(),cAdd.Data());
378
379
380 bool fFIND[6][1000];
381 double maxpt[2][1000];for(int i=0;i<1000;i++){maxpt[0][i]=0;maxpt[1][i]=0;}
382 int nearrecID[1000]; for(int i=0;i<1000;i++){nearrecID[i]=99999;}
383 AliAODJet* jetsAOD;
384 for(int algorithm=0;algorithm<2;algorithm++){
385 //for LHC11a1 LHC11a2 official
386 if((!IsMC&&algorithm==0))continue;
387
388 if(algorithm==0)fJetBranch = Branchname_gen.Data();
389 if(algorithm==1)fJetBranch = Branchname_rec.Data();
390
391 TClonesArray* jets = dynamic_cast <TClonesArray*> (fAODIn->FindListObject(fJetBranch.Data()));
652845bc 392 if(!jets){
393 printf(" Tere are no Branch named %s \n",fJetBranch.Data());
394 continue;
395 }
635693ee 396 Int_t nj = jets->GetEntriesFast();
397 if (fDebug) printf("There are %5d jets in the event \n", nj);
398
399
400 Jet_n[algorithm] = nj;
401 for(int njet =0;njet<nj;njet++){
402 jetsAOD = (AliAODJet*) (jets->At(njet));
403 Jet_pt [algorithm][njet] = jetsAOD->Pt();
404 Jet_phi [algorithm][njet] = jetsAOD->Phi();
405 Jet_eta [algorithm][njet] = jetsAOD->Eta();
406 double eta_cut_Jet=0.5;
407 TRefArray *reftracks = jetsAOD->GetRefTracks();
408 int ntracks=reftracks->GetEntriesFast();
409 if(TMath::Abs(Jet_eta[algorithm][njet])<eta_cut_Jet){
410 //------------calc max pt in Jet----------------
411 double trackpt;
412 double sumtrackpt=0;//test
413 for(int ntr=0;ntr<ntracks;ntr++){// calc. max pt of track which is in Jet
414 AliAODTrack *AODtrack = dynamic_cast<AliAODTrack*>(reftracks->At(ntr));
415 if(AODtrack){
00a3f17b 416 Bool_t bgoodT=false;
417 if(Filtermask!=768){
418 if(AODtrack->TestFilterMask(Filtermask))bgoodT=true;
419 }
420 else{
421 if(AODtrack->IsHybridGlobalConstrainedGlobal())bgoodT=true; //for hybrid Track cuts
422 }
423 if(!bgoodT)continue;
424 trackpt = AODtrack->Pt();
425 sumtrackpt += trackpt;
426 if(trackpt>maxpt[algorithm][njet]){
427 maxpt[algorithm][njet] = trackpt;
635693ee 428 }
429 }
430 }// track Loop
431
432 if(algorithm==0){
433 fH1jetMCAKT04_pt[0]->Fill(Jet_pt[algorithm][njet]);
434 if(maxpt[algorithm][njet]<40)fH1jetMCAKT04_pt[1]->Fill(Jet_pt[algorithm][njet]);
435 if(ntracks>1) fH1jetMCAKT04_pt[2]->Fill(Jet_pt[algorithm][njet]);
436 if(ntracks>2) fH1jetMCAKT04_pt[3]->Fill(Jet_pt[algorithm][njet]);
437 if(ntracks>3) fH1jetMCAKT04_pt[4]->Fill(Jet_pt[algorithm][njet]);
438 if((maxpt[algorithm][njet]<40)&&(ntracks>1))fH1jetMCAKT04_pt[5]->Fill(Jet_pt[algorithm][njet]);
439 fH2jetMCAKT04_Jetpt_maxpt->Fill(maxpt[algorithm][njet],Jet_pt[algorithm][njet]);
440 }
441 if(algorithm==1){
442 fH1jetAKT04_pt[0]->Fill(Jet_pt[algorithm][njet]);
443 if(maxpt[algorithm][njet]<40)fH1jetAKT04_pt[1]->Fill(Jet_pt[algorithm][njet]);
444 if(ntracks>1) fH1jetAKT04_pt[2]->Fill(Jet_pt[algorithm][njet]);
445 if(ntracks>2) fH1jetAKT04_pt[3]->Fill(Jet_pt[algorithm][njet]);
446 if(ntracks>3) fH1jetAKT04_pt[4]->Fill(Jet_pt[algorithm][njet]);
447 if((maxpt[algorithm][njet]<40)&&(ntracks>1))fH1jetAKT04_pt[5]->Fill(Jet_pt[algorithm][njet]);
448 fH2jetAKT04_Jetpt_maxpt->Fill(maxpt[algorithm][njet],Jet_pt[algorithm][njet]);
449 }
450 }
451 }
452 if(!(IsMC&&algorithm==1))continue;
453 for(int njetMC =0;njetMC<Jet_n[0];njetMC++){
454 double eta_cut_Jet=0.5;
455 if(TMath::Abs(Jet_eta[0][njetMC])<eta_cut_Jet){
456 //Find muched jet pare=====================================
457 for(int cut=0;cut<6;cut++){
458 double min_R=10.;
459 for(int njetAOD=0;njetAOD<Jet_n[1];njetAOD++){
f21feb54 460 jets = dynamic_cast <TClonesArray*> (fAODIn->FindListObject(fJetBranch.Data()));
461 if(!jets)continue;
462 jetsAOD = (AliAODJet*) (jets->At(njetAOD));
463 TRefArray *reftracks = jetsAOD->GetRefTracks();
464 int ntracks=reftracks->GetEntriesFast();
635693ee 465 if(cut==1){if(maxpt[1][njetAOD]>=40.)continue;}
466 if(cut==2){if(ntracks==1)continue;}
467 if(cut==3){if(ntracks<=2)continue;}
468 if(cut==4){if(ntracks<=3)continue;}
469 if(cut==5){if(maxpt[1][njetAOD]>=40.)continue;if(ntracks==1)continue;}
470 double DelR = DeltaR(Jet_phi[0][njetMC],Jet_phi[1][njetAOD],Jet_eta[0][njetMC],Jet_eta[1][njetAOD]);
471 if(DelR<min_R){
472 nearrecID[njetMC]=njetAOD;
473 min_R=DelR;
474 }
475 }
476 if(min_R<0.4){
477 min_R=10.;
478 int neargenID=99999;
479 for(int njet =0;njet<Jet_n[0];njet++){
480 double DelR = DeltaR(Jet_phi[0][njet],Jet_phi[1][nearrecID[njetMC]],Jet_eta[0][njet],Jet_eta[1][nearrecID[njetMC]]);
481 if(DelR<min_R){
482 neargenID=njet;
483 min_R=DelR;
484 }
485 }
486 if((min_R<0.4)&&(neargenID==njetMC))fFIND[cut][njetMC]=true;
487 }
488 }
489 //======================================================
490 }
491 }
492 }//algorithm
493 if(IsMC){
494 for(int njetMC =0;njetMC<Jet_n[0];njetMC++){
495 double eta_cut_Jet=0.5;
496 if(TMath::Abs(Jet_eta[0][njetMC])<eta_cut_Jet){
497 for(int cut=0;cut<6;cut++){
498 if(fFIND[cut][njetMC]==true){
499 fH1jetMCAKT04_match[cut]->Fill(Jet_pt[0][njetMC]);
500 fH2jetMCAKT04_Eratio[cut]->Fill(Jet_pt[0][njetMC],Jet_pt[1][nearrecID[njetMC]]/Jet_pt[0][njetMC]);
501 }
502 }
503 }//eta window
504 }//njet loop
505 }
506
507 // End of jet analysis ---------
508 // Post output data.
509 PostData(1, fHistList);
510 return;
511}
512
513//________________________________________________________________________
514void AliAnalysisTaskCheckSingleTrackJetRejection::Terminate(Option_t *)
515{
516 // Terminate analysis
517 //
518 if (fDebug) printf("AnalysisTaskPt: Terminate() \n");
519
520}
521
522
523Bool_t AliAnalysisTaskCheckSingleTrackJetRejection::JetSelected(AliAODJet *jet){
524 Bool_t selected = false;
525
526 if(!jet)return selected;
527
528 if(fabs(jet->Eta())<fJetRecEtaWindow&&jet->Pt()>fMinJetPt){
529 selected = kTRUE;
530 }
531 return selected;
532
533}
534
535
536Double_t AliAnalysisTaskCheckSingleTrackJetRejection::DeltaPhi(Double_t phi1,Double_t phi2){
537 Float_t pi=TMath::Pi();
538 Double_t dphi = phi1-phi2;
539 if (dphi<(-1./2*pi))dphi = dphi +2*pi;
540 else if(dphi>( 3./2*pi))dphi = dphi -2*pi;
541 return dphi;
542}
543Double_t AliAnalysisTaskCheckSingleTrackJetRejection::DeltaR(Double_t phi1,Double_t phi2,Double_t eta1,Double_t eta2){
544 Float_t pi=TMath::Pi();
545 Double_t dphi = DeltaPhi(phi1,phi2);
546 if (dphi<(-1./2*pi))dphi = dphi +2*pi;
547 else if(dphi>( 3./2*pi))dphi = dphi -2*pi;
548 Double_t Deta = eta1 - eta2;
549 Double_t deltaR = sqrt(pow(dphi,2)+pow(Deta,2));
550 return deltaR;
551}
552
553Bool_t AliAnalysisTaskCheckSingleTrackJetRejection::PythiaInfoFromFile(const char* currFile,Float_t &fXsec,Float_t &fTrials){
554 //
555 // get the cross section and the trails either from pyxsec.root or from pysec_hists.root
556 // This is to called in Notify and should provide the path to the AOD/ESD file
557 // Copied from AliAnalysisTaskJetSpectrum2
558 //
559
560 TString file(currFile);
561 fXsec = 0;
562 fTrials = 1;
563
564 if(file.Contains("root_archive.zip#")){
565 Ssiz_t pos1 = file.Index("root_archive",12,TString::kExact);
566 Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
567 file.Replace(pos+1,20,"");
568 }
569 else {
570 // not an archive take the basename....
571 file.ReplaceAll(gSystem->BaseName(file.Data()),"");
572 }
573 // Printf("%s",file.Data());
574
575
576 TFile *filexsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")); // problem that we cannot really
577 //test the existance of a file in a archive so we have to lvie with open error message from root
578 if(!filexsec){
579 // next trial fetch the histgram file
580 filexsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
581 if(!filexsec){
582 // not a severe condition but inciate that we have no information
583 return kFALSE;
584 }
585 else{
586 // find the tlist we want to be independtent of the name so use the Tkey
587 TKey* key = (TKey*)filexsec->GetListOfKeys()->At(0);
588 if(!key){
589 filexsec->Close();
590 return kFALSE;
591 }
592 TList *list = dynamic_cast<TList*>(key->ReadObj());
593 if(!list){
594 filexsec->Close();
595 return kFALSE;
596 }
597 fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
598 fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
599 filexsec->Close();
600 }
601 } // no tree pyxsec.root
602 else {
603 TTree *xtree = (TTree*)filexsec->Get("Xsection");
604 if(!xtree){
605 filexsec->Close();
606 return kFALSE;
607 }
608 UInt_t ntrials = 0;
609 Double_t xsection = 0;
610 xtree->SetBranchAddress("xsection",&xsection);
611 xtree->SetBranchAddress("ntrials",&ntrials);
612 xtree->GetEntry(0);
613 fTrials = ntrials;
614 fXsec = xsection;
615 filexsec->Close();
616 }
617 return kTRUE;
618}
619
620
f4b5b842 621//Float_t AliAnalysisTaskCheckSingleTrackJetRejection::GetTotalEvents(const char* currFile){
622// Float_t totalevent;
623// TString file_es(currFile);
624// if(file_es.Contains("root_archive.zip#")){
625// Ssiz_t pos1 = file_es.Index("root_archive",12,TString::kExact);
626// Ssiz_t pos = file_es.Index("#",1,pos1,TString::kExact);
627// file_es.Replace(pos+1,20,"");
628// }
629// else {
630// // not an archive take the basename....
631// file_es.ReplaceAll(gSystem->BaseName(file_es.Data()),"");
632// }
633//
634// TString cAdd = "";
635// cAdd += Form("%02d_",(int)((Radius+0.01)*10.));
636// cAdd += Form("B%d",(int)BackM);
637// cAdd += Form("_Filter%05d",Filtermask);
638// cAdd += Form("_Cut%05d",(int)(1000.*TrackPtcut));
639// cAdd += Form("_Skip%02d",SkipCone);
640// TString Dirname,Listname;
641// Dirname = Form("PWG4_cluster_AOD__%s%s",JFAlg.Data(),cAdd.Data());
642// Listname = Form("pwg4cluster_AOD__%s%s",JFAlg.Data(),cAdd.Data());
643//
644// TFile *feventstat = TFile::Open(Form("%s%s",file_es.Data(),"JetTasksOutput.root"));
645// gROOT->Cd(Dirname.Data());
646// TList *templist = (TList*)gROOT->FindObject(Listname.Data());
647// TH1F* temphist = (TH1F*)templist->FindObject("fh1Trials");
648// totalevent = temphist->Integral();
649// //cout<<temphist->Integral()<<endl;
650// delete feventstat;
651// return totalevent;
652//
653//}
635693ee 654
655