AliPhysics  608b256 (608b256)
AliAnalysisTaskEmcalJetCDF.h
Go to the documentation of this file.
1 #ifndef ALIANALYSISTASKEMCALJETCDF_H
2 #define ALIANALYSISTASKEMCALJETCDF_H
3 
9 /* Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
10  * See cxx source for full Copyright notice */
11 
12 #include <TFile.h>
13 #include <TChain.h>
14 #include <TH1D.h>
15 #include <TH2D.h>
16 #include <TArrayD.h>
17 #include <TVector2.h>
18 #include <TString.h>
19 #include "THistManager.h"
20 #include "AliLog.h"
21 #include "AliAnalysisManager.h"
23 #include "AliEmcalEmbeddingQA.h"
24 
31 public:
32 
34  AliAnalysisTaskEmcalJetCDF ( const char *name );
36 
38  void Terminate ( Option_t *option );
39 
40 protected:
41  void ExecOnce();
42  Bool_t Run() ;
43 
47 
51  TObject* GetHistogram ( const char* histName );
52 
54 
56 
57  // Event selection
59  AliEventCuts fEventCuts;
62 
63  // MC options
65 
67 
68 
69 private:
70  AliAnalysisTaskEmcalJetCDF ( const AliAnalysisTaskEmcalJetCDF& ); // not implemented
72 
74  ClassDef ( AliAnalysisTaskEmcalJetCDF, 8 );
76 
77 };
78 
79 namespace PWGJE {
80 namespace EMCALJetTasks {
81 namespace AliAnalysisTaskEmcalJetCDF_NS {
82 
84  typedef std::pair<Double_t, Int_t> ptidx_pair;
85 
87  struct sort_descend
88  {
89  bool operator () ( const ptidx_pair &p1, const ptidx_pair &p2 ) { return p1.first > p2.first ; }
90  };
91 
95  inline Double_t Mag2 (const AliVParticle& trk)
96  { return trk.Px()*trk.Px() + trk.Py()*trk.Py() + trk.Pz()*trk.Pz(); }
97 
101  inline Double_t Mag(const AliVParticle& trk) { return TMath::Sqrt(Mag2(trk)); }
102 
107  inline Double_t Dot (const AliVParticle& trk1, const AliVParticle& trk2 )
108  { return trk1.Px()*trk2.Px() + trk1.Py()*trk2.Py() + trk1.Pz()*trk2.Pz(); }
109 
114  inline Double_t Perp2( const AliVParticle& trk1, const AliVParticle& trk2) {
115  Double_t mag1 = Mag2(trk1);
116  Double_t mag2 = Mag2(trk2);
117  Double_t dotp = Dot(trk1,trk2);
118  if (mag2 > 0.0) { mag1 -= dotp*dotp/mag2; }
119  if (mag1 <= 0) { mag1 = 0; }
120  return mag1;
121  }
122 
127  inline Double_t Perp (const AliVParticle& trk1, const AliVParticle& trk2 ) { return TMath::Sqrt(Perp2(trk1,trk2)); }
128 
132  std::vector<Int_t> SortTracksPt ( AliVEvent* event );
133 
137  std::vector<Int_t> SortTracksPt ( AliParticleContainer* track_container );
138 
143  inline Double_t Z_ptot ( const AliEmcalJet* jet, const AliVParticle* trk ) // Get Z of constituent trk ; p total
144  {
145  if (trk->P() < 1e-6) return 0.;
146  return (trk != 0) ? trk->P()/ jet->P() : 0.;
147  }
148 
153  inline Double_t Z_pt ( const AliEmcalJet* jet, const AliVParticle* trk ) // Get Z of constituent trk ; pt
154  {
155  if (trk->P() < 1e-6) return 0.;
156  return (trk != 0) ? trk->Pt() / jet->Pt() : 0.;
157  }
158 
161  inline Double_t Xi ( Double_t z ) { return TMath::Log ( 1/z ); } // Get Xi of value z
162 
167  inline Double_t DeltaR ( const AliVParticle* part1, const AliVParticle* part2 )
168  {
169  Double_t dPhi = part1->Phi() - part2->Phi();
170  Double_t dEta = part1->Eta() - part2->Eta();
171  dPhi = TVector2::Phi_mpi_pi ( dPhi );
172  return TMath::Sqrt ( dPhi * dPhi + dEta * dEta );
173  }
174 
179  inline Int_t IdxInArray ( Int_t index, TArrayI &array )
180  {
181  for ( Int_t i = 0; i < array.GetSize(); i++ ) { if ( index == array[i] ) { return i; } }
182  return -1;
183  }
184 
192  const char* ntracks = "usedefault",
193  const char* nclusters = "usedefault",
194  const char* ncells = "usedefault",
195  const char* ntracks_mc = "",
196  const char* tag = "CDF"
197  );
198 
210  AliJetContainer* jetCont,
211  Float_t jetptmin = 1.,
212  Float_t jetptmax = 500.,
213  Float_t jetareacutperc = 0.,
214  Int_t leadhadtype = 2,
215  Int_t nLeadJets = 1,
216  Float_t mintrackpt = 0.15,
217  Float_t maxtrackpt = 1000.
218  );
219 
220 
233 TChain* CreateChain ( const char* filelist = "filelist.txt",
234  const char* cTreeNameArg = "auto",
235  const char* friends = "",
236  UInt_t iNumFiles = -1,
237  UInt_t iStartWithFile = 1
238  );
239 
240 
241 
245 inline bool PeriodIsMC (const char* str) {
246  TString period (str);
247  if (!period.IsNull()) {
248  period.ToLower();
249  if ( period.BeginsWith("lhc") && (period.Length() > 6) ) { return true; }
250  }
251  return false;
252  }
253 
254 
258 inline TString GetPeriod (const char* file_path) {
259  TString period = "";
260  TString sFile(file_path);
261  sFile.ToLower();
262 
263  if (!sFile.IsNull()) {
264  // split string in tokens (libs)
265  TObjArray* tokens_list = sFile.Tokenize("/");
266  tokens_list->SetOwner(kTRUE);
267  TIter next_str(tokens_list);
268  TObjString* token = NULL;
269  while ((token=(TObjString*)next_str())) {
270  TString token_str = token->GetString();
271  if ( token_str.BeginsWith("lhc") ) { period = token_str; break; }
272  }
273  delete tokens_list;
274  }
275  return period;
276  }
277 
278 
282 inline TString GetPass ( const char* file_path) {
283  TString pass = "";
284  TString sFile (file_path);
285  sFile.ToLower();
286 
287  if (!sFile.IsNull()) {
288  // split string in tokens (libs)
289  TObjArray* tokens_list = sFile.Tokenize("/");
290  tokens_list->SetOwner(kTRUE);
291  TIter next_str(tokens_list);
292  TObjString* token = NULL;
293  while ((token=(TObjString*)next_str())) {
294  TString token_str = token->GetString();
295  if ( token_str.BeginsWith("pass") ) { pass = token_str; break; }
296  }
297  delete tokens_list;
298  }
299 
300  return pass;
301  }
302 
303 
307 inline TString GetFileFromPath ( const char* file_path = "" ) {
308 TString file (file_path);
309 TString file_name ("");
310 TObjArray* token_list = file.Tokenize("/");
311 token_list->SetOwner(kTRUE);
312 
313 if ( token_list->GetEntries() > 0 )
314  { file_name = ((TObjString*)token_list->At(token_list->GetLast()))->GetString(); }
315 
316 delete token_list;
317 return file_name;
318 }
319 
320 
324 inline bool SaveManager ( const char* file_name) {
325  AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
326  if ( !mgr ) { ::Error ( "SaveManager", "No analysis manager to connect to." ); return kFALSE; }
327 
328  TFile pOutFile (file_name,"RECREATE");
329  if ( ! pOutFile.cd() ) { ::Error ( "SaveManager", "Could not use the new created file" ); return kFALSE; }
330  Int_t written_bytes = mgr->Write();
331  pOutFile.Close();
332  if (written_bytes == 0 ) { ::Error ( "SaveManager", "0 bytes written saving manager to file" ); return kFALSE; }
333  return kTRUE;
334  }
335 
340 inline Double_t JetPtRho(const AliEmcalJet* jet, Double_t rho) {
341  return jet->Pt() - jet->Area() * rho;
342  }
343 
344 } // namespace AliAnalysisTaskEmcalJetCDF_NS
345 } // namespace EMCALJetTasks
346 } // namespace PWGJE
347 
348 
349 #endif // end of #ifndef ALIANALYSISTASKEMCALJETCDF_H
350 
351 // kate: indent-mode none; indent-width 2; replace-tabs on;
std::vector< Int_t > SortTracksPt(AliParticleContainer *track_container)
Bool_t fUseManualEventCuts
Flag to use manual event cuts.
THistManager fHistManager
Histogram manager.
Double_t Area() const
Definition: AliEmcalJet.h:130
Bool_t fUseAliEventCuts
Flag to use AliEventCuts (otherwise AliAnalysisTaskEmcal will be used)
double Double_t
Definition: External.C:58
Bool_t IsEventSelected()
Performing event selection.
std::pair< Double_t, Int_t > ptidx_pair
(pt,index) pair
void ExecOnce()
Perform steps needed to initialize the analysis.
Container for particles within the EMCAL framework.
QA Class for EMCal Embedding Framework.
Double_t Z_ptot(const AliEmcalJet *jet, const AliVParticle *trk)
Double_t Perp2(const AliVParticle &trk1, const AliVParticle &trk2)
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
Double_t DeltaR(const AliVParticle *part1, const AliVParticle *part2)
TObject * GetHistogram(const char *histName)
Double_t JetPtRho(const AliEmcalJet *jet, Double_t rho)
Double_t Z_pt(const AliEmcalJet *jet, const AliVParticle *trk)
AliEventCuts fEventCuts
event selection utility
functional for sorting pair by first element - descending
Double_t Perp(const AliVParticle &trk1, const AliVParticle &trk2)
Double_t Pt() const
Definition: AliEmcalJet.h:109
AliJetContainer * jetContSetParams(AliJetContainer *jetCont, Float_t jetptmin=1., Float_t jetptmax=500., Float_t jetareacutperc=0., Int_t leadhadtype=2, Int_t nLeadJets=1, Float_t mintrackpt=0.15, Float_t maxtrackpt=1000.)
AliMCParticleContainer * fGeneratorLevel
! generator level container
virtual ~AliAnalysisTaskEmcalJetCDF()
Destructor.
Double_t P() const
Definition: AliEmcalJet.h:110
AliAnalysisTaskEmcalJetCDF * AddTaskEmcalJetCDF(const char *ntracks="usedefault", const char *nclusters="usedefault", const char *ncells="usedefault", const char *ntracks_mc="", const char *tag="CDF")
TFile * file
TList with histograms for a given trigger.
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
Container class for histograms.
Definition: THistManager.h:99
AliEmcalEmbeddingQA fEmbeddingQA
! QA hists for embedding (will only be added if embedding)
TChain * CreateChain(const char *filelist="filelist.txt", const char *cTreeNameArg="auto", const char *friends="", UInt_t iNumFiles=-1, UInt_t iStartWithFile=1)
const char Option_t
Definition: External.C:48
Double_t Dot(const AliVParticle &trk1, const AliVParticle &trk2)
bool Bool_t
Definition: External.C:53
AliAnalysisTaskEmcalJetCDF & operator=(const AliAnalysisTaskEmcalJetCDF &)
Container for MC-true particles within the EMCAL framework.
Container for jet within the EMCAL jet framework.
Analysis of jet shapes and FF of all jets and leading jets.
TList * fEventCutList
! Output list for event cut histograms