AliPhysics  2ad5f07 (2ad5f07)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalJetExtractor.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Jet extractor task
4 // Stores jets (e.g. as AliEmcalJet objects) in trees
5 
6 // Author: R. Haake
7 
8 #include "AliParticleContainer.h"
9 #include "AliEmcalJet.h"
10 #include "AliAODTrack.h"
11 #include "TRandom3.h"
13 
14 
18 
19 //________________________________________________________________________
21 {
22 // dummy destructor
23 }
24 
25 //________________________________________________________________________
27 {
28 // dummy destructor
29 }
30 
31 
32 //________________________________________________________________________
34  AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetExtractor", kTRUE),
35  fJetsCont(0),
36  fTracksCont(0),
37  fJetBuffer(0),
38  fJetsOutput(0),
39  fCounter(0),
40  fRandom(0),
41  fIsExtractionDefined(0),
42  fExtractionType(0),
43  fExtractionCriterium(0),
44  fExtractionMinPt(0),
45  fExtractionMaxPt(0),
46  fExtractionPercentage(1.0)
47 {
48  // Default constructor.
50 }
51 
52 //________________________________________________________________________
54  AliAnalysisTaskEmcalJet(name, kTRUE),
55  fJetsCont(0),
56  fTracksCont(0),
57  fJetBuffer(0),
58  fJetsOutput(0),
59  fCounter(0),
60  fRandom(0),
61  fIsExtractionDefined(0),
62  fExtractionType(0),
63  fExtractionCriterium(0),
64  fExtractionMinPt(0),
65  fExtractionMaxPt(0),
66  fExtractionPercentage(1.0)
67 {
69 }
70 
71 //________________________________________________________________________
73 {
74  if(fRandom) delete fRandom;
75  if(fJetsOutput) delete fJetsOutput;
76 }
77 
78 //________________________________________________________________________
79 void AliAnalysisTaskEmcalJetExtractor::DefineExtraction(Int_t type, Int_t criterium, Double_t minPt, Double_t maxPt, Double_t percentage)
80 {
82  {
83  fIsExtractionDefined = kTRUE;
84  fExtractionType = type;
85  fExtractionCriterium = criterium;
86  fExtractionMinPt = minPt;
87  fExtractionMaxPt = maxPt;
88  fExtractionPercentage = percentage;
89  }
90  else
91  AliFatal("Tried to define extraction twice -- aborting.");
92 }
93 
94 //________________________________________________________________________
96 {
98  fRandom = new TRandom3(0);
99  // ### Basic container settings
101  if(fJetsCont) { //get particles connected to jets
103  } else { //no jets, just analysis tracks
104  fTracksCont = dynamic_cast<AliTrackContainer*>(GetTrackContainer(0));
105  }
106  if(fTracksCont) fTracksCont->SetClassName("AliAODTrack");
107 
108  // Histograms
109  AddHistogram2D<TH2D>("hTrackPhiEta", "Track angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Tracks}/d#phi d#eta");
110  AddHistogram2D<TH2D>("hJetPhiEta", "Jet angular distribution #phi/#eta", "COLZ", 180, 0., 2*TMath::Pi(), 100, -2.5, 2.5, "#phi", "#eta", "dN^{Jets}/d#phi d#eta");
111 
112  PrintSettings();
113 
114  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
115 }
116 
117 //________________________________________________________________________
119 {
120  std::cout << std::endl << Form("########### Printing settings for task %s", GetName()) << std::endl;
121  std::cout << "Jet container cuts:\n";
122  fJetsCont->PrintCuts();
123  std::cout << "---------------------------------------\n";
124  std::cout << "Particle container: " << fTracksCont->GetName() << std::endl;
125  std::cout << "Track filter: " << fTracksCont->GetTrackFilterType() << " (0 - no filter, 1 - custom, 2 - hybrid, 3 - TPC only)\n";
126  std::cout << "---------------------------------------\n";
127  std::cout << "Extraction type: " << fExtractionType << " (0 - AliEmcalJet, 1 - AliBasicJet, 2 - AliBasicJet w/constituents)\n";
128  std::cout << "Extraction criterium: " << fExtractionCriterium << " (0 - MinBias, 1 - Signal jets, 2 - Background jets)\n";
129  std::cout << "Extraction min pT (after background correction if rho given): " << fExtractionMinPt << "\n";
130  std::cout << "Extraction max pT (after background correction if rho given): " << fExtractionMaxPt << "\n";
131  std::cout << "Extraction percentage (rest is thrown away): " << fExtractionPercentage << "\n";
132  std::cout << "###########" << std::endl << std::endl;
133 }
134 
135 //________________________________________________________________________
137 {
138  return kTRUE;
139 }
140 
141 //________________________________________________________________________
144  if (fJetsCont && fJetsCont->GetArray() == 0)
145  {
146  fJetsCont = 0;
147  AliFatal("No jet container found.");
148  }
149  if (fTracksCont && fTracksCont->GetArray() == 0)
150  {
151  fTracksCont = 0;
152  AliFatal("No track container found.");
153  }
155  {
156  AliFatal("Extraction conditions not defined. Use DefineExtraction().");
157  }
158 
159 
160  fJetsOutput = new TTree("ExtractedJets", "ExtractedJets");
161  if(fExtractionType==0) // AliEmcalJet
162  fJetsOutput->Branch("Jets", "AliEmcalJet", &fJetBuffer, 1000);
163  else if(fExtractionType==1) // AliBasicJet w/o constituents
164  fJetsOutput->Branch("Jets", "AliBasicJet", &fJetBuffer, 1000);
165  else if(fExtractionType==2) // AliBasicJet w/ constituents
166  fJetsOutput->Branch("Jets", "AliBasicJet", &fJetBuffer, 1000);
167 
168  TList* list = dynamic_cast<TList*> (GetOutputData(1));
169  list->Add(fJetsOutput);
170 
171 }
172 
173 //________________________________________________________________________
175 {
176  Long64_t eventID = InputEvent()->GetHeader()->GetEventIdAsLong();
177 
179  AliAODTrack *track = static_cast<AliAODTrack*>(fTracksCont->GetNextAcceptParticle());
180  // All tracks plots
181  while(track) {
182  FillHistogram("hTrackPhiEta", track->Phi(), track->Eta());
183  track = static_cast<AliAODTrack*>(fTracksCont->GetNextAcceptParticle());
184  }
185 
188  while(jet) {
189  // Select jets according to extraction criterium
190  if( ((jet->Pt()-jet->Area()*fJetsCont->GetRhoVal()) < fExtractionMinPt) || ((jet->Pt()-jet->Area()*fJetsCont->GetRhoVal()) >= fExtractionMaxPt) )
191  {
192  jet = fJetsCont->GetNextAcceptJet();
193  continue;
194  }
195 
196  // "Minimum bias" criterium
197  if(fExtractionCriterium==0)
198  {
199  // do nothing
200  }
201 
202  // Discard jets statistically
203  if(fRandom->Rndm() >= fExtractionPercentage)
204  {
205  jet = fJetsCont->GetNextAcceptJet();
206  continue;
207  }
208 
209  FillHistogram("hJetPhiEta", jet->Phi(), jet->Eta());
210 
211  // Create the jet object that will be saved to the tree
212  if(fExtractionType==0) // AliEmcalJet
213  {
214  fJetBuffer = jet;
215  fJetsOutput->Fill();
216  }
217  else if( (fExtractionType==1) || (fExtractionType==2) )// AliBasicJet
218  {
219  AliBasicJet basicJet(jet->Eta(), jet->Phi(), jet->Pt(), jet->Charge(), fJetsCont->GetJetRadius(), jet->Area(), fJetsCont->GetRhoVal(), eventID, fCent);
220  if(fExtractionType==1)
221  fJetBuffer = &basicJet;
222  else if(fExtractionType==2) // AliBasicJet w/ constituents
223  {
224  for(Int_t i = 0; i < jet->GetNumberOfTracks(); i++)
225  {
226  AliVParticle* particle = static_cast<AliVParticle*>(jet->TrackAt(i, fTracksCont->GetArray()));
227  basicJet.AddJetConstituent(particle->Eta(), particle->Phi(), particle->Pt(), particle->Charge());
228  fJetBuffer = &basicJet;
229  }
230  }
231  fJetsOutput->Fill();
232  }
233  jet = fJetsCont->GetNextAcceptJet();
234  }
235 
236  return kTRUE;
237 
238 }
239 
240 //########################################################################
241 // HISTOGRAM HELPER FUNCTIONS
242 //########################################################################
243 
244 //________________________________________________________________________
245 inline void AliAnalysisTaskEmcalJetExtractor::FillHistogram(const char * key, Double_t x)
246 {
247  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
248  if(!tmpHist)
249  {
250  AliError(Form("Cannot find histogram <%s> ",key)) ;
251  return;
252  }
253 
254  tmpHist->Fill(x);
255 }
256 
257 //________________________________________________________________________
258 inline void AliAnalysisTaskEmcalJetExtractor::FillHistogram(const char * key, Double_t x, Double_t y)
259 {
260  TH1* tmpHist = static_cast<TH1*>(fOutput->FindObject(key));
261  if(!tmpHist)
262  {
263  AliError(Form("Cannot find histogram <%s> ",key));
264  return;
265  }
266 
267  if (tmpHist->IsA()->GetBaseClass("TH1"))
268  static_cast<TH1*>(tmpHist)->Fill(x,y); // Fill x with y
269  else if (tmpHist->IsA()->GetBaseClass("TH2"))
270  static_cast<TH2*>(tmpHist)->Fill(x,y); // Fill x,y with 1
271 }
272 
273 //________________________________________________________________________
274 inline void AliAnalysisTaskEmcalJetExtractor::FillHistogram(const char * key, Double_t x, Double_t y, Double_t add)
275 {
276  TH2* tmpHist = static_cast<TH2*>(fOutput->FindObject(key));
277  if(!tmpHist)
278  {
279  AliError(Form("Cannot find histogram <%s> ",key));
280  return;
281  }
282 
283  tmpHist->Fill(x,y,add);
284 }
285 
286 //________________________________________________________________________
287 template <class T> T* AliAnalysisTaskEmcalJetExtractor::AddHistogram1D(const char* name, const char* title, const char* options, Int_t xBins, Double_t xMin, Double_t xMax, const char* xTitle, const char* yTitle)
288 {
289  T* tmpHist = new T(name, title, xBins, xMin, xMax);
290 
291  tmpHist->GetXaxis()->SetTitle(xTitle);
292  tmpHist->GetYaxis()->SetTitle(yTitle);
293  tmpHist->SetOption(options);
294  tmpHist->SetMarkerStyle(kFullCircle);
295  tmpHist->Sumw2();
296 
297  fOutput->Add(tmpHist);
298 
299  return tmpHist;
300 }
301 
302 //________________________________________________________________________
303 template <class T> T* AliAnalysisTaskEmcalJetExtractor::AddHistogram2D(const char* name, const char* title, const char* options, Int_t xBins, Double_t xMin, Double_t xMax, Int_t yBins, Double_t yMin, Double_t yMax, const char* xTitle, const char* yTitle, const char* zTitle)
304 {
305  T* tmpHist = new T(name, title, xBins, xMin, xMax, yBins, yMin, yMax);
306  tmpHist->GetXaxis()->SetTitle(xTitle);
307  tmpHist->GetYaxis()->SetTitle(yTitle);
308  tmpHist->GetZaxis()->SetTitle(zTitle);
309  tmpHist->SetOption(options);
310  tmpHist->SetMarkerStyle(kFullCircle);
311  tmpHist->Sumw2();
312 
313  fOutput->Add(tmpHist);
314 
315  return tmpHist;
316 }
Short_t Charge() const
Definition: AliEmcalJet.h:62
Double_t Area() const
Definition: AliEmcalJet.h:69
Double_t GetRhoVal() const
const char * title
Definition: MakeQAPdf.C:26
AliJetContainer * GetJetContainer(Int_t i=0) const
Double_t Eta() const
Definition: AliEmcalJet.h:60
Double_t Phi() const
Definition: AliEmcalJet.h:55
Container with name, TClonesArray and cuts for particles.
TList * list
TTree * fJetsOutput
buffer for one jet (that will be saved to the tree)
TList * fOutput
!output list
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:83
AliParticleContainer * GetParticleContainer() const
Double_t fCent
!event centrality
TClonesArray * GetArray() const
ETrackFilterType_t GetTrackFilterType() const
void SetClassName(const char *clname)
AliEmcalJet * GetNextAcceptJet()
Double_t Pt() const
Definition: AliEmcalJet.h:47
const char * GetName() const
Float_t GetJetRadius() const
Short_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:107
T * AddHistogram1D(const char *name="CustomHistogram", const char *title="NO_TITLE", const char *options="", Int_t xBins=100, Double_t xMin=0.0, Double_t xMax=20.0, const char *xTitle="x axis", const char *yTitle="y axis")
AliTrackContainer * GetTrackContainer(Int_t i=0) const
void SetMakeGeneralHistograms(Bool_t g)
void DefineExtraction(Int_t type, Int_t criterium, Double_t minPt, Double_t maxPt, Double_t percentage)
ClassImp(AliAnalysisTaskEmcalJetExtractor) ClassImp(AliBasicJet) ClassImp(AliBasicJetConstituent) AliBasicJet
virtual AliVParticle * GetNextAcceptParticle()
T * AddHistogram2D(const char *name="CustomHistogram", const char *title="NO_TITLE", const char *options="", Int_t xBins=100, Double_t xMin=0.0, Double_t xMax=20.0, Int_t yBins=100, Double_t yMin=0.0, Double_t yMax=20.0, const char *xTitle="x axis", const char *yTitle="y axis", const char *zTitle="z axis")
void ResetCurrentID(Int_t i=-1)
void FillHistogram(const char *key, Double_t x)