AliPhysics  master (3d17d9d)
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 <TSystem.h>
13 #include <TFile.h>
14 #include <TChain.h>
15 #include <TH1D.h>
16 #include <TH2D.h>
17 #include <TArrayD.h>
18 #include <TVector2.h>
19 #include <TString.h>
20 #include <TPRegexp.h>
21 #include "THistManager.h"
22 #include "AliLog.h"
23 #include "AliAnalysisManager.h"
25 #include "AliEmcalEmbeddingQA.h"
26 #include "AliEmcalJetTask.h"
27 
34 public:
35 
37  AliAnalysisTaskEmcalJetCDF ( const char *name );
39 
41  void Terminate ( Option_t *option );
42 
45 
46 
47 
48 protected:
49  void ExecOnce();
50  Bool_t Run() ;
51 
55 
59  TObject* GetHistogram ( const char* histName );
60 
62 
64 
65  // Event selection
67  AliEventCuts fEventCuts;
70 
71  // MC options
73 
75 
76 
77 private:
78  AliAnalysisTaskEmcalJetCDF ( const AliAnalysisTaskEmcalJetCDF& ); // not implemented
80 
82  ClassDef ( AliAnalysisTaskEmcalJetCDF, 8 );
84 
85 };
86 
87 namespace PWGJE {
88 namespace EMCALJetTasks {
89 namespace AliAnalysisTaskEmcalJetCDF_NS {
90 
92  typedef std::pair<Double_t, Int_t> ptidx_pair;
93 
95  struct sort_descend
96  {
97  bool operator () ( const ptidx_pair &p1, const ptidx_pair &p2 ) { return p1.first > p2.first ; }
98  };
99 
103  inline Double_t Mag2 (const AliVParticle& trk)
104  { return trk.Px()*trk.Px() + trk.Py()*trk.Py() + trk.Pz()*trk.Pz(); }
105 
109  inline Double_t Mag(const AliVParticle& trk) { return TMath::Sqrt(Mag2(trk)); }
110 
115  inline Double_t Dot (const AliVParticle& trk1, const AliVParticle& trk2 )
116  { return trk1.Px()*trk2.Px() + trk1.Py()*trk2.Py() + trk1.Pz()*trk2.Pz(); }
117 
122  inline Double_t Perp2( const AliVParticle& trk1, const AliVParticle& trk2) {
123  Double_t mag1 = Mag2(trk1);
124  Double_t mag2 = Mag2(trk2);
125  Double_t dotp = Dot(trk1,trk2);
126  if (mag2 > 0.0) { mag1 -= dotp*dotp/mag2; }
127  if (mag1 <= 0) { mag1 = 0; }
128  return mag1;
129  }
130 
135  inline Double_t Perp (const AliVParticle& trk1, const AliVParticle& trk2 ) { return TMath::Sqrt(Perp2(trk1,trk2)); }
136 
140  std::vector<Int_t> SortTracksPt ( AliVEvent* event );
141 
145  std::vector<Int_t> SortTracksPt ( AliParticleContainer* track_container );
146 
151  inline Double_t Z_ptot ( const AliEmcalJet* jet, const AliVParticle* trk ) // Get Z of constituent trk ; p total
152  {
153  if (trk->P() < 1e-6) return 0.;
154  return (trk != 0) ? trk->P()/ jet->P() : 0.;
155  }
156 
161  inline Double_t Z_pt ( const AliEmcalJet* jet, const AliVParticle* trk ) // Get Z of constituent trk ; pt
162  {
163  if (trk->P() < 1e-6) return 0.;
164  return (trk != 0) ? trk->Pt() / jet->Pt() : 0.;
165  }
166 
169  inline Double_t Xi ( Double_t z ) { return TMath::Log ( 1/z ); } // Get Xi of value z
170 
175  inline Double_t DeltaR ( const AliVParticle* part1, const AliVParticle* part2 )
176  {
177  Double_t dPhi = part1->Phi() - part2->Phi();
178  Double_t dEta = part1->Eta() - part2->Eta();
179  dPhi = TVector2::Phi_mpi_pi ( dPhi );
180  return TMath::Sqrt ( dPhi * dPhi + dEta * dEta );
181  }
182 
187  inline Int_t IdxInArray ( Int_t index, TArrayI &array )
188  {
189  for ( Int_t i = 0; i < array.GetSize(); i++ ) { if ( index == array[i] ) { return i; } }
190  return -1;
191  }
192 
200  const char* ntracks = "usedefault",
201  const char* nclusters = "usedefault",
202  const char* ncells = "usedefault",
203  const char* ntracks_mc = "",
204  const char* tag = "CDF"
205  );
206 
218  AliJetContainer* jetCont,
219  Float_t jetptmin = 1.,
220  Float_t jetptmax = 500.,
221  Float_t jetareacutperc = 0.,
222  Int_t leadhadtype = 2,
223  Int_t nLeadJets = 1,
224  Float_t mintrackpt = 0.15,
225  Float_t maxtrackpt = 1000.
226  );
227 
228 
241 TChain* CreateChain ( const char* filelist = "filelist.txt",
242  const char* cTreeNameArg = "auto",
243  const char* friends = "",
244  UInt_t iNumFiles = -1,
245  UInt_t iStartWithFile = 1
246  );
247 
248 
249 
253 inline bool PeriodIsMC (const char* str) {
254  TString period (str);
255  if (!period.IsNull()) {
256  period.ToLower();
257  if ( period.BeginsWith("lhc") && (period.Length() > 6) ) { return true; }
258  }
259  return false;
260  }
261 
262 
266 inline TString GetPeriod (const char* file_path) {
267  TString period = "";
268  TString sFile(file_path);
269  sFile.ToLower();
270 
271  if (!sFile.IsNull()) {
272  // split string in tokens (libs)
273  TObjArray* tokens_list = sFile.Tokenize("/");
274  tokens_list->SetOwner(kTRUE);
275  TIter next_str(tokens_list);
276  TObjString* token = NULL;
277  while ((token=(TObjString*)next_str())) {
278  TString token_str = token->GetString();
279  if ( token_str.BeginsWith("lhc") ) { period = token_str; break; }
280  }
281  delete tokens_list;
282  }
283  return period;
284  }
285 
286 
290 inline TString GetPass ( const char* file_path) {
291  TString pass = "";
292  TString sFile (file_path);
293  sFile.ToLower();
294 
295  if (!sFile.IsNull()) {
296  // split string in tokens (libs)
297  TObjArray* tokens_list = sFile.Tokenize("/");
298  tokens_list->SetOwner(kTRUE);
299  TIter next_str(tokens_list);
300  TObjString* token = NULL;
301  while ((token=(TObjString*)next_str())) {
302  TString token_str = token->GetString();
303  if ( token_str.BeginsWith("pass") ) { pass = token_str; break; }
304  }
305  delete tokens_list;
306  }
307 
308  return pass;
309  }
310 
311 
315 inline TString GetFileFromPath ( const char* file_path = "" ) {
316 TString file (file_path);
317 TString file_name ("");
318 TObjArray* token_list = file.Tokenize("/");
319 token_list->SetOwner(kTRUE);
320 
321 if ( token_list->GetEntries() > 0 )
322  { file_name = ((TObjString*)token_list->At(token_list->GetLast()))->GetString(); }
323 
324 delete token_list;
325 return file_name;
326 }
327 
328 
332 inline bool SaveManager ( const char* file_name) {
333  AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
334  if ( !mgr ) { ::Error ( "SaveManager", "No analysis manager to connect to." ); return kFALSE; }
335 
336  TFile pOutFile (file_name,"RECREATE");
337  if ( ! pOutFile.cd() ) { ::Error ( "SaveManager", "Could not use the new created file" ); return kFALSE; }
338  Int_t written_bytes = mgr->Write();
339  pOutFile.Close();
340  if (written_bytes == 0 ) { ::Error ( "SaveManager", "0 bytes written saving manager to file" ); return kFALSE; }
341  return kTRUE;
342  }
343 
348 inline Double_t JetPtRho(const AliEmcalJet* jet, Double_t rho) {
349  return jet->Pt() - jet->Area() * rho;
350  }
351 
358 return task->AddJetContainer ( static_cast<AliAnalysisTaskEmcalJet::EJetType_t>(jf->GetJetType()),
359  static_cast<AliAnalysisTaskEmcalJet::EJetAlgo_t>(jf->GetJetAlgo()),
360  static_cast<AliAnalysisTaskEmcalJet::ERecoScheme_t>(jf->GetRecombScheme()),
361  jf->GetRadius(),
362  static_cast<UInt_t>(acc),
363  jf->GetJetsTag() );
364 }
365 
374 return task->AddJetContainer ( static_cast<AliAnalysisTaskEmcalJet::EJetType_t>(jf->GetJetType()),
375  static_cast<AliAnalysisTaskEmcalJet::EJetAlgo_t>(jf->GetJetAlgo()),
376  static_cast<AliAnalysisTaskEmcalJet::ERecoScheme_t>(jf->GetRecombScheme()),
377  jf->GetRadius(),
378  static_cast<UInt_t>(acc),
379  partcont,
380  cluscont,
381  jf->GetJetsTag() );
382 }
383 
386 inline void load_config(const char* file) {
387 std::ifstream filestream(file);
388 std::string str;
389 TPRegexp regex ("[a-zA-Z0-9_]+=[^\r^\n^\t^\f^\v^ ]+");
390 
391 while (std::getline(filestream, str)) {
392  TString line (str);
393  TString pair = line(regex);
394  if (!pair.IsNull()) {
395  TObjArray* decl = pair.Tokenize("=");
396  TString key = ((TObjString*)decl->At(0))->GetString();
397  TString value = ((TObjString*)decl->At(1))->GetString();
398  value = value.Strip(TString::EStripType::kBoth, '\"');
399  value = value.Strip(TString::EStripType::kBoth, '\'');
400  gSystem->Setenv(key.Data(), value.Data());
401  delete decl;
402  }
403  }
404 }
405 
406 } // namespace AliAnalysisTaskEmcalJetCDF_NS
407 } // namespace EMCALJetTasks
408 } // namespace PWGJE
409 
410 
411 #endif // end of #ifndef ALIANALYSISTASKEMCALJETCDF_H
412 
413 // 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
AliJetContainer * AddJetContainerJetTask(AliAnalysisTaskEmcalJet *task, AliEmcalJetTask *jf, AliEmcalJet::JetAcceptanceType acc)
Bool_t IsEventSelected()
Performing event selection.
std::pair< Double_t, Int_t > ptidx_pair
(pt,index) pair
JetAcceptanceType
Bit definition for jet geometry acceptance. Cut implemented in AliJetContainer by comparing jet&#39;s bit...
Definition: AliEmcalJet.h:66
TSystem * gSystem
void ExecOnce()
Perform steps needed to initialize the analysis.
AliJetContainer * AddJetContainer(const char *n, TString defaultCutType, Float_t jetRadius=0.4)
Double_t GetRadius()
General jet finder task implementing a wrapper for FastJet.
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)
AliJetContainer * AddJetContainerJetTaskCustomPartClus(AliAnalysisTaskEmcalJet *task, AliEmcalJetTask *jf, AliEmcalJet::JetAcceptanceType acc, AliParticleContainer *partcont, AliClusterContainer *cluscont)
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 structure for EMCAL clusters.
Container for MC-true particles within the EMCAL framework.
Container for jet within the EMCAL jet framework.
const char * GetJetsTag()
Analysis of jet shapes and FF of all jets and leading jets.
TList * fEventCutList
! Output list for event cut histograms