1 // **************************************
2 // Task used for jet finding in the Kine Train (generation and analysis on the fly, no detector effects)
3 // Output is stored in an exchance container
4 // *******************************************
7 /**************************************************************************
8 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
10 * Author: The ALICE Off-line Project. *
11 * Contributors are mentioned in the code where appropriate. *
13 * Permission to use, copy, modify and distribute this software and its *
14 * documentation strictly for non-commercial purposes is hereby granted *
15 * without fee, provided that the above copyright notice appears in all *
16 * copies and that both the copyright notice and this permission notice *
17 * appear in the supporting documentation. The authors make no claims *
18 * about the suitability of this software for any purpose. It is *
19 * provided "as is" without express or implied warranty. *
20 **************************************************************************/
26 #include <TInterpreter.h>
28 #include <TRefArray.h>
37 #include <TLorentzVector.h>
38 #include <TClonesArray.h>
39 #include "TDatabasePDG.h"
42 #include "AliAnalysisTaskJetClusterKine.h"
43 #include "AliAnalysisManager.h"
44 #include "AliJetFinder.h"
45 #include "AliJetHeader.h"
46 #include "AliJetReader.h"
47 #include "AliESDEvent.h"
48 #include "AliAODEvent.h"
49 #include "AliAODHandler.h"
50 #include "AliAODExtension.h"
51 #include "AliAODTrack.h"
52 #include "AliAODJet.h"
53 #include "AliVParticle.h" //FK//
54 #include "AliAODMCParticle.h"
55 #include "AliMCEventHandler.h"
56 #include "AliMCEvent.h"
58 #include "AliGenEventHeader.h" //FK//
59 #include "AliGenPythiaEventHeader.h"
60 #include "AliJetKineReaderHeader.h"
61 #include "AliGenCocktailEventHeader.h"
62 #include "AliInputEventHandler.h"
63 #include "AliAODJetEventBackground.h"
65 #include "fastjet/PseudoJet.hh"
66 #include "fastjet/ClusterSequenceArea.hh"
67 #include "fastjet/AreaDefinition.hh"
68 #include "fastjet/JetDefinition.hh"
69 // get info on how fastjet was configured
70 #include "fastjet/config.h"
74 ClassImp(AliAnalysisTaskJetClusterKine)
76 AliAnalysisTaskJetClusterKine::~AliAnalysisTaskJetClusterKine(){
83 if(fTCAJetsOut)fTCAJetsOut->Delete();
88 //_____________________________________________________________________
90 AliAnalysisTaskJetClusterKine::AliAnalysisTaskJetClusterKine():
95 fTrackTypeGen(kTrackKineCharged), //Kine charged?
99 fJetOutputMinPt(0.150),
100 fMaxTrackPtInJet(100.),
103 fOutContainer(kNoOutput), //FF//
106 fAlgorithm(fastjet::kt_algorithm),
107 fStrategy(fastjet::Best),
108 fRecombScheme(fastjet::BIpt_scheme),
109 fAreaType(fastjet::active_area),
111 fActiveAreaRepeats(1),
118 fh1PtHardTrials(0x0),
121 fh1NConstLeadingGen(0x0),
123 fh1PtJetsLeadingGenIn(0x0),
124 fh1PtJetConstGen(0x0),
125 fh1PtJetConstLeadingGen(0x0),
126 fh1PtTracksGenIn(0x0),
130 fh2NConstLeadingPt(0x0),
132 fh2LeadingJetPhiEta(0x0),
134 fh2LeadingJetEtaPt(0x0),
136 fh2JetsLeadingPhiEta(0x0),
137 fh2JetsLeadingPhiPt(0x0),
138 fh2JetsLeadingPhiPtW(0x0),
146 //_____________________________________________________________________
148 AliAnalysisTaskJetClusterKine::AliAnalysisTaskJetClusterKine(const char* name):
149 AliAnalysisTaskSE(name),
153 fTrackTypeGen(kTrackKineCharged), //kine charged?
155 fTrackEtaWindow(0.9),
157 fJetOutputMinPt(0.150),
158 fMaxTrackPtInJet(100.),
161 fOutContainer(kNoOutput),//FF//
164 fAlgorithm(fastjet::kt_algorithm),
165 fStrategy(fastjet::Best),
166 fRecombScheme(fastjet::BIpt_scheme),
167 fAreaType(fastjet::active_area),
169 fActiveAreaRepeats(1),
176 fh1PtHardTrials(0x0),
179 fh1NConstLeadingGen(0x0),
181 fh1PtJetsLeadingGenIn(0x0),
182 fh1PtJetConstGen(0x0),
183 fh1PtJetConstLeadingGen(0x0),
184 fh1PtTracksGenIn(0x0),
188 fh2NConstLeadingPt(0x0),
190 fh2LeadingJetPhiEta(0x0),
192 fh2LeadingJetEtaPt(0x0),
194 fh2JetsLeadingPhiEta(0x0),
195 fh2JetsLeadingPhiPt(0x0),
196 fh2JetsLeadingPhiPtW(0x0),
202 DefineOutput(1, TList::Class());
203 DefineOutput(2, TClonesArray::Class());
207 //_____________________________________________________________________
209 Bool_t AliAnalysisTaskJetClusterKine::Notify(){
211 // Implemented Notify() to read the cross sections
212 // and number of trials from pyxsec.root
217 //_____________________________________________________________________
219 void AliAnalysisTaskJetClusterKine::UserCreateOutputObjects(){
222 // Create the output container
228 if(fDebug > 1) printf("AnalysisTaskJetCluster::UserCreateOutputObjects() \n");
231 if(fNonStdBranch.Length()!=0){
232 // only create the output branch if we have a name
233 // Create a new branch for jets...
234 // -> cleared in the UserExec....
235 // here we can also have the case that the brnaches are written to a separate file
236 fTCAJetsOut = new TClonesArray("AliAODJet", 0);
237 fTCAJetsOut->SetName(fNonStdBranch.Data());
238 if(fOutContainer==kAODBranch){ //FF//
239 AddAODBranch("TClonesArray",&fTCAJetsOut,fNonStdFile.Data());
243 //if(fNonStdFile.Length()!=0){
245 // case that we have an AOD extension we need to fetch the jets from the extended output
246 // we identify the extension aod event by looking for the branchname
247 //AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
248 // case that we have an AOD extension we need can fetch the background maybe from the extended output
249 //fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
254 if(!fHistList) fHistList = new TList();
255 fHistList->SetOwner(kTRUE);
256 PostData(1, fHistList); // post data in any case once
259 if(fOutContainer==kExchCont){
260 fTCAJetsOut->SetOwner();
261 PostData(2, fTCAJetsOut); //FF// post data in any case once
264 Bool_t oldStatus = TH1::AddDirectoryStatus();
265 TH1::AddDirectory(kFALSE);
271 const Int_t nBinPt = 100;
272 Double_t binLimitsPt[nBinPt+1];
273 for(Int_t iPt = 0;iPt <= nBinPt;iPt++){
275 binLimitsPt[iPt] = 0.0;
277 binLimitsPt[iPt] = binLimitsPt[iPt-1] + 2.0;
281 const Int_t nBinPhi = 90;
282 Double_t binLimitsPhi[nBinPhi+1];
283 for(Int_t iPhi = 0; iPhi<=nBinPhi; iPhi++){
285 binLimitsPhi[iPhi] = -1.*TMath::Pi();
287 binLimitsPhi[iPhi] = binLimitsPhi[iPhi-1] + 1/(Float_t)nBinPhi * TMath::Pi()*2;
291 const Int_t nBinEta = 40;
292 Double_t binLimitsEta[nBinEta+1];
293 for(Int_t iEta = 0;iEta<=nBinEta;iEta++){
295 binLimitsEta[iEta] = -2.0;
297 binLimitsEta[iEta] = binLimitsEta[iEta-1] + 0.1;
301 const int nChMax = 5000;
303 fh1Xsec = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
304 fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
306 fh1Trials = new TH1F("fh1Trials","trials root file",1,0,1);
307 fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
310 fh1NJetsGen = new TH1F("fh1NJetsGen","N reconstructed jets",120,-0.5,119.5);
312 fh1NConstGen = new TH1F("fh1NConstGen","# jet constituents",120,-0.5,119.5);
313 fh1NConstLeadingGen = new TH1F("fh1NConstLeadingGen","jet constituents",120,-0.5,119.5);
316 fh1PtHard = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",nBinPt,binLimitsPt);
317 fh1PtHardNoW = new TH1F("fh1PtHardNoW","PYTHIA Pt hard no weight;p_{T,hard}",nBinPt,binLimitsPt);
318 fh1PtHardTrials = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",nBinPt,binLimitsPt);
320 fh1PtJetsGenIn = new TH1F("fh1PtJetsGenIn","Gen jets P_T;p_{T} (GeV/c)",nBinPt,binLimitsPt);
321 fh1PtJetsLeadingGenIn = new TH1F("fh1PtJetsLeadingGenIn","Gen jets P_T;p_{T} (GeV/c)",nBinPt,binLimitsPt);
322 fh1PtJetConstGen = new TH1F("fh1PtJetsConstGen","Rec jets constituents P_T;p_{T} (GeV/c)",nBinPt,binLimitsPt);
323 fh1PtJetConstLeadingGen = new TH1F("fh1PtJetsConstLeadingGen","Gen jets constituents P_T;p_{T} (GeV/c)",nBinPt,binLimitsPt);
324 fh1PtTracksGenIn = new TH1F("fh1PtTracksGenIn",Form("Gen tracks P_T #eta < %1.2f ;p_{T} (GeV/c)",fTrackEtaWindow),nBinPt,binLimitsPt);
325 fh1Nch = new TH1F("fh1Nch","charged multiplicity; N_{ch}",nChMax,-0.5,nChMax-0.5);
328 fh1Z = new TH1F("fh1Z",";zvtx",100,-25,25);
332 fh2NConstPt = new TH2F("fh2NConstPt","Number of constituents ;p_{T} (GeV/c);N",nBinPt,binLimitsPt,50,-0.5,49.5);
333 fh2NConstLeadingPt = new TH2F("fh2NConstLeadingPt","Number of constituents ;p_{T} (GeV/c);N",nBinPt,binLimitsPt,50,-0.5,49.5);
336 fh2JetPhiEta = new TH2F("fh2JetPhiEta","eta vs phi all jets;#phi;#eta",
337 nBinPhi,0.,2.*TMath::Pi(),nBinEta,binLimitsEta);
338 fh2LeadingJetPhiEta = new TH2F("fh2LeadingJetPhiEta","eta vs phi leading jets;#phi;#eta",
339 nBinPhi,0.,2.*TMath::Pi(),nBinEta,binLimitsEta);
341 fh2JetEtaPt = new TH2F("fh2JetEtaPt","pt vs eta all jets;#eta;p_{T}",
342 nBinEta,binLimitsEta,nBinPt,binLimitsPt);
343 fh2LeadingJetEtaPt = new TH2F("fh2LeadingJetEtaPt","pT vs eta leading jets;#eta;p_{T}",
344 nBinEta,binLimitsEta,nBinPt,binLimitsPt);
346 fh2TrackEtaPt = new TH2F("fh2TrackEtaPt","pt vs eta all jets;#eta;p_{T}",
347 nBinEta,binLimitsEta,nBinPt,binLimitsPt);
351 fh2JetsLeadingPhiEta = new TH2F("fh2JetsLeadingPhiEta","delta eta vs delta phi to leading jet;#Delta#phi;#Delta#eta",
352 nBinPhi,binLimitsPhi,nBinEta,binLimitsEta);
353 fh2JetsLeadingPhiPt = new TH2F("fh2JetsLeadingPhiPt","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
354 nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
356 fh2JetsLeadingPhiPtW = new TH2F("fh2JetsLeadingPhiPtW","leading p_T vs delta phi p_T weigted to leading jet;#Delta#phi;p_{T} (GeV/c)",
357 nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
362 const Int_t saveLevel = 3; // large save level more histos
364 fHistList->Add(fh1Xsec);
365 fHistList->Add(fh1Trials);
367 fHistList->Add(fh1NJetsGen);
368 fHistList->Add(fh1NConstGen);
369 fHistList->Add(fh1NConstLeadingGen);
370 fHistList->Add(fh1PtJetsGenIn);
371 fHistList->Add(fh1PtJetsLeadingGenIn);
372 fHistList->Add(fh1PtTracksGenIn);
373 fHistList->Add(fh1PtJetConstGen);
374 fHistList->Add(fh1PtJetConstLeadingGen);
375 fHistList->Add(fh1Nch);
376 fHistList->Add(fh1Z);
378 fHistList->Add(fh2NConstPt);
379 fHistList->Add(fh2NConstLeadingPt);
380 fHistList->Add(fh2JetPhiEta);
381 fHistList->Add(fh2LeadingJetPhiEta);
382 fHistList->Add(fh2JetEtaPt);
383 fHistList->Add(fh2LeadingJetEtaPt);
384 fHistList->Add(fh2TrackEtaPt);
385 fHistList->Add(fh2JetsLeadingPhiEta);
386 fHistList->Add(fh2JetsLeadingPhiPt);
387 fHistList->Add(fh2JetsLeadingPhiPtW);
390 // =========== Switch on Sumw2 for all histos ===========
391 for(Int_t i=0; i<fHistList->GetEntries(); ++i){
392 TH1 *h1 = dynamic_cast<TH1*>(fHistList->At(i));
397 THnSparse *hn = dynamic_cast<THnSparse*>(fHistList->At(i));
400 TH1::AddDirectory(oldStatus);
402 //_________________________________________________________________________
404 void AliAnalysisTaskJetClusterKine::Init(){
406 if(fDebug > 1) printf("AnalysisTaskJetClusterKine::Init() \n");
407 fMcHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); //FK//
410 //_________________________________________________________________________
411 void AliAnalysisTaskJetClusterKine::LocalInit(){
417 if(fDebug > 1) printf("AnalysisTaskJetClusterKine::LocalInit() \n");
421 //_________________________________________________________________________
422 void AliAnalysisTaskJetClusterKine::UserExec(Option_t* /*option*/){
423 // handle and reset the output jet branch
425 if(fDebug > 1) printf("AliAnalysisTaskJetClusterKine::UserExec \n");
426 if(fTCAJetsOut) fTCAJetsOut->Delete(); //clean TClonesArray
429 // Execute analysis for current event
433 fMcEvent = fMcHandler->MCEvent();
435 if(fDebug > 1) printf("AnalysisTaskJetClusterKine::Handler() fMcHandler=NULL\n");
436 PostData(1, fHistList);
438 if(fOutContainer==kExchCont){
439 PostData(2, fTCAJetsOut); //FF//
445 if(fDebug > 1) printf("AnalysisTaskJetClusterKine::Exec() fMcEvent=NULL \n");
446 PostData(1, fHistList);
448 if(fOutContainer==kExchCont){
449 PostData(2, fTCAJetsOut); //FF//
455 const AliVVertex *vtxMC = fMcEvent->GetPrimaryVertex();
456 Float_t zVtx = vtxMC->GetZ();
457 if(TMath::Abs(zVtx) > fVtxZCut){ //vertex cut
458 PostData(1, fHistList);
460 if(fOutContainer==kExchCont){
461 PostData(2, fTCAJetsOut); //FF//
469 fh1Trials->Fill("#sum{ntrials}",1);
470 if(fDebug > 10) Printf("%s:%d",(char*)__FILE__,__LINE__);
472 // we simply fetch the mc particles as a list of AliVParticles
474 TList genParticles; //list of particles with above min pT and good eta range
476 Int_t nParticles = GetListOfTracks(&genParticles, fTrackTypeGen); //fill the list
477 fh1Nch->Fill(nParticles);
479 if(fDebug>2) Printf("%s:%d Selected Gen tracks: %d %d",(char*)__FILE__,__LINE__, nParticles, genParticles.GetEntries());
483 vector<fastjet::PseudoJet> inputParticlesKine;
485 for(int i = 0; i < genParticles.GetEntries(); i++){ //loop over generated particles
486 AliVParticle *vp = (AliVParticle*) genParticles.At(i);
487 fastjet::PseudoJet jInp(vp->Px(),vp->Py(),vp->Pz(),vp->P());
488 jInp.set_user_index(i);
489 inputParticlesKine.push_back(jInp);
493 fRef->Delete(); // make sure to delete before placement new...
494 new(fRef) TRefArray(TProcessID::GetProcessWithUID(vp)); //TRefArray does not work with toy model ...
496 fRef->Add(vp); //TRefArray does not work with toy model ...
498 }//generated particles
500 if(inputParticlesKine.size()==0){ //FK//
501 if(fDebug)Printf("%s:%d No input particles found, skipping event",(char*)__FILE__,__LINE__);
502 PostData(1, fHistList);
504 if(fOutContainer==kExchCont){
505 PostData(2, fTCAJetsOut); //FF//
512 // employ setters for these...
513 // now create the object that holds info about ghosts
515 fastjet::GhostedAreaSpec ghostSpec(fGhostEtamax, fActiveAreaRepeats, fGhostArea);
516 fastjet::AreaType areaType = fastjet::active_area;
517 fastjet::AreaDefinition areaDef = fastjet::AreaDefinition(areaType,ghostSpec);
518 fastjet::JetDefinition jetDef(fAlgorithm, fRparam, fRecombScheme, fStrategy);
519 fastjet::ClusterSequenceArea clustSeq(inputParticlesKine, jetDef, areaDef); //FK//
521 //range where to compute background
522 Double_t phiMin = 0, phiMax = 0, rapMin = 0, rapMax = 0;
524 phiMax = 2*TMath::Pi();
525 rapMax = fGhostEtamax - fRparam;
526 rapMin = - fGhostEtamax + fRparam;
527 fastjet::RangeDefinition range(rapMin,rapMax, phiMin, phiMax);
530 const vector <fastjet::PseudoJet> &inclusiveJets = clustSeq.inclusive_jets();
531 const vector <fastjet::PseudoJet> &sortedJets = sorted_by_pt(inclusiveJets);
534 fh1NJetsGen->Fill(sortedJets.size()); //the number of found jets
536 // loop over all jets an fill information, first on is the leading jet
538 Int_t nGen = inclusiveJets.size();
539 if(inclusiveJets.size()>0){
542 AliAODJet leadingJet (sortedJets[0].px(), sortedJets[0].py(), sortedJets[0].pz(), sortedJets[0].E());
543 Double_t area = clustSeq.area(sortedJets[0]);
544 leadingJet.SetEffArea(area,0);
545 Float_t pt = leadingJet.Pt();
546 Float_t phi = leadingJet.Phi();
547 if(phi<0) phi += TMath::Pi()*2.;
548 Float_t eta = leadingJet.Eta();
551 Int_t nAodOutJets = 0;
552 AliAODJet *aodOutJet = NULL;
555 // correlation of leading jet with tracks
556 //FK// TIterator *particleIter = genParticles.MakeIterator();//FK//
557 //FK// particleIter->Reset();
558 //FK// AliVParticle *tmpParticle = NULL;
559 //FK// while((tmpParticle = (AliVParticle*)(particleIter->Next()))){
560 //FK// Float_t tmpPt = tmpParticle->Pt();
561 //FK// // correlation
562 //FK// Float_t tmpPhi = tmpParticle->Phi();
563 //FK// if(tmpPhi<0) tmpPhi+= TMath::Pi()*2.;
564 //FK// Float_t dPhi = phi - tmpPhi;
565 //FK// if(dPhi > TMath::Pi()) dPhi = dPhi - 2.*TMath::Pi(); //-pi,pi
566 //FK// if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi(); //-pi,pi
567 //FK// fh2TracksLeadingJetPhiPt->Fill(dPhi,pt);
568 //FK// fh2TracksLeadingJetPhiPtW->Fill(dPhi,pt,tmpPt);
571 TLorentzVector vecareab;
572 for(int j = 0; j < nGen;j++){ //loop over inclusive jets
573 AliAODJet tmpGenJet (sortedJets[j].px(), sortedJets[j].py(), sortedJets[j].pz(), sortedJets[j].E());
575 Float_t tmpPt = tmpGenJet.Pt(); //incl jet pt
577 if((tmpPt > fJetOutputMinPt) && fTCAJetsOut){// cut on the non-background subtracted...
578 aodOutJet = new ((*fTCAJetsOut)[nAodOutJets++]) AliAODJet(tmpGenJet);
579 aodOutJet->GetRefTracks()->Clear();
580 Double_t area1 = clustSeq.area(sortedJets[j]);
581 aodOutJet->SetEffArea(area1,0);
582 fastjet::PseudoJet vecarea=clustSeq.area_4vector(sortedJets[j]);
583 vecareab.SetPxPyPzE(vecarea.px(),vecarea.py(),vecarea.pz(),vecarea.e());
584 aodOutJet->SetVectorAreaCharged(&vecareab);
587 fh1PtJetsGenIn->Fill(tmpPt); //incl jet Pt
588 // Fill Spectra with constituentsemacs
589 const vector<fastjet::PseudoJet> &constituents = clustSeq.constituents(sortedJets[j]);
591 fh1NConstGen->Fill(constituents.size()); //number of constituents
592 fh2NConstPt->Fill(tmpPt,constituents.size()); //number of constituents vs jet pt
594 // loop over constiutents and fill spectrum
595 AliVParticle *partLead = 0;
598 for(unsigned int ic = 0; ic < constituents.size(); ic++){
599 AliVParticle *part = (AliVParticle*) genParticles.At(constituents[ic].user_index());
601 fh1PtJetConstGen->Fill(part->Pt()); //pt of constituents
604 aodOutJet->AddTrack(fRef->At(constituents[ic].user_index()));//FK//
606 if(part->Pt()>fMaxTrackPtInJet){
607 aodOutJet->SetTrigger(AliAODJet::kHighTrackPtTriggered);
611 if(part->Pt()>ptLead){
616 if(j==0) fh1PtJetConstLeadingGen->Fill(part->Pt()); //pt of leading jet constituents
621 //set pT of leading constituent of jet
622 aodOutJet->SetPtLeading(partLead->Pt());
627 Float_t tmpPhi = tmpGenJet.Phi(); //incl jet phi
628 Float_t tmpEta = tmpGenJet.Eta(); //incl jet eta
630 if(tmpPhi<0) tmpPhi += TMath::Pi()*2.;
631 fh2JetPhiEta->Fill(tmpGenJet.Phi(),tmpEta);
632 fh2JetEtaPt->Fill(tmpEta,tmpPt);
634 if(j==0){ //leading jet
635 fh1PtJetsLeadingGenIn->Fill(tmpPt); //leading jet pt
636 fh2LeadingJetPhiEta->Fill(tmpPhi,tmpEta);
637 fh2LeadingJetEtaPt->Fill(tmpEta,tmpPt);
638 fh1NConstLeadingGen->Fill(constituents.size()); //number of constituents in leading jet
639 fh2NConstLeadingPt->Fill(tmpPt,constituents.size());
643 Float_t dPhi = phi - tmpPhi;
644 if(dPhi > TMath::Pi()) dPhi = dPhi - 2.*TMath::Pi(); //-pi,pi
645 if(dPhi<(-1.*TMath::Pi())) dPhi = dPhi + 2.*TMath::Pi(); //-pi.pi
646 Float_t dEta = eta - tmpGenJet.Eta();
647 fh2JetsLeadingPhiEta->Fill(dPhi,dEta);
648 fh2JetsLeadingPhiPt->Fill(dPhi,pt);
649 fh2JetsLeadingPhiPtW->Fill(dPhi,pt,tmpPt);
650 }// loop over reconstructed jets
651 //delete particleIter;
655 // fill track information
656 Int_t nTrackOver = genParticles.GetSize();
657 // do the same for tracks and jets
660 TIterator *particleIter = genParticles.MakeIterator();
661 AliVParticle *tmpGen = 0;
662 particleIter->Reset();
664 while((tmpGen = (AliVParticle*)(particleIter->Next()))){
665 Float_t tmpPt = tmpGen->Pt();
666 Float_t tmpEta = tmpGen->Eta();
667 fh1PtTracksGenIn->Fill(tmpPt);
668 fh2TrackEtaPt->Fill(tmpEta,tmpPt);
678 if(fTCAJetsOut)Printf("%s:%d Rec Jets %d",(char*)__FILE__,__LINE__,fTCAJetsOut->GetEntriesFast());
680 PostData(1, fHistList);
682 if(fOutContainer==kExchCont){
683 PostData(2, fTCAJetsOut); //FF//
688 //____________________________________________________________________________
690 void AliAnalysisTaskJetClusterKine::Terminate(Option_t* /*option*/){
692 // Terminate analysis
694 if (fDebug > 1) printf("AnalysisJetCluster: Terminate() \n");
699 //_____________________________________________________________________________________
701 Int_t AliAnalysisTaskJetClusterKine::GetListOfTracks(TList *list, Int_t type){
704 // get list of tracks/particles for different types
707 if(fDebug>2) Printf("%s:%d Selecting tracks with %d",(char*)__FILE__,__LINE__,type);
709 Int_t iCount = 0; //number of particles
711 if(type == kTrackKineAll || type == kTrackKineCharged){
712 if(! fMcEvent) return iCount;
714 //we want to have alivpartilces so use get track
715 for(int it = 0; it < fMcEvent->GetNumberOfTracks(); ++it){
716 if(!fMcEvent->IsPhysicalPrimary(it)) continue;
718 AliMCParticle* part = (AliMCParticle*)fMcEvent->GetTrack(it);
719 if(TMath::Abs(part->Eta()) > fTrackEtaWindow) continue; //apply eta cut
720 if(part->Pt() < fTrackPtCut) continue; //apply pT cut
722 if(type == kTrackKineAll){ // full jets
726 }else if(type == kTrackKineCharged){ //charged jets
727 if(part->Particle()->GetPDG()->Charge()==0) continue;
736 return iCount; //number of particles in the list