AliPhysics  b43479f (b43479f)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliFJWrapper.h
Go to the documentation of this file.
1 #ifndef AliFJWrapper_H
2 #define AliFJWrapper_H
3 
4 #if !defined(__CINT__)
5 
6 #include <vector>
7 #include <TString.h>
8 #include "AliLog.h"
9 #include "FJ_includes.h"
10 #include "AliJetShape.h"
11 
13 {
14  public:
15  AliFJWrapper(const char *name, const char *title);
16  virtual ~AliFJWrapper();
17 
18  virtual void AddInputVector (Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index = -99999);
19  virtual void AddInputVector (const fastjet::PseudoJet& vec, Int_t index = -99999);
20  virtual void AddInputVectors(const std::vector<fastjet::PseudoJet>& vecs, Int_t offsetIndex = -99999);
21  virtual void AddInputGhost (Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index = -99999);
22  virtual const char *ClassName() const { return "AliFJWrapper"; }
23  virtual void Clear(const Option_t* /*opt*/ = "");
24  virtual void ClearMemory();
25  virtual void CopySettingsFrom (const AliFJWrapper& wrapper);
26  virtual void GetMedianAndSigma(Double_t& median, Double_t& sigma, Int_t remove = 0) const;
27  fastjet::ClusterSequenceArea* GetClusterSequence() const { return fClustSeq; }
28  fastjet::ClusterSequence* GetClusterSequenceSA() const { return fClustSeqSA; }
29  fastjet::ClusterSequenceActiveAreaExplicitGhosts* GetClusterSequenceGhosts() const { return fClustSeqActGhosts; }
30  const std::vector<fastjet::PseudoJet>& GetInputVectors() const { return fInputVectors; }
31  const std::vector<fastjet::PseudoJet>& GetInputGhosts() const { return fInputGhosts; }
32  const std::vector<fastjet::PseudoJet>& GetInclusiveJets() const { return fInclusiveJets; }
33  const std::vector<fastjet::PseudoJet>& GetFilteredJets() const { return fFilteredJets; }
34  std::vector<fastjet::PseudoJet> GetJetConstituents(UInt_t idx) const;
35  std::vector<fastjet::PseudoJet> GetFilteredJetConstituents(UInt_t idx) const;
36  Double_t GetMedianUsedForBgSubtraction() const { return fMedUsedForBgSub; }
37  const char* GetName() const { return fName; }
38  const char* GetTitle() const { return fTitle; }
39  Double_t GetJetArea (UInt_t idx) const;
40  fastjet::PseudoJet GetJetAreaVector (UInt_t idx) const;
41  Double_t GetFilteredJetArea (UInt_t idx) const;
42  fastjet::PseudoJet GetFilteredJetAreaVector(UInt_t idx) const;
43  Double_t GetJetSubtractedPt (UInt_t idx) const;
44  virtual std::vector<double> GetSubtractedJetsPts(Double_t median_pt = -1, Bool_t sorted = kFALSE);
45  Bool_t GetLegacyMode() { return fLegacyMode; }
46  Bool_t GetDoFilterArea() { return fDoFilterArea; }
47 #ifdef FASTJET_VERSION
48  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetMass() const {return fGenSubtractorInfoJetMass ; }
49  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetAngularity() const {return fGenSubtractorInfoJetAngularity ; }
50  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetpTD() const {return fGenSubtractorInfoJetpTD ; }
51  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetCircularity() const {return fGenSubtractorInfoJetCircularity ; }
52  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetSigma2() const {return fGenSubtractorInfoJetSigma2 ; }
53  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetConstituent() const {return fGenSubtractorInfoJetConstituent ; }
54  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetLeSub() const {return fGenSubtractorInfoJetLeSub ; }
55  const std::vector<fastjet::PseudoJet> GetConstituentSubtrJets() const {return fConstituentSubtrJets ; }
56  Int_t CreateGenSub(); // fastjet::contrib::GenericSubtractor
57  Int_t CreateConstituentSub(); // fastjet::contrib::ConstituentSubtractor
58 #endif
59  virtual std::vector<double> GetGRNumerator() const { return fGRNumerator ; }
60  virtual std::vector<double> GetGRDenominator() const { return fGRDenominator ; }
61  virtual std::vector<double> GetGRNumeratorSub() const { return fGRNumeratorSub ; }
62  virtual std::vector<double> GetGRDenominatorSub() const { return fGRDenominatorSub ; }
63 
64  virtual void RemoveLastInputVector();
65 
66  virtual Int_t Run();
67  virtual Int_t Filter();
68  virtual Int_t DoGenericSubtractionJetMass();
69  virtual Int_t DoGenericSubtractionGR(Int_t ijet);
70  virtual Int_t DoGenericSubtractionJetAngularity();
71  virtual Int_t DoGenericSubtractionJetpTD();
72  virtual Int_t DoGenericSubtractionJetCircularity();
73  virtual Int_t DoGenericSubtractionJetSigma2();
74  virtual Int_t DoGenericSubtractionJetConstituent();
75  virtual Int_t DoGenericSubtractionJetLeSub();
76  virtual Int_t DoConstituentSubtraction();
77 
78  void SetName(const char* name) { fName = name; }
79  void SetTitle(const char* title) { fTitle = title; }
80  void SetStrategy(const fastjet::Strategy &strat) { fStrategy = strat; }
81  void SetAlgorithm(const fastjet::JetAlgorithm &algor) { fAlgor = algor; }
82  void SetRecombScheme(const fastjet::RecombinationScheme &scheme) { fScheme = scheme; }
83  void SetAreaType(const fastjet::AreaType &atype) { fAreaType = atype; }
84  void SetNRepeats(Int_t nrepeat) { fNGhostRepeats = nrepeat; }
85  void SetGhostArea(Double_t gharea) { fGhostArea = gharea; }
86  void SetMaxRap(Double_t maxrap) { fMaxRap = maxrap; }
87  void SetR(Double_t r) { fR = r; }
88  void SetGridScatter(Double_t gridSc) { fGridScatter = gridSc; }
89  void SetKtScatter(Double_t ktSc) { fKtScatter = ktSc; }
90  void SetMeanGhostKt(Double_t meankt) { fMeanGhostKt = meankt; }
91  void SetPluginAlgor(Int_t plugin) { fPluginAlgor = plugin; }
92  void SetUseArea4Vector(Bool_t useA4v) { fUseArea4Vector = useA4v; }
93  void SetupAlgorithmfromOpt(const char *option);
94  void SetupAreaTypefromOpt(const char *option);
95  void SetupSchemefromOpt(const char *option);
96  void SetupStrategyfromOpt(const char *option);
97  void SetLegacyMode (Bool_t mode) { fLegacyMode ^= mode; }
98  void SetLegacyFJ();
99  void SetUseExternalBkg(Bool_t b, Double_t rho, Double_t rhom) { fUseExternalBkg = b; fRho = rho; fRhom = rhom;}
100  void SetRMaxAndStep(Double_t rmax, Double_t dr) {fRMax = rmax; fDRStep = dr; }
101  void SetRhoRhom (Double_t rho, Double_t rhom) { fUseExternalBkg = kTRUE; fRho = rho; fRhom = rhom;} // if using rho,rhom then fUseExternalBkg is true
102 
103  protected:
104  TString fName;
105  TString fTitle;
106  std::vector<fastjet::PseudoJet> fInputVectors;
107  std::vector<fastjet::PseudoJet> fInputGhosts;
108  std::vector<fastjet::PseudoJet> fInclusiveJets;
109  std::vector<fastjet::PseudoJet> fFilteredJets;
110  std::vector<double> fSubtractedJetsPt;
111  std::vector<fastjet::PseudoJet> fConstituentSubtrJets;
112  fastjet::AreaDefinition *fAreaDef;
113  fastjet::VoronoiAreaSpec *fVorAreaSpec;
114  fastjet::GhostedAreaSpec *fGhostedAreaSpec;
115  fastjet::JetDefinition *fJetDef;
116  fastjet::JetDefinition::Plugin *fPlugin;
117 #ifndef FASTJET_VERSION
118  fastjet::RangeDefinition *fRange;
119 #else
120  fastjet::Selector *fRange;
121 #endif
122  fastjet::ClusterSequenceArea *fClustSeq;
123  fastjet::ClusterSequence *fClustSeqSA;
124  fastjet::ClusterSequenceActiveAreaExplicitGhosts *fClustSeqActGhosts;
125  fastjet::Strategy fStrategy;
126  fastjet::JetAlgorithm fAlgor;
127  fastjet::RecombinationScheme fScheme;
128  fastjet::AreaType fAreaType;
130  Double_t fGhostArea;
131  Double_t fMaxRap;
132  Double_t fR;
133  // no setters for the moment - used default values in the constructor
134  Double_t fGridScatter;
135  Double_t fKtScatter;
136  Double_t fMeanGhostKt;
137  Int_t fPluginAlgor;
138  // extra parameters
139  Double_t fMedUsedForBgSub;
141 #ifdef FASTJET_VERSION
142  fastjet::JetMedianBackgroundEstimator *fBkrdEstimator;
143  //from contrib package
144  fastjet::contrib::GenericSubtractor *fGenSubtractor;
145  fastjet::contrib::ConstituentSubtractor *fConstituentSubtractor;
146  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetMass;
147  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoGRNum;
148  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoGRDen;
149  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetAngularity;
150  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetpTD;
151  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetCircularity;
152  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetSigma2;
153  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetConstituent;
154  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetLeSub;
155 #endif
156  Bool_t fDoFilterArea;
157  Bool_t fLegacyMode;
159  Double_t fRho; // pT background density
160  Double_t fRhom; // mT background density
161  Double_t fRMax;
162  Double_t fDRStep;
163  std::vector<double> fGRNumerator;
164  std::vector<double> fGRDenominator;
165  std::vector<double> fGRNumeratorSub;
166  std::vector<double> fGRDenominatorSub;
167 
168  virtual void SubtractBackground(const Double_t median_pt = -1);
169 
170  private:
171  AliFJWrapper();
172  AliFJWrapper(const AliFJWrapper& wrapper);
173  AliFJWrapper& operator = (const AliFJWrapper& wrapper);
174 };
175 #endif
176 #endif
177 
178 #ifdef AliFJWrapper_CXX
179 #undef AliFJWrapper_CXX
180 
181 #if defined __GNUC__
182 #pragma GCC system_header
183 #endif
184 
185 namespace fj = fastjet;
186 
187 //_________________________________________________________________________________________________
188 AliFJWrapper::AliFJWrapper(const char *name, const char *title)
189  :
190  fName (name)
191  , fTitle (title)
192  , fInputVectors ( )
193  , fInputGhosts ( )
194  , fInclusiveJets ( )
195  , fFilteredJets ( )
196  , fSubtractedJetsPt ( )
197  , fConstituentSubtrJets ( )
198  , fAreaDef (0)
199  , fVorAreaSpec (0)
200  , fGhostedAreaSpec (0)
201  , fJetDef (0)
202  , fPlugin (0)
203  , fRange (0)
204  , fClustSeq (0)
205  , fClustSeqSA (0)
206  , fClustSeqActGhosts (0)
207  , fStrategy (fj::Best)
208  , fAlgor (fj::kt_algorithm)
209  , fScheme (fj::BIpt_scheme)
210  , fAreaType (fj::active_area)
211  , fNGhostRepeats (1)
212  , fGhostArea (0.005)
213  , fMaxRap (1.)
214  , fR (0.4)
215  , fGridScatter (1.0)
216  , fKtScatter (0.1)
217  , fMeanGhostKt (1e-100)
218  , fPluginAlgor (0)
219  , fMedUsedForBgSub (0)
220  , fUseArea4Vector (kFALSE)
221 #ifdef FASTJET_VERSION
222  , fBkrdEstimator (0)
223  , fGenSubtractor (0)
224  , fConstituentSubtractor (0)
225  , fGenSubtractorInfoJetMass ( )
226  , fGenSubtractorInfoGRNum ( )
227  , fGenSubtractorInfoGRDen ( )
228  , fGenSubtractorInfoJetAngularity ( )
229  , fGenSubtractorInfoJetpTD ( )
230  , fGenSubtractorInfoJetCircularity( )
231  , fGenSubtractorInfoJetSigma2()
232  , fGenSubtractorInfoJetConstituent ( )
233  , fGenSubtractorInfoJetLeSub ( )
234 #endif
235  , fDoFilterArea (false)
236  , fLegacyMode (false)
237  , fUseExternalBkg (false)
238  , fRho (0)
239  , fRhom (0)
240  , fRMax(2.)
241  , fDRStep(0.04)
242  , fGRNumerator()
243  , fGRDenominator()
244  , fGRNumeratorSub()
245  , fGRDenominatorSub()
246 {
247  // Constructor.
248 }
249 
250 //_________________________________________________________________________________________________
252 {
253  // Destructor.
254  ClearMemory();
255 }
256 
257 //_________________________________________________________________________________________________
259 {
260  // Destructor.
261  if (fAreaDef) { delete fAreaDef; fAreaDef = NULL; }
262  if (fVorAreaSpec) { delete fVorAreaSpec; fVorAreaSpec = NULL; }
263  if (fGhostedAreaSpec) { delete fGhostedAreaSpec; fGhostedAreaSpec = NULL; }
264  if (fJetDef) { delete fJetDef; fJetDef = NULL; }
265  if (fPlugin) { delete fPlugin; fPlugin = NULL; }
266  if (fRange) { delete fRange; fRange = NULL; }
267  if (fClustSeq) { delete fClustSeq; fClustSeq = NULL; }
268  if (fClustSeqSA) { delete fClustSeqSA; fClustSeqSA = NULL; }
270  #ifdef FASTJET_VERSION
271  if (fBkrdEstimator) { delete fBkrdEstimator; fBkrdEstimator = NULL; }
272  if (fGenSubtractor) { delete fGenSubtractor; fGenSubtractor = NULL; }
273  if (fConstituentSubtractor) { delete fConstituentSubtractor; fConstituentSubtractor = NULL; }
274  #endif
275 }
276 
277 //_________________________________________________________________________________________________
278 void AliFJWrapper::CopySettingsFrom(const AliFJWrapper& wrapper)
279 {
280  // Copy some settings.
281  // You very often want to keep most of the settings
282  // but change only the algorithm or R - do it after call to this function
283 
284  fStrategy = wrapper.fStrategy;
285  fAlgor = wrapper.fAlgor;
286  fScheme = wrapper.fScheme;
287  fAreaType = wrapper.fAreaType;
288  fNGhostRepeats = wrapper.fNGhostRepeats;
289  fGhostArea = wrapper.fGhostArea;
290  fMaxRap = wrapper.fMaxRap;
291  fR = wrapper.fR;
292  fGridScatter = wrapper.fGridScatter;
293  fKtScatter = wrapper.fKtScatter;
294  fMeanGhostKt = wrapper.fMeanGhostKt;
295  fPluginAlgor = wrapper.fPluginAlgor;
297  fLegacyMode = wrapper.fLegacyMode;
299  fRho = wrapper.fRho;
300  fRhom = wrapper.fRhom;
301 }
302 
303 //_________________________________________________________________________________________________
304 void AliFJWrapper::Clear(const Option_t */*opt*/)
305 {
306  // Simply clear the input vectors.
307  // Make sure done on every event if the instance is reused
308  // Reset the median to zero.
309 
310  fInputVectors.clear();
311  fInputGhosts.clear();
312  fMedUsedForBgSub = 0;
313 
314  // for the moment brute force delete everything
315  ClearMemory();
316 }
317 
318 //_________________________________________________________________________________________________
320 {
321  // Remove last input vector
322 
323  fInputVectors.pop_back();
324 }
325 
326 //_________________________________________________________________________________________________
327 void AliFJWrapper::AddInputVector(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index)
328 {
329  // Make the input pseudojet.
330 
331  fastjet::PseudoJet inVec(px, py, pz, E);
332 
333  // Salvatore Aiola: not sure why this was done...
334  //if (index > -99999) {
335  inVec.set_user_index(index);
336  //} else {
337  //inVec.set_user_index(fInputVectors.size());
338  //}
339 
340  // add to the fj container of input vectors
341  fInputVectors.push_back(inVec);
342 }
343 
344 //_________________________________________________________________________________________________
345 void AliFJWrapper::AddInputVector(const fj::PseudoJet& vec, Int_t index)
346 {
347  // Add an input pseudojet.
348 
349  fj::PseudoJet inVec = vec;
350 
351  // Salvatore Aiola: not sure why this was done...
353  inVec.set_user_index(index);
354  //} else {
355  //inVec.set_user_index(fInputVectors.size());
356  //}
357 
358  // add to the fj container of input vectors
359  fInputVectors.push_back(inVec);
360 }
361 
362 //_________________________________________________________________________________________________
363 void AliFJWrapper::AddInputVectors(const std::vector<fj::PseudoJet>& vecs, Int_t offsetIndex)
364 {
365  // Add the input from vector of pseudojets.
366 
367  for (UInt_t i = 0; i < vecs.size(); ++i) {
368  fj::PseudoJet inVec = vecs[i];
369  if (offsetIndex > -99999)
370  inVec.set_user_index(fInputVectors.size() + offsetIndex);
371  // add to the fj container of input vectors
372  fInputVectors.push_back(inVec);
373  }
374 }
375 
376 //_________________________________________________________________________________________________
377 void AliFJWrapper::AddInputGhost(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index)
378 {
379  // Make the input pseudojet.
380 
381  fastjet::PseudoJet inVec(px, py, pz, E);
382 
383  if (index > -99999) {
384  inVec.set_user_index(index);
385  } else {
386  inVec.set_user_index(fInputGhosts.size());
387  }
388 
389  // add to the fj container of input vectors
390  fInputGhosts.push_back(inVec);
391  if (!fDoFilterArea) fDoFilterArea = kTRUE;
392 }
393 
394 //_________________________________________________________________________________________________
395 Double_t AliFJWrapper::GetJetArea(UInt_t idx) const
396 {
397  // Get the jet area.
398 
399  Double_t retval = -1; // really wrong area..
400  if ( idx < fInclusiveJets.size() ) {
401  retval = fClustSeq->area(fInclusiveJets[idx]);
402  } else {
403  AliError(Form("[e] ::GetJetArea wrong index: %d",idx));
404  }
405  return retval;
406 }
407 
408 //_________________________________________________________________________________________________
409 Double_t AliFJWrapper::GetFilteredJetArea(UInt_t idx) const
410 {
411  // Get the filtered jet area.
412 
413  Double_t retval = -1; // really wrong area..
414  if (fDoFilterArea && fClustSeqActGhosts && (idx<fFilteredJets.size())) {
415  retval = fClustSeqActGhosts->area(fFilteredJets[idx]);
416  } else {
417  AliError(Form("[e] ::GetFilteredJetArea wrong index: %d",idx));
418  }
419  return retval;
420 }
421 
422 //_________________________________________________________________________________________________
423 fastjet::PseudoJet AliFJWrapper::GetJetAreaVector(UInt_t idx) const
424 {
425  // Get the jet area as vector.
426  fastjet::PseudoJet retval;
427  if ( idx < fInclusiveJets.size() ) {
428  retval = fClustSeq->area_4vector(fInclusiveJets[idx]);
429  } else {
430  AliError(Form("[e] ::GetJetArea wrong index: %d",idx));
431  }
432  return retval;
433 }
434 
435 //_________________________________________________________________________________________________
436 fastjet::PseudoJet AliFJWrapper::GetFilteredJetAreaVector(UInt_t idx) const
437 {
438  // Get the jet area as vector.
439  fastjet::PseudoJet retval;
440  if (fDoFilterArea && fClustSeqActGhosts && (idx<fFilteredJets.size())) {
441  retval = fClustSeqActGhosts->area_4vector(fFilteredJets[idx]);
442  } else {
443  AliError(Form("[e] ::GetFilteredJetArea wrong index: %d",idx));
444  }
445  return retval;
446 }
447 
448 //_________________________________________________________________________________________________
449 std::vector<double> AliFJWrapper::GetSubtractedJetsPts(Double_t median_pt, Bool_t sorted)
450 {
451  // Get subtracted jets pTs, returns vector.
452 
453  SubtractBackground(median_pt);
454 
455  if (kTRUE == sorted) {
456  std::sort(fSubtractedJetsPt.begin(), fSubtractedJetsPt.begin());
457  }
458  return fSubtractedJetsPt;
459 }
460 
461 //_________________________________________________________________________________________________
462 Double_t AliFJWrapper::GetJetSubtractedPt(UInt_t idx) const
463 {
464  // Get subtracted jets pTs, returns Double_t.
465 
466  Double_t retval = -99999.; // really wrong pt..
467  if ( idx < fSubtractedJetsPt.size() ) {
468  retval = fSubtractedJetsPt[idx];
469  }
470  return retval;
471 }
472 
473 //_________________________________________________________________________________________________
474 std::vector<fastjet::PseudoJet>
475 AliFJWrapper::GetJetConstituents(UInt_t idx) const
476 {
477  // Get jets constituents.
478 
479  std::vector<fastjet::PseudoJet> retval;
480 
481  if ( idx < fInclusiveJets.size() ) {
482  retval = fClustSeq->constituents(fInclusiveJets[idx]);
483  } else {
484  AliError(Form("[e] ::GetJetConstituents wrong index: %d",idx));
485  }
486 
487  return retval;
488 }
489 
490 //_________________________________________________________________________________________________
491 std::vector<fastjet::PseudoJet>
493 {
494  // Get jets constituents.
495 
496  std::vector<fastjet::PseudoJet> retval;
497 
498  if ( idx < fFilteredJets.size() ) {
499  if (fClustSeqSA) retval = fClustSeqSA->constituents(fFilteredJets[idx]);
500  if (fClustSeqActGhosts) retval = fClustSeqActGhosts->constituents(fFilteredJets[idx]);
501  } else {
502  AliError(Form("[e] ::GetFilteredJetConstituents wrong index: %d",idx));
503  }
504 
505  return retval;
506 }
507 
508 //_________________________________________________________________________________________________
509 void AliFJWrapper::GetMedianAndSigma(Double_t &median, Double_t &sigma, Int_t remove) const
510 {
511  // Get the median and sigma from fastjet.
512  // User can also do it on his own because the cluster sequence is exposed (via a getter)
513 
514  if (!fClustSeq) {
515  AliError("[e] Run the jfinder first.");
516  return;
517  }
518 
519  Double_t mean_area = 0;
520  try {
521  if(0 == remove) {
522  fClustSeq->get_median_rho_and_sigma(*fRange, fUseArea4Vector, median, sigma, mean_area);
523  } else {
524  std::vector<fastjet::PseudoJet> input_jets = sorted_by_pt(fClustSeq->inclusive_jets());
525  input_jets.erase(input_jets.begin(), input_jets.begin() + remove);
526  fClustSeq->get_median_rho_and_sigma(input_jets, *fRange, fUseArea4Vector, median, sigma, mean_area);
527  input_jets.clear();
528  }
529  } catch (fj::Error) {
530  AliError(" [w] FJ Exception caught.");
531  median = -1.;
532  sigma = -1;
533  }
534 }
535 
536 //_________________________________________________________________________________________________
537 Int_t AliFJWrapper::Run()
538 {
539  // Run the actual jet finder.
540 
541  if (fAreaType == fj::voronoi_area) {
542  // Rfact - check dependence - default is 1.
543  // NOTE: hardcoded variable!
544  fVorAreaSpec = new fj::VoronoiAreaSpec(1.);
545  fAreaDef = new fj::AreaDefinition(*fVorAreaSpec);
546  } else {
547  fGhostedAreaSpec = new fj::GhostedAreaSpec(fMaxRap,
549  fGhostArea,
550  fGridScatter,
551  fKtScatter,
552  fMeanGhostKt);
553 
554  fAreaDef = new fj::AreaDefinition(*fGhostedAreaSpec, fAreaType);
555  }
556 
557  // this is acceptable by fastjet:
558 #ifndef FASTJET_VERSION
559  fRange = new fj::RangeDefinition(fMaxRap - 0.95 * fR);
560 #else
561  fRange = new fj::Selector(fj::SelectorAbsRapMax(fMaxRap - 0.95 * fR));
562 #endif
563 
564  if (fAlgor == fj::plugin_algorithm) {
565  if (fPluginAlgor == 0) {
566  // SIS CONE ALGOR
567  // NOTE: hardcoded split parameter
568  Double_t overlap_threshold = 0.75; // NOTE: this actually splits a lot: thr/min(pt1,pt2)
569  fPlugin = new fj::SISConePlugin(fR,
570  overlap_threshold,
571  0, //search of stable cones - zero = until no more
572  1.0); // this should be seed effectively for proto jets
573  fJetDef = new fastjet::JetDefinition(fPlugin);
574  } else if (fPluginAlgor == 1) {
575  // CDF cone
576  // NOTE: hardcoded split parameter
577  Double_t overlap_threshold = 0.75; // NOTE: this actually splits a lot: thr/min(pt1,pt2)
578  fPlugin = new fj::CDFMidPointPlugin(fR,
579  overlap_threshold,
580  1.0, //search of stable cones - zero = until no more
581  1.0); // this should be seed effectively for proto jets
582  fJetDef = new fastjet::JetDefinition(fPlugin);
583  } else {
584  AliError("[e] Unrecognized plugin number!");
585  }
586  } else {
587  fJetDef = new fj::JetDefinition(fAlgor, fR, fScheme, fStrategy);
588  }
589 
590  try {
591  fClustSeq = new fj::ClusterSequenceArea(fInputVectors, *fJetDef, *fAreaDef);
592  } catch (fj::Error) {
593  AliError(" [w] FJ Exception caught.");
594  return -1;
595  }
596 
597  // FJ3 :: Define an JetMedianBackgroundEstimator just in case it will be used
598 #ifdef FASTJET_VERSION
599  fBkrdEstimator = new fj::JetMedianBackgroundEstimator(fj::SelectorAbsRapMax(fMaxRap));
600 #endif
601 
602  if (fLegacyMode) { SetLegacyFJ(); } // for FJ 2.x even if fLegacyMode is set, SetLegacyFJ is dummy
603 
604  // inclusive jets:
605  fInclusiveJets.clear();
606  fInclusiveJets = fClustSeq->inclusive_jets(0.0);
607 
608  return 0;
609 }
610 
611 //_________________________________________________________________________________________________
612 Int_t AliFJWrapper::Filter()
613 {
614 //
615 // AliFJWrapper::Filter
616 //
617 
618  fJetDef = new fj::JetDefinition(fAlgor, fR, fScheme, fStrategy);
619 
620  if (fDoFilterArea) {
621  if (fInputGhosts.size()>0) {
622  try {
623  fClustSeqActGhosts = new fj::ClusterSequenceActiveAreaExplicitGhosts(fInputVectors,
624  *fJetDef,
625  fInputGhosts,
626  fGhostArea);
627  } catch (fj::Error) {
628  AliError(" [w] FJ Exception caught.");
629  return -1;
630  }
631 
632  fFilteredJets.clear();
633  fFilteredJets = fClustSeqActGhosts->inclusive_jets(0.0);
634  } else {
635  return -1;
636  }
637  } else {
638  try {
639  fClustSeqSA = new fastjet::ClusterSequence(fInputVectors, *fJetDef);
640  } catch (fj::Error) {
641  AliError(" [w] FJ Exception caught.");
642  return -1;
643  }
644 
645  fFilteredJets.clear();
646  fFilteredJets = fClustSeqSA->inclusive_jets(0.0);
647  }
648 
649  return 0;
650 }
651 
652 //_________________________________________________________________________________________________
654 {
655  // This methods enable legacy behaviour (FastJet 2.x) when AliROOT is compiled with FastJet 3.x
656 #ifdef FASTJET_VERSION
657  std::cout << "WARNING! Setting FastJet in legacy mode" << std::endl;
658  if (fGhostedAreaSpec) { fGhostedAreaSpec->set_fj2_placement(kTRUE); }
659  if (fBkrdEstimator) {
660  fBkrdEstimator->set_provide_fj2_sigma(kTRUE);
661  fBkrdEstimator->set_use_area_4vector(kFALSE);
662  }
663 #endif
664 }
665 
666 //_________________________________________________________________________________________________
667 void AliFJWrapper::SubtractBackground(Double_t median_pt)
668 {
669  // Subtract the background (specify the value - see below the meaning).
670  // Negative argument means the bg will be determined with the current algorithm
671  // this is the default behavior. Zero is allowed
672  // Note: user may set the switch for area4vector based subtraction.
673 
674  Double_t median = 0;
675  Double_t sigma = 0;
676  Double_t mean_area = 0;
677 
678  // clear the subtracted jet pt's vector<double>
679  fSubtractedJetsPt.clear();
680 
681  // check what was specified (default is -1)
682  if (median_pt < 0) {
683  try {
684  fClustSeq->get_median_rho_and_sigma(*fRange, fUseArea4Vector, median, sigma, mean_area);
685  }
686 
687  catch (fj::Error) {
688  AliError(" [w] FJ Exception caught.");
689  median = -9999.;
690  sigma = -1;
691  fMedUsedForBgSub = median;
692  return;
693  }
694  fMedUsedForBgSub = median;
695  } else {
696  // we do not know the sigma in this case
697  sigma = -1;
698  if (0.0 == median_pt) {
699  // AliWarning(Form(" [w] Median specified for bg subtraction is ZERO: %f", median_pt )) << endl;
700  fMedUsedForBgSub = 0.;
701  } else {
702  fMedUsedForBgSub = median_pt;
703  }
704  }
705 
706  // subtract:
707  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
708  if ( fUseArea4Vector ) {
709  // subtract the background using the area4vector
710  fj::PseudoJet area4v = fClustSeq->area_4vector(fInclusiveJets[i]);
711  fj::PseudoJet jet_sub = fInclusiveJets[i] - area4v * fMedUsedForBgSub;
712  fSubtractedJetsPt.push_back(jet_sub.perp()); // here we put only the pt of the jet - note: this can be negative
713  } else {
714  // subtract the background using scalars
715  // fj::PseudoJet jet_sub = fInclusiveJets[i] - area * fMedUsedForBgSub_;
716  Double_t area = fClustSeq->area(fInclusiveJets[i]);
717  // standard subtraction
718  Double_t pt_sub = fInclusiveJets[i].perp() - fMedUsedForBgSub * area;
719  fSubtractedJetsPt.push_back(pt_sub); // here we put only the pt of the jet - note: this can be negative
720  }
721  }
722 }
723 
724 //_________________________________________________________________________________________________
726  //Do generic subtraction for jet mass
727 #ifdef FASTJET_VERSION
728  CreateGenSub();
729 
730  // Define jet shape
731  AliJetShapeMass shapeMass;
732 
733  // clear the generic subtractor info vector
734  fGenSubtractorInfoJetMass.clear();
735  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
736  fj::contrib::GenericSubtractorInfo info;
737  if(fInclusiveJets[i].perp()>1.e-4)
738  double subtracted_shape = (*fGenSubtractor)(shapeMass, fInclusiveJets[i], info);
739  fGenSubtractorInfoJetMass.push_back(info);
740  }
741 #endif
742  return 0;
743 }
744 
745 //_________________________________________________________________________________________________
746 Int_t AliFJWrapper::DoGenericSubtractionGR(Int_t ijet) {
747  //Do generic subtraction for jet mass
748 #ifdef FASTJET_VERSION
749  CreateGenSub();
750 
751  if(ijet>fInclusiveJets.size()) return 0;
752 
753  fGRNumerator.clear();
754  fGRDenominator.clear();
755  fGRNumeratorSub.clear();
756  fGRDenominatorSub.clear();
757 
758  // Define jet shape
759  for(Double_t r = 0.; r<fRMax; r+=fDRStep) {
760  AliJetShapeGRNum shapeGRNum(r,fDRStep);
761  AliJetShapeGRDen shapeGRDen(r,fDRStep);
762 
763  // clear the generic subtractor info vector
764  fGenSubtractorInfoGRNum.clear();
765  fGenSubtractorInfoGRDen.clear();
766  fj::contrib::GenericSubtractorInfo infoNum;
767  fj::contrib::GenericSubtractorInfo infoDen;
768  if(fInclusiveJets[ijet].perp()>1.e-4) {
769  double sub_num = (*fGenSubtractor)(shapeGRNum, fInclusiveJets[ijet], infoNum);
770  double sub_den = (*fGenSubtractor)(shapeGRDen, fInclusiveJets[ijet], infoDen);
771  }
772  fGenSubtractorInfoGRNum.push_back(infoNum);
773  fGenSubtractorInfoGRDen.push_back(infoDen);
774  fGRNumerator.push_back(infoNum.unsubtracted());
775  fGRDenominator.push_back(infoDen.unsubtracted());
776  fGRNumeratorSub.push_back(infoNum.second_order_subtracted());
777  fGRDenominatorSub.push_back(infoDen.second_order_subtracted());
778  }
779 #endif
780  return 0;
781 }
782 //_________________________________________________________________________________________________
784  //Do generic subtraction for jet mass
785 #ifdef FASTJET_VERSION
786  CreateGenSub();
787 
788  // Define jet shape
789  AliJetShapeAngularity shapeAngularity;
790 
791  // clear the generic subtractor info vector
792  fGenSubtractorInfoJetAngularity.clear();
793  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
794  fj::contrib::GenericSubtractorInfo infoAng;
795  if(fInclusiveJets[i].perp()>1.e-4)
796  double subtracted_shape = (*fGenSubtractor)(shapeAngularity, fInclusiveJets[i], infoAng);
797  fGenSubtractorInfoJetAngularity.push_back(infoAng);
798  }
799 #endif
800  return 0;
801 }
802 //_________________________________________________________________________________________________
804  //Do generic subtraction for jet mass
805 #ifdef FASTJET_VERSION
806  CreateGenSub();
807 
808  // Define jet shape
809  AliJetShapepTD shapepTD;
810 
811  // clear the generic subtractor info vector
812  fGenSubtractorInfoJetpTD.clear();
813  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
814  fj::contrib::GenericSubtractorInfo infopTD;
815  if(fInclusiveJets[i].perp()>1.e-4)
816  double subtracted_shape = (*fGenSubtractor)(shapepTD, fInclusiveJets[i], infopTD);
817  fGenSubtractorInfoJetpTD.push_back(infopTD);
818  }
819 #endif
820  return 0;
821 }
822 //_________________________________________________________________________________________________
824  //Do generic subtraction for jet mass
825 #ifdef FASTJET_VERSION
826  CreateGenSub();
827 
828  // Define jet shape
829  AliJetShapeCircularity shapecircularity;
830 
831  // clear the generic subtractor info vector
832  fGenSubtractorInfoJetCircularity.clear();
833  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
834  fj::contrib::GenericSubtractorInfo infoCirc;
835  if(fInclusiveJets[i].perp()>1.e-4)
836  double subtracted_shape = (*fGenSubtractor)(shapecircularity, fInclusiveJets[i], infoCirc);
837  fGenSubtractorInfoJetCircularity.push_back(infoCirc);
838  }
839 #endif
840  return 0;
841 }
842 //_________________________________________________________________________________________________
844  //Do generic subtraction for jet mass
845 #ifdef FASTJET_VERSION
846  CreateGenSub();
847 
848  // Define jet shape
849  AliJetShapeSigma2 shapesigma2;
850 
851  // clear the generic subtractor info vector
852  fGenSubtractorInfoJetSigma2.clear();
853  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
854  fj::contrib::GenericSubtractorInfo infoSigma;
855  if(fInclusiveJets[i].perp()>1.e-4)
856  double subtracted_shape = (*fGenSubtractor)(shapesigma2, fInclusiveJets[i], infoSigma);
857  fGenSubtractorInfoJetSigma2.push_back(infoSigma);
858  }
859 #endif
860  return 0;
861 }
862 //_________________________________________________________________________________________________
864  //Do generic subtraction for jet mass
865 #ifdef FASTJET_VERSION
866  CreateGenSub();
867 
868  // Define jet shape
869  AliJetShapeConstituent shapeconst;
870 
871  // clear the generic subtractor info vector
872  fGenSubtractorInfoJetConstituent.clear();
873  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
874  fj::contrib::GenericSubtractorInfo infoConst;
875  if(fInclusiveJets[i].perp()>1.e-4)
876  double subtracted_shape = (*fGenSubtractor)(shapeconst, fInclusiveJets[i], infoConst);
877  fGenSubtractorInfoJetConstituent.push_back(infoConst);
878  }
879 #endif
880  return 0;
881 }
882 
883 //_________________________________________________________________________________________________
885  //Do generic subtraction for jet mass
886 #ifdef FASTJET_VERSION
887  CreateGenSub();
888 
889  // Define jet shape
890  AliJetShapeLeSub shapeLeSub;
891 
892  // clear the generic subtractor info vector
893  fGenSubtractorInfoJetLeSub.clear();
894  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
895  fj::contrib::GenericSubtractorInfo infoLeSub;
896  if(fInclusiveJets[i].perp()>1.e-4)
897  double subtracted_shape = (*fGenSubtractor)(shapeLeSub, fInclusiveJets[i], infoLeSub);
898  fGenSubtractorInfoJetLeSub.push_back(infoLeSub);
899  }
900 #endif
901  return 0;
902 }
903 
904 //_________________________________________________________________________________________________
906  //Do constituent subtraction
907 #ifdef FASTJET_VERSION
908  CreateConstituentSub();
909  // fConstituentSubtractor->set_alpha(/* double alpha */);
910  // fConstituentSubtractor->set_max_deltaR(/* double max_deltaR */);
911 
912  //clear constituent subtracted jets
913  fConstituentSubtrJets.clear();
914  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
915  fj::PseudoJet subtracted_jet(0.,0.,0.,0.);
916  if(fInclusiveJets[i].perp()>0.)
917  subtracted_jet = (*fConstituentSubtractor)(fInclusiveJets[i]);
918  fConstituentSubtrJets.push_back(subtracted_jet);
919  }
920  if(fConstituentSubtractor) { delete fConstituentSubtractor; fConstituentSubtractor = NULL; }
921 
922 #endif
923  return 0;
924 }
925 
926 //_________________________________________________________________________________________________
927 Int_t AliFJWrapper::CreateGenSub() {
928  //Do generic subtraction for jet mass
929  #ifdef FASTJET_VERSION
930  if (fGenSubtractor) { delete fGenSubtractor; } // protect against memory leaks
931 
932  if (fUseExternalBkg)
933  { fGenSubtractor = new fj::contrib::GenericSubtractor(fRho,fRhom); }
934  else
935  {
936  fGenSubtractor = new fj::contrib::GenericSubtractor(fBkrdEstimator);
937  #if FASTJET_VERSION_NUMBER >= 30100
938  fGenSubtractor->set_common_bge_for_rho_and_rhom(); // see contrib 1.020 GenericSubtractor.hh line 62
939  #endif
940  }
941 
942  #endif
943  return 0;
944 }
945 
946 //_________________________________________________________________________________________________
947 Int_t AliFJWrapper::CreateConstituentSub() {
948  //Do generic subtraction for jet mass
949  #ifdef FASTJET_VERSION
950  if (fConstituentSubtractor) { delete fConstituentSubtractor; } // protect against memory leaks
951 
952  // see ConstituentSubtractor.hh signatures
953  // ConstituentSubtractor(double rho, double rhom=0, double alpha=0, double maxDeltaR=-1)
954  if (fUseExternalBkg) { fConstituentSubtractor = new fj::contrib::ConstituentSubtractor(fRho,fRhom); }
955  else { fConstituentSubtractor = new fj::contrib::ConstituentSubtractor(fBkrdEstimator); }
956 
957  #endif
958  return 0;
959 }
960 
961 //_________________________________________________________________________________________________
962 void AliFJWrapper::SetupAlgorithmfromOpt(const char *option)
963 {
964  // Setup algorithm from char.
965 
966  std::string opt(option);
967 
968  if (!opt.compare("kt")) fAlgor = fj::kt_algorithm;
969  if (!opt.compare("antikt")) fAlgor = fj::antikt_algorithm;
970  if (!opt.compare("cambridge")) fAlgor = fj::cambridge_algorithm;
971  if (!opt.compare("genkt")) fAlgor = fj::genkt_algorithm;
972  if (!opt.compare("cambridge_passive")) fAlgor = fj::cambridge_for_passive_algorithm;
973  if (!opt.compare("genkt_passive")) fAlgor = fj::genkt_for_passive_algorithm;
974  if (!opt.compare("ee_kt")) fAlgor = fj::ee_kt_algorithm;
975  if (!opt.compare("ee_genkt")) fAlgor = fj::ee_genkt_algorithm;
976  if (!opt.compare("plugin")) fAlgor = fj::plugin_algorithm;
977 }
978 
979 //_________________________________________________________________________________________________
980 void AliFJWrapper::SetupAreaTypefromOpt(const char *option)
981 {
982  // Setup area type from char.
983 
984  std::string opt(option);
985 
986  if (!opt.compare("active")) fAreaType = fj::active_area;
987  if (!opt.compare("invalid")) fAreaType = fj::invalid_area;
988  if (!opt.compare("active_area_explicit_ghosts")) fAreaType = fj::active_area_explicit_ghosts;
989  if (!opt.compare("one_ghost_passive")) fAreaType = fj::one_ghost_passive_area;
990  if (!opt.compare("passive")) fAreaType = fj::passive_area;
991  if (!opt.compare("voronoi")) fAreaType = fj::voronoi_area;
992 }
993 
994 //_________________________________________________________________________________________________
995 void AliFJWrapper::SetupSchemefromOpt(const char *option)
996 {
997  //
998  // setup scheme from char
999  //
1000 
1001  std::string opt(option);
1002 
1003  if (!opt.compare("BIpt")) fScheme = fj::BIpt_scheme;
1004  if (!opt.compare("BIpt2")) fScheme = fj::BIpt2_scheme;
1005  if (!opt.compare("E")) fScheme = fj::E_scheme;
1006  if (!opt.compare("pt")) fScheme = fj::pt_scheme;
1007  if (!opt.compare("pt2")) fScheme = fj::pt2_scheme;
1008  if (!opt.compare("Et")) fScheme = fj::Et_scheme;
1009  if (!opt.compare("Et2")) fScheme = fj::Et2_scheme;
1010 }
1011 
1012 //_________________________________________________________________________________________________
1013 void AliFJWrapper::SetupStrategyfromOpt(const char *option)
1014 {
1015  // Setup strategy from char.
1016 
1017  std::string opt(option);
1018 
1019  if (!opt.compare("Best")) fStrategy = fj::Best;
1020  if (!opt.compare("N2MinHeapTiled")) fStrategy = fj::N2MinHeapTiled;
1021  if (!opt.compare("N2Tiled")) fStrategy = fj::N2Tiled;
1022  if (!opt.compare("N2PoorTiled")) fStrategy = fj::N2PoorTiled;
1023  if (!opt.compare("N2Plain")) fStrategy = fj::N2Plain;
1024  if (!opt.compare("N3Dumb")) fStrategy = fj::N3Dumb;
1025  if (!opt.compare("NlnN")) fStrategy = fj::NlnN;
1026  if (!opt.compare("NlnN3pi")) fStrategy = fj::NlnN3pi;
1027  if (!opt.compare("NlnN4pi")) fStrategy = fj::NlnN4pi;
1028  if (!opt.compare("NlnNCam4pi")) fStrategy = fj::NlnNCam4pi;
1029  if (!opt.compare("NlnNCam2pi2R")) fStrategy = fj::NlnNCam2pi2R;
1030  if (!opt.compare("NlnNCam")) fStrategy = fj::NlnNCam;
1031  if (!opt.compare("plugin")) fStrategy = fj::plugin_strategy;
1032 }
1033 #endif
fastjet::PseudoJet GetJetAreaVector(UInt_t idx) const
std::vector< double > fGRNumeratorSub
Definition: AliFJWrapper.h:165
void SetRMaxAndStep(Double_t rmax, Double_t dr)
Definition: AliFJWrapper.h:100
void SetUseExternalBkg(Bool_t b, Double_t rho, Double_t rhom)
Definition: AliFJWrapper.h:99
Double_t GetFilteredJetArea(UInt_t idx) const
std::vector< fastjet::PseudoJet > fInputVectors
Definition: AliFJWrapper.h:106
Double_t fR
Definition: AliFJWrapper.h:132
virtual Int_t Run()
virtual Int_t Filter()
const char * title
Definition: MakeQAPdf.C:26
virtual std::vector< double > GetSubtractedJetsPts(Double_t median_pt=-1, Bool_t sorted=kFALSE)
fastjet::ClusterSequence * fClustSeqSA
Definition: AliFJWrapper.h:123
Int_t fPluginAlgor
Definition: AliFJWrapper.h:137
fastjet::ClusterSequenceActiveAreaExplicitGhosts * fClustSeqActGhosts
Definition: AliFJWrapper.h:124
Double_t fKtScatter
Definition: AliFJWrapper.h:135
Int_t fNGhostRepeats
Definition: AliFJWrapper.h:129
virtual Int_t DoGenericSubtractionJetpTD()
Double_t GetMedianUsedForBgSubtraction() const
Definition: AliFJWrapper.h:36
virtual void AddInputVectors(const std::vector< fastjet::PseudoJet > &vecs, Int_t offsetIndex=-99999)
std::vector< fastjet::PseudoJet > fInclusiveJets
Definition: AliFJWrapper.h:108
fastjet::JetAlgorithm fAlgor
Definition: AliFJWrapper.h:126
const char * GetTitle() const
Definition: AliFJWrapper.h:38
fastjet::RecombinationScheme fScheme
Definition: AliFJWrapper.h:127
Double_t fRMax
Definition: AliFJWrapper.h:161
virtual void CopySettingsFrom(const AliFJWrapper &wrapper)
void SetupAlgorithmfromOpt(const char *option)
virtual Int_t DoGenericSubtractionJetAngularity()
virtual Int_t DoGenericSubtractionGR(Int_t ijet)
virtual Int_t DoGenericSubtractionJetMass()
void SetMaxRap(Double_t maxrap)
Definition: AliFJWrapper.h:86
Double_t GetJetSubtractedPt(UInt_t idx) const
fastjet::PseudoJet GetFilteredJetAreaVector(UInt_t idx) const
void SetRecombScheme(const fastjet::RecombinationScheme &scheme)
Definition: AliFJWrapper.h:82
fastjet::RangeDefinition * fRange
Definition: AliFJWrapper.h:118
fastjet::ClusterSequenceActiveAreaExplicitGhosts * GetClusterSequenceGhosts() const
Definition: AliFJWrapper.h:29
fastjet::ClusterSequence * GetClusterSequenceSA() const
Definition: AliFJWrapper.h:28
void SetStrategy(const fastjet::Strategy &strat)
Definition: AliFJWrapper.h:80
std::vector< fastjet::PseudoJet > fFilteredJets
Definition: AliFJWrapper.h:109
void SetupAreaTypefromOpt(const char *option)
virtual void ClearMemory()
fastjet::AreaType fAreaType
Definition: AliFJWrapper.h:128
void SetR(Double_t r)
Definition: AliFJWrapper.h:87
Double_t fMaxRap
Definition: AliFJWrapper.h:131
fastjet::ClusterSequenceArea * GetClusterSequence() const
Definition: AliFJWrapper.h:27
Double_t * sigma
const std::vector< fastjet::PseudoJet > & GetInclusiveJets() const
Definition: AliFJWrapper.h:32
AliFJWrapper & operator=(const AliFJWrapper &wrapper)
Double_t GetJetArea(UInt_t idx) const
std::vector< double > fSubtractedJetsPt
Definition: AliFJWrapper.h:110
fastjet::JetDefinition::Plugin * fPlugin
Definition: AliFJWrapper.h:116
void SetLegacyMode(Bool_t mode)
Definition: AliFJWrapper.h:97
Bool_t fDoFilterArea
Definition: AliFJWrapper.h:156
void SetAlgorithm(const fastjet::JetAlgorithm &algor)
Definition: AliFJWrapper.h:81
virtual Int_t DoGenericSubtractionJetCircularity()
void SetupSchemefromOpt(const char *option)
const std::vector< fastjet::PseudoJet > & GetInputGhosts() const
Definition: AliFJWrapper.h:31
void SetUseArea4Vector(Bool_t useA4v)
Definition: AliFJWrapper.h:92
fastjet::ClusterSequenceArea * fClustSeq
Definition: AliFJWrapper.h:122
void SetName(const char *name)
Definition: AliFJWrapper.h:78
virtual void SubtractBackground(const Double_t median_pt=-1)
std::vector< fastjet::PseudoJet > fInputGhosts
Definition: AliFJWrapper.h:107
Double_t fGridScatter
Definition: AliFJWrapper.h:134
virtual void Clear(const Option_t *="")
Bool_t fLegacyMode
Definition: AliFJWrapper.h:157
fastjet::VoronoiAreaSpec * fVorAreaSpec
Definition: AliFJWrapper.h:113
Int_t mode
Definition: anaM.C:40
Double_t fDRStep
Definition: AliFJWrapper.h:162
std::vector< double > fGRDenominatorSub
Definition: AliFJWrapper.h:166
Bool_t GetLegacyMode()
Definition: AliFJWrapper.h:45
fastjet::GhostedAreaSpec * fGhostedAreaSpec
Definition: AliFJWrapper.h:114
Bool_t fUseArea4Vector
Definition: AliFJWrapper.h:140
void SetNRepeats(Int_t nrepeat)
Definition: AliFJWrapper.h:84
virtual std::vector< double > GetGRNumeratorSub() const
Definition: AliFJWrapper.h:61
virtual ~AliFJWrapper()
std::vector< double > fGRDenominator
Definition: AliFJWrapper.h:164
Double_t fGhostArea
Definition: AliFJWrapper.h:130
virtual void GetMedianAndSigma(Double_t &median, Double_t &sigma, Int_t remove=0) const
TString fName
Definition: AliFJWrapper.h:104
void SetGhostArea(Double_t gharea)
Definition: AliFJWrapper.h:85
Double_t fMedUsedForBgSub
Definition: AliFJWrapper.h:139
std::vector< fastjet::PseudoJet > GetJetConstituents(UInt_t idx) const
fastjet::JetDefinition * fJetDef
Definition: AliFJWrapper.h:115
virtual Int_t DoGenericSubtractionJetSigma2()
void SetupStrategyfromOpt(const char *option)
const std::vector< fastjet::PseudoJet > & GetInputVectors() const
Definition: AliFJWrapper.h:30
fastjet::AreaDefinition * fAreaDef
Definition: AliFJWrapper.h:112
virtual std::vector< double > GetGRDenominatorSub() const
Definition: AliFJWrapper.h:62
Bool_t fUseExternalBkg
Definition: AliFJWrapper.h:158
void SetTitle(const char *title)
Definition: AliFJWrapper.h:79
void SetPluginAlgor(Int_t plugin)
Definition: AliFJWrapper.h:91
virtual Int_t DoConstituentSubtraction()
Double_t fRhom
Definition: AliFJWrapper.h:160
virtual std::vector< double > GetGRDenominator() const
Definition: AliFJWrapper.h:60
void SetRhoRhom(Double_t rho, Double_t rhom)
Definition: AliFJWrapper.h:101
void SetGridScatter(Double_t gridSc)
Definition: AliFJWrapper.h:88
std::vector< fastjet::PseudoJet > fConstituentSubtrJets
Definition: AliFJWrapper.h:111
void SetMeanGhostKt(Double_t meankt)
Definition: AliFJWrapper.h:90
virtual void AddInputVector(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index=-99999)
virtual void AddInputGhost(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index=-99999)
void SetLegacyFJ()
std::vector< fastjet::PseudoJet > GetFilteredJetConstituents(UInt_t idx) const
void SetKtScatter(Double_t ktSc)
Definition: AliFJWrapper.h:89
virtual std::vector< double > GetGRNumerator() const
Definition: AliFJWrapper.h:59
virtual const char * ClassName() const
Definition: AliFJWrapper.h:22
virtual void RemoveLastInputVector()
void SetAreaType(const fastjet::AreaType &atype)
Definition: AliFJWrapper.h:83
Double_t fRho
Definition: AliFJWrapper.h:159
Bool_t GetDoFilterArea()
Definition: AliFJWrapper.h:46
virtual Int_t DoGenericSubtractionJetLeSub()
virtual Int_t DoGenericSubtractionJetConstituent()
const char * GetName() const
Definition: AliFJWrapper.h:37
const std::vector< fastjet::PseudoJet > & GetFilteredJets() const
Definition: AliFJWrapper.h:33
std::vector< double > fGRNumerator
Definition: AliFJWrapper.h:163
Double_t fMeanGhostKt
Definition: AliFJWrapper.h:136
fastjet::Strategy fStrategy
Definition: AliFJWrapper.h:125
TString fTitle
Definition: AliFJWrapper.h:105