AliPhysics  a9863a5 (a9863a5)
 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 
12 
14 {
15  public:
16  AliFJWrapper(const char *name, const char *title);
17  virtual ~AliFJWrapper();
18 
19  virtual void AddInputVector (Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index = -99999);
20  virtual void AddInputVector (const fastjet::PseudoJet& vec, Int_t index = -99999);
21  virtual void AddInputVectors(const std::vector<fastjet::PseudoJet>& vecs, Int_t offsetIndex = -99999);
22  virtual void AddInputGhost (Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index = -99999);
23  virtual const char *ClassName() const { return "AliFJWrapper"; }
24  virtual void Clear(const Option_t* /*opt*/ = "");
25  virtual void ClearMemory();
26  virtual void CopySettingsFrom (const AliFJWrapper& wrapper);
27  virtual void GetMedianAndSigma(Double_t& median, Double_t& sigma, Int_t remove = 0) const;
28  fastjet::ClusterSequenceArea* GetClusterSequence() const { return fClustSeq; }
29  fastjet::ClusterSequence* GetClusterSequenceSA() const { return fClustSeqSA; }
30  fastjet::ClusterSequenceActiveAreaExplicitGhosts* GetClusterSequenceGhosts() const { return fClustSeqActGhosts; }
31  const std::vector<fastjet::PseudoJet>& GetInputVectors() const { return fInputVectors; }
32  const std::vector<fastjet::PseudoJet>& GetInputGhosts() const { return fInputGhosts; }
33  const std::vector<fastjet::PseudoJet>& GetInclusiveJets() const { return fInclusiveJets; }
34  const std::vector<fastjet::PseudoJet>& GetFilteredJets() const { return fFilteredJets; }
35  std::vector<fastjet::PseudoJet> GetJetConstituents(UInt_t idx) const;
36  std::vector<fastjet::PseudoJet> GetFilteredJetConstituents(UInt_t idx) const;
37  Double_t GetMedianUsedForBgSubtraction() const { return fMedUsedForBgSub; }
38  const char* GetName() const { return fName; }
39  const char* GetTitle() const { return fTitle; }
40  Double_t GetJetArea (UInt_t idx) const;
41  fastjet::PseudoJet GetJetAreaVector (UInt_t idx) const;
42  Double_t GetFilteredJetArea (UInt_t idx) const;
43  fastjet::PseudoJet GetFilteredJetAreaVector(UInt_t idx) const;
44  Double_t GetJetSubtractedPt (UInt_t idx) const;
45  virtual std::vector<double> GetSubtractedJetsPts(Double_t median_pt = -1, Bool_t sorted = kFALSE);
46  Bool_t GetLegacyMode() { return fLegacyMode; }
47  Bool_t GetDoFilterArea() { return fDoFilterArea; }
48  Double_t NSubjettiness(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Int_t Option=0);
49  Double32_t NSubjettinessDerivativeSub(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Double_t JetR, fastjet::PseudoJet jet, Int_t Option=0);
50 #ifdef FASTJET_VERSION
51  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetMass() const {return fGenSubtractorInfoJetMass ; }
52  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetAngularity() const {return fGenSubtractorInfoJetAngularity ; }
53  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetpTD() const {return fGenSubtractorInfoJetpTD ; }
54  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetCircularity() const {return fGenSubtractorInfoJetCircularity ; }
55  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetSigma2() const {return fGenSubtractorInfoJetSigma2 ; }
56  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetConstituent() const {return fGenSubtractorInfoJetConstituent ; }
57  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetLeSub() const {return fGenSubtractorInfoJetLeSub ; }
58  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJet1subjettiness_kt() const {return fGenSubtractorInfoJet1subjettiness_kt ; }
59  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJet2subjettiness_kt() const {return fGenSubtractorInfoJet2subjettiness_kt ; }
60  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJet3subjettiness_kt() const {return fGenSubtractorInfoJet3subjettiness_kt ; }
61  const std::vector<fastjet::contrib::GenericSubtractorInfo> GetGenSubtractorInfoJetOpeningAngle_kt() const {return fGenSubtractorInfoJetOpeningAngle_kt ; }
62  const std::vector<fastjet::PseudoJet> GetConstituentSubtrJets() const {return fConstituentSubtrJets ; }
63  Int_t CreateGenSub(); // fastjet::contrib::GenericSubtractor
64  Int_t CreateConstituentSub(); // fastjet::contrib::ConstituentSubtractor
65 #endif
66  virtual std::vector<double> GetGRNumerator() const { return fGRNumerator ; }
67  virtual std::vector<double> GetGRDenominator() const { return fGRDenominator ; }
68  virtual std::vector<double> GetGRNumeratorSub() const { return fGRNumeratorSub ; }
69  virtual std::vector<double> GetGRDenominatorSub() const { return fGRDenominatorSub ; }
70 
71  virtual void RemoveLastInputVector();
72 
73  virtual Int_t Run();
74  virtual Int_t Filter();
75  virtual Int_t DoGenericSubtractionJetMass();
76  virtual Int_t DoGenericSubtractionGR(Int_t ijet);
77  virtual Int_t DoGenericSubtractionJetAngularity();
78  virtual Int_t DoGenericSubtractionJetpTD();
79  virtual Int_t DoGenericSubtractionJetCircularity();
80  virtual Int_t DoGenericSubtractionJetSigma2();
81  virtual Int_t DoGenericSubtractionJetConstituent();
82  virtual Int_t DoGenericSubtractionJetLeSub();
87  virtual Int_t DoConstituentSubtraction();
88 
89  void SetName(const char* name) { fName = name; }
90  void SetTitle(const char* title) { fTitle = title; }
91  void SetStrategy(const fastjet::Strategy &strat) { fStrategy = strat; }
92  void SetAlgorithm(const fastjet::JetAlgorithm &algor) { fAlgor = algor; }
93  void SetRecombScheme(const fastjet::RecombinationScheme &scheme) { fScheme = scheme; }
94  void SetAreaType(const fastjet::AreaType &atype) { fAreaType = atype; }
95  void SetNRepeats(Int_t nrepeat) { fNGhostRepeats = nrepeat; }
96  void SetGhostArea(Double_t gharea) { fGhostArea = gharea; }
97  void SetMaxRap(Double_t maxrap) { fMaxRap = maxrap; }
98  void SetR(Double_t r) { fR = r; }
99  void SetGridScatter(Double_t gridSc) { fGridScatter = gridSc; }
100  void SetKtScatter(Double_t ktSc) { fKtScatter = ktSc; }
101  void SetMeanGhostKt(Double_t meankt) { fMeanGhostKt = meankt; }
102  void SetPluginAlgor(Int_t plugin) { fPluginAlgor = plugin; }
103  void SetUseArea4Vector(Bool_t useA4v) { fUseArea4Vector = useA4v; }
104  void SetupAlgorithmfromOpt(const char *option);
105  void SetupAreaTypefromOpt(const char *option);
106  void SetupSchemefromOpt(const char *option);
107  void SetupStrategyfromOpt(const char *option);
108  void SetLegacyMode (Bool_t mode) { fLegacyMode ^= mode; }
109  void SetLegacyFJ();
110  void SetUseExternalBkg(Bool_t b, Double_t rho, Double_t rhom) { fUseExternalBkg = b; fRho = rho; fRhom = rhom;}
111  void SetRMaxAndStep(Double_t rmax, Double_t dr) {fRMax = rmax; fDRStep = dr; }
112  void SetRhoRhom (Double_t rho, Double_t rhom) { fUseExternalBkg = kTRUE; fRho = rho; fRhom = rhom;} // if using rho,rhom then fUseExternalBkg is true
113  void SetMinJetPt(Double_t MinPt) {fMinJetPt=MinPt;}
114 
115  protected:
116  TString fName;
117  TString fTitle;
118  std::vector<fastjet::PseudoJet> fInputVectors;
119  std::vector<fastjet::PseudoJet> fInputGhosts;
120  std::vector<fastjet::PseudoJet> fInclusiveJets;
121  std::vector<fastjet::PseudoJet> fFilteredJets;
122  std::vector<double> fSubtractedJetsPt;
123  std::vector<fastjet::PseudoJet> fConstituentSubtrJets;
124  fastjet::AreaDefinition *fAreaDef;
125  fastjet::VoronoiAreaSpec *fVorAreaSpec;
126  fastjet::GhostedAreaSpec *fGhostedAreaSpec;
127  fastjet::JetDefinition *fJetDef;
128  fastjet::JetDefinition::Plugin *fPlugin;
129 #ifndef FASTJET_VERSION
130  fastjet::RangeDefinition *fRange;
131 #else
132  fastjet::Selector *fRange;
133 #endif
134  fastjet::ClusterSequenceArea *fClustSeq;
135  fastjet::ClusterSequence *fClustSeqSA;
136  fastjet::ClusterSequenceActiveAreaExplicitGhosts *fClustSeqActGhosts;
137  fastjet::Strategy fStrategy;
138  fastjet::JetAlgorithm fAlgor;
139  fastjet::RecombinationScheme fScheme;
140  fastjet::AreaType fAreaType;
142  Double_t fGhostArea;
143  Double_t fMaxRap;
144  Double_t fR;
145  Double_t fMinJetPt;
146  // no setters for the moment - used default values in the constructor
147  Double_t fGridScatter;
148  Double_t fKtScatter;
149  Double_t fMeanGhostKt;
150  Int_t fPluginAlgor;
151  // extra parameters
152  Double_t fMedUsedForBgSub;
154 #ifdef FASTJET_VERSION
155  fastjet::JetMedianBackgroundEstimator *fBkrdEstimator;
156  //from contrib package
157  fastjet::contrib::GenericSubtractor *fGenSubtractor;
158  fastjet::contrib::ConstituentSubtractor *fConstituentSubtractor;
159  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetMass;
160  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoGRNum;
161  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoGRDen;
162  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetAngularity;
163  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetpTD;
164  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetCircularity;
165  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetSigma2;
166  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetConstituent;
167  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetLeSub;
168  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJet1subjettiness_kt;
169  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJet2subjettiness_kt;
170  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJet3subjettiness_kt;
171  std::vector<fastjet::contrib::GenericSubtractorInfo> fGenSubtractorInfoJetOpeningAngle_kt;
172 #endif
173  Bool_t fDoFilterArea;
174  Bool_t fLegacyMode;
176  Double_t fRho; // pT background density
177  Double_t fRhom; // mT background density
178  Double_t fRMax;
179  Double_t fDRStep;
180  std::vector<double> fGRNumerator;
181  std::vector<double> fGRDenominator;
182  std::vector<double> fGRNumeratorSub;
183  std::vector<double> fGRDenominatorSub;
184 
185  virtual void SubtractBackground(const Double_t median_pt = -1);
186 
187  private:
188  AliFJWrapper();
189  AliFJWrapper(const AliFJWrapper& wrapper);
190  AliFJWrapper& operator = (const AliFJWrapper& wrapper);
191 };
192 #endif
193 #endif
194 
195 #ifdef AliFJWrapper_CXX
196 #undef AliFJWrapper_CXX
197 
198 #if defined __GNUC__
199 #pragma GCC system_header
200 #endif
201 
202 namespace fj = fastjet;
203 
204 //_________________________________________________________________________________________________
205 AliFJWrapper::AliFJWrapper(const char *name, const char *title)
206  :
207  fName (name)
208  , fTitle (title)
209  , fInputVectors ( )
210  , fInputGhosts ( )
211  , fInclusiveJets ( )
212  , fFilteredJets ( )
213  , fSubtractedJetsPt ( )
214  , fConstituentSubtrJets ( )
215  , fAreaDef (0)
216  , fVorAreaSpec (0)
217  , fGhostedAreaSpec (0)
218  , fJetDef (0)
219  , fPlugin (0)
220  , fRange (0)
221  , fClustSeq (0)
222  , fClustSeqSA (0)
223  , fClustSeqActGhosts (0)
224  , fStrategy (fj::Best)
225  , fAlgor (fj::kt_algorithm)
226  , fScheme (fj::BIpt_scheme)
227  , fAreaType (fj::active_area)
228  , fNGhostRepeats (1)
229  , fGhostArea (0.005)
230  , fMaxRap (1.)
231  , fR (0.4)
232  , fGridScatter (1.0)
233  , fKtScatter (0.1)
234  , fMeanGhostKt (1e-100)
235  , fPluginAlgor (0)
236  , fMedUsedForBgSub (0)
237  , fUseArea4Vector (kFALSE)
238 #ifdef FASTJET_VERSION
239  , fBkrdEstimator (0)
240  , fGenSubtractor (0)
241  , fConstituentSubtractor (0)
242  , fGenSubtractorInfoJetMass ( )
243  , fGenSubtractorInfoGRNum ( )
244  , fGenSubtractorInfoGRDen ( )
245  , fGenSubtractorInfoJetAngularity ( )
246  , fGenSubtractorInfoJetpTD ( )
247  , fGenSubtractorInfoJetCircularity( )
248  , fGenSubtractorInfoJetSigma2()
249  , fGenSubtractorInfoJetConstituent ( )
250  , fGenSubtractorInfoJetLeSub ( )
251  , fGenSubtractorInfoJet1subjettiness_kt ( )
252  , fGenSubtractorInfoJet2subjettiness_kt ( )
253  , fGenSubtractorInfoJet3subjettiness_kt ( )
254  , fGenSubtractorInfoJetOpeningAngle_kt ( )
255 #endif
256  , fDoFilterArea (false)
257  , fLegacyMode (false)
258  , fUseExternalBkg (false)
259  , fRho (0)
260  , fRhom (0)
261  , fRMax(2.)
262  , fDRStep(0.04)
263  , fGRNumerator()
264  , fGRDenominator()
265  , fGRNumeratorSub()
266  , fGRDenominatorSub()
267 {
268  // Constructor.
269 }
270 
271 //_________________________________________________________________________________________________
273 {
274  // Destructor.
275  ClearMemory();
276 }
277 
278 //_________________________________________________________________________________________________
280 {
281  // Destructor.
282  if (fAreaDef) { delete fAreaDef; fAreaDef = NULL; }
283  if (fVorAreaSpec) { delete fVorAreaSpec; fVorAreaSpec = NULL; }
284  if (fGhostedAreaSpec) { delete fGhostedAreaSpec; fGhostedAreaSpec = NULL; }
285  if (fJetDef) { delete fJetDef; fJetDef = NULL; }
286  if (fPlugin) { delete fPlugin; fPlugin = NULL; }
287  if (fRange) { delete fRange; fRange = NULL; }
288  if (fClustSeq) { delete fClustSeq; fClustSeq = NULL; }
289  if (fClustSeqSA) { delete fClustSeqSA; fClustSeqSA = NULL; }
291  #ifdef FASTJET_VERSION
292  if (fBkrdEstimator) { delete fBkrdEstimator; fBkrdEstimator = NULL; }
293  if (fGenSubtractor) { delete fGenSubtractor; fGenSubtractor = NULL; }
294  if (fConstituentSubtractor) { delete fConstituentSubtractor; fConstituentSubtractor = NULL; }
295  #endif
296 }
297 
298 //_________________________________________________________________________________________________
299 void AliFJWrapper::CopySettingsFrom(const AliFJWrapper& wrapper)
300 {
301  // Copy some settings.
302  // You very often want to keep most of the settings
303  // but change only the algorithm or R - do it after call to this function
304 
305  fStrategy = wrapper.fStrategy;
306  fAlgor = wrapper.fAlgor;
307  fScheme = wrapper.fScheme;
308  fAreaType = wrapper.fAreaType;
309  fNGhostRepeats = wrapper.fNGhostRepeats;
310  fGhostArea = wrapper.fGhostArea;
311  fMaxRap = wrapper.fMaxRap;
312  fR = wrapper.fR;
313  fGridScatter = wrapper.fGridScatter;
314  fKtScatter = wrapper.fKtScatter;
315  fMeanGhostKt = wrapper.fMeanGhostKt;
316  fPluginAlgor = wrapper.fPluginAlgor;
318  fLegacyMode = wrapper.fLegacyMode;
320  fRho = wrapper.fRho;
321  fRhom = wrapper.fRhom;
322 }
323 
324 //_________________________________________________________________________________________________
325 void AliFJWrapper::Clear(const Option_t */*opt*/)
326 {
327  // Simply clear the input vectors.
328  // Make sure done on every event if the instance is reused
329  // Reset the median to zero.
330 
331  fInputVectors.clear();
332  fInputGhosts.clear();
333  fMedUsedForBgSub = 0;
334 
335  // for the moment brute force delete everything
336  ClearMemory();
337 }
338 
339 //_________________________________________________________________________________________________
341 {
342  // Remove last input vector
343 
344  fInputVectors.pop_back();
345 }
346 
347 //_________________________________________________________________________________________________
348 void AliFJWrapper::AddInputVector(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index)
349 {
350  // Make the input pseudojet.
351 
352  fastjet::PseudoJet inVec(px, py, pz, E);
353 
354  // Salvatore Aiola: not sure why this was done...
355  //if (index > -99999) {
356  inVec.set_user_index(index);
357  //} else {
358  //inVec.set_user_index(fInputVectors.size());
359  //}
360 
361  // add to the fj container of input vectors
362  fInputVectors.push_back(inVec);
363 }
364 
365 //_________________________________________________________________________________________________
366 void AliFJWrapper::AddInputVector(const fj::PseudoJet& vec, Int_t index)
367 {
368  // Add an input pseudojet.
369 
370  fj::PseudoJet inVec = vec;
371 
372  // Salvatore Aiola: not sure why this was done...
374  inVec.set_user_index(index);
375  //} else {
376  //inVec.set_user_index(fInputVectors.size());
377  //}
378 
379  // add to the fj container of input vectors
380  fInputVectors.push_back(inVec);
381 }
382 
383 //_________________________________________________________________________________________________
384 void AliFJWrapper::AddInputVectors(const std::vector<fj::PseudoJet>& vecs, Int_t offsetIndex)
385 {
386  // Add the input from vector of pseudojets.
387 
388  for (UInt_t i = 0; i < vecs.size(); ++i) {
389  fj::PseudoJet inVec = vecs[i];
390  if (offsetIndex > -99999)
391  inVec.set_user_index(fInputVectors.size() + offsetIndex);
392  // add to the fj container of input vectors
393  fInputVectors.push_back(inVec);
394  }
395 }
396 
397 //_________________________________________________________________________________________________
398 void AliFJWrapper::AddInputGhost(Double_t px, Double_t py, Double_t pz, Double_t E, Int_t index)
399 {
400  // Make the input pseudojet.
401 
402  fastjet::PseudoJet inVec(px, py, pz, E);
403 
404  if (index > -99999) {
405  inVec.set_user_index(index);
406  } else {
407  inVec.set_user_index(fInputGhosts.size());
408  }
409 
410  // add to the fj container of input vectors
411  fInputGhosts.push_back(inVec);
412  if (!fDoFilterArea) fDoFilterArea = kTRUE;
413 }
414 
415 //_________________________________________________________________________________________________
416 Double_t AliFJWrapper::GetJetArea(UInt_t idx) const
417 {
418  // Get the jet area.
419 
420  Double_t retval = -1; // really wrong area..
421  if ( idx < fInclusiveJets.size() ) {
422  retval = fClustSeq->area(fInclusiveJets[idx]);
423  } else {
424  AliError(Form("[e] ::GetJetArea wrong index: %d",idx));
425  }
426  return retval;
427 }
428 
429 //_________________________________________________________________________________________________
430 Double_t AliFJWrapper::GetFilteredJetArea(UInt_t idx) const
431 {
432  // Get the filtered jet area.
433 
434  Double_t retval = -1; // really wrong area..
435  if (fDoFilterArea && fClustSeqActGhosts && (idx<fFilteredJets.size())) {
436  retval = fClustSeqActGhosts->area(fFilteredJets[idx]);
437  } else {
438  AliError(Form("[e] ::GetFilteredJetArea wrong index: %d",idx));
439  }
440  return retval;
441 }
442 
443 //_________________________________________________________________________________________________
444 fastjet::PseudoJet AliFJWrapper::GetJetAreaVector(UInt_t idx) const
445 {
446  // Get the jet area as vector.
447  fastjet::PseudoJet retval;
448  if ( idx < fInclusiveJets.size() ) {
449  retval = fClustSeq->area_4vector(fInclusiveJets[idx]);
450  } else {
451  AliError(Form("[e] ::GetJetArea wrong index: %d",idx));
452  }
453  return retval;
454 }
455 
456 //_________________________________________________________________________________________________
457 fastjet::PseudoJet AliFJWrapper::GetFilteredJetAreaVector(UInt_t idx) const
458 {
459  // Get the jet area as vector.
460  fastjet::PseudoJet retval;
461  if (fDoFilterArea && fClustSeqActGhosts && (idx<fFilteredJets.size())) {
462  retval = fClustSeqActGhosts->area_4vector(fFilteredJets[idx]);
463  } else {
464  AliError(Form("[e] ::GetFilteredJetArea wrong index: %d",idx));
465  }
466  return retval;
467 }
468 
469 //_________________________________________________________________________________________________
470 std::vector<double> AliFJWrapper::GetSubtractedJetsPts(Double_t median_pt, Bool_t sorted)
471 {
472  // Get subtracted jets pTs, returns vector.
473 
474  SubtractBackground(median_pt);
475 
476  if (kTRUE == sorted) {
477  std::sort(fSubtractedJetsPt.begin(), fSubtractedJetsPt.begin());
478  }
479  return fSubtractedJetsPt;
480 }
481 
482 //_________________________________________________________________________________________________
483 Double_t AliFJWrapper::GetJetSubtractedPt(UInt_t idx) const
484 {
485  // Get subtracted jets pTs, returns Double_t.
486 
487  Double_t retval = -99999.; // really wrong pt..
488  if ( idx < fSubtractedJetsPt.size() ) {
489  retval = fSubtractedJetsPt[idx];
490  }
491  return retval;
492 }
493 
494 //_________________________________________________________________________________________________
495 std::vector<fastjet::PseudoJet>
496 AliFJWrapper::GetJetConstituents(UInt_t idx) const
497 {
498  // Get jets constituents.
499 
500  std::vector<fastjet::PseudoJet> retval;
501 
502  if ( idx < fInclusiveJets.size() ) {
503  retval = fClustSeq->constituents(fInclusiveJets[idx]);
504  } else {
505  AliError(Form("[e] ::GetJetConstituents wrong index: %d",idx));
506  }
507 
508  return retval;
509 }
510 
511 //_________________________________________________________________________________________________
512 std::vector<fastjet::PseudoJet>
514 {
515  // Get jets constituents.
516 
517  std::vector<fastjet::PseudoJet> retval;
518 
519  if ( idx < fFilteredJets.size() ) {
520  if (fClustSeqSA) retval = fClustSeqSA->constituents(fFilteredJets[idx]);
521  if (fClustSeqActGhosts) retval = fClustSeqActGhosts->constituents(fFilteredJets[idx]);
522  } else {
523  AliError(Form("[e] ::GetFilteredJetConstituents wrong index: %d",idx));
524  }
525 
526  return retval;
527 }
528 
529 //_________________________________________________________________________________________________
530 void AliFJWrapper::GetMedianAndSigma(Double_t &median, Double_t &sigma, Int_t remove) const
531 {
532  // Get the median and sigma from fastjet.
533  // User can also do it on his own because the cluster sequence is exposed (via a getter)
534 
535  if (!fClustSeq) {
536  AliError("[e] Run the jfinder first.");
537  return;
538  }
539 
540  Double_t mean_area = 0;
541  try {
542  if(0 == remove) {
543  fClustSeq->get_median_rho_and_sigma(*fRange, fUseArea4Vector, median, sigma, mean_area);
544  } else {
545  std::vector<fastjet::PseudoJet> input_jets = sorted_by_pt(fClustSeq->inclusive_jets());
546  input_jets.erase(input_jets.begin(), input_jets.begin() + remove);
547  fClustSeq->get_median_rho_and_sigma(input_jets, *fRange, fUseArea4Vector, median, sigma, mean_area);
548  input_jets.clear();
549  }
550  } catch (fj::Error) {
551  AliError(" [w] FJ Exception caught.");
552  median = -1.;
553  sigma = -1;
554  }
555 }
556 
557 //_________________________________________________________________________________________________
558 Int_t AliFJWrapper::Run()
559 {
560  // Run the actual jet finder.
561 
562  if (fAreaType == fj::voronoi_area) {
563  // Rfact - check dependence - default is 1.
564  // NOTE: hardcoded variable!
565  fVorAreaSpec = new fj::VoronoiAreaSpec(1.);
566  fAreaDef = new fj::AreaDefinition(*fVorAreaSpec);
567  } else {
568  fGhostedAreaSpec = new fj::GhostedAreaSpec(fMaxRap,
570  fGhostArea,
571  fGridScatter,
572  fKtScatter,
573  fMeanGhostKt);
574 
575  fAreaDef = new fj::AreaDefinition(*fGhostedAreaSpec, fAreaType);
576  }
577 
578  // this is acceptable by fastjet:
579 #ifndef FASTJET_VERSION
580  fRange = new fj::RangeDefinition(fMaxRap - 0.95 * fR);
581 #else
582  fRange = new fj::Selector(fj::SelectorAbsRapMax(fMaxRap - 0.95 * fR));
583 #endif
584 
585  if (fAlgor == fj::plugin_algorithm) {
586  if (fPluginAlgor == 0) {
587  // SIS CONE ALGOR
588  // NOTE: hardcoded split parameter
589  Double_t overlap_threshold = 0.75; // NOTE: this actually splits a lot: thr/min(pt1,pt2)
590  fPlugin = new fj::SISConePlugin(fR,
591  overlap_threshold,
592  0, //search of stable cones - zero = until no more
593  1.0); // this should be seed effectively for proto jets
594  fJetDef = new fastjet::JetDefinition(fPlugin);
595  } else if (fPluginAlgor == 1) {
596  // CDF cone
597  // NOTE: hardcoded split parameter
598  Double_t overlap_threshold = 0.75; // NOTE: this actually splits a lot: thr/min(pt1,pt2)
599  fPlugin = new fj::CDFMidPointPlugin(fR,
600  overlap_threshold,
601  1.0, //search of stable cones - zero = until no more
602  1.0); // this should be seed effectively for proto jets
603  fJetDef = new fastjet::JetDefinition(fPlugin);
604  } else {
605  AliError("[e] Unrecognized plugin number!");
606  }
607  } else {
608  fJetDef = new fj::JetDefinition(fAlgor, fR, fScheme, fStrategy);
609  }
610 
611  try {
612  fClustSeq = new fj::ClusterSequenceArea(fInputVectors, *fJetDef, *fAreaDef);
613  } catch (fj::Error) {
614  AliError(" [w] FJ Exception caught.");
615  return -1;
616  }
617 
618  // FJ3 :: Define an JetMedianBackgroundEstimator just in case it will be used
619 #ifdef FASTJET_VERSION
620  fBkrdEstimator = new fj::JetMedianBackgroundEstimator(fj::SelectorAbsRapMax(fMaxRap));
621 #endif
622 
623  if (fLegacyMode) { SetLegacyFJ(); } // for FJ 2.x even if fLegacyMode is set, SetLegacyFJ is dummy
624 
625  // inclusive jets:
626  fInclusiveJets.clear();
627  fInclusiveJets = fClustSeq->inclusive_jets(0.0);
628 
629  return 0;
630 }
631 
632 //_________________________________________________________________________________________________
633 Int_t AliFJWrapper::Filter()
634 {
635 //
636 // AliFJWrapper::Filter
637 //
638 
639  fJetDef = new fj::JetDefinition(fAlgor, fR, fScheme, fStrategy);
640 
641  if (fDoFilterArea) {
642  if (fInputGhosts.size()>0) {
643  try {
644  fClustSeqActGhosts = new fj::ClusterSequenceActiveAreaExplicitGhosts(fInputVectors,
645  *fJetDef,
646  fInputGhosts,
647  fGhostArea);
648  } catch (fj::Error) {
649  AliError(" [w] FJ Exception caught.");
650  return -1;
651  }
652 
653  fFilteredJets.clear();
654  fFilteredJets = fClustSeqActGhosts->inclusive_jets(0.0);
655  } else {
656  return -1;
657  }
658  } else {
659  try {
660  fClustSeqSA = new fastjet::ClusterSequence(fInputVectors, *fJetDef);
661  } catch (fj::Error) {
662  AliError(" [w] FJ Exception caught.");
663  return -1;
664  }
665 
666  fFilteredJets.clear();
667  fFilteredJets = fClustSeqSA->inclusive_jets(0.0);
668  }
669 
670  return 0;
671 }
672 
673 //_________________________________________________________________________________________________
675 {
676  // This methods enable legacy behaviour (FastJet 2.x) when AliROOT is compiled with FastJet 3.x
677 #ifdef FASTJET_VERSION
678  std::cout << "WARNING! Setting FastJet in legacy mode" << std::endl;
679  if (fGhostedAreaSpec) { fGhostedAreaSpec->set_fj2_placement(kTRUE); }
680  if (fBkrdEstimator) {
681  fBkrdEstimator->set_provide_fj2_sigma(kTRUE);
682  fBkrdEstimator->set_use_area_4vector(kFALSE);
683  }
684 #endif
685 }
686 
687 //_________________________________________________________________________________________________
688 void AliFJWrapper::SubtractBackground(Double_t median_pt)
689 {
690  // Subtract the background (specify the value - see below the meaning).
691  // Negative argument means the bg will be determined with the current algorithm
692  // this is the default behavior. Zero is allowed
693  // Note: user may set the switch for area4vector based subtraction.
694 
695  Double_t median = 0;
696  Double_t sigma = 0;
697  Double_t mean_area = 0;
698 
699  // clear the subtracted jet pt's vector<double>
700  fSubtractedJetsPt.clear();
701 
702  // check what was specified (default is -1)
703  if (median_pt < 0) {
704  try {
705  fClustSeq->get_median_rho_and_sigma(*fRange, fUseArea4Vector, median, sigma, mean_area);
706  }
707 
708  catch (fj::Error) {
709  AliError(" [w] FJ Exception caught.");
710  median = -9999.;
711  sigma = -1;
712  fMedUsedForBgSub = median;
713  return;
714  }
715  fMedUsedForBgSub = median;
716  } else {
717  // we do not know the sigma in this case
718  sigma = -1;
719  if (0.0 == median_pt) {
720  // AliWarning(Form(" [w] Median specified for bg subtraction is ZERO: %f", median_pt )) << endl;
721  fMedUsedForBgSub = 0.;
722  } else {
723  fMedUsedForBgSub = median_pt;
724  }
725  }
726 
727  // subtract:
728  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
729  if ( fUseArea4Vector ) {
730  // subtract the background using the area4vector
731  fj::PseudoJet area4v = fClustSeq->area_4vector(fInclusiveJets[i]);
732  fj::PseudoJet jet_sub = fInclusiveJets[i] - area4v * fMedUsedForBgSub;
733  fSubtractedJetsPt.push_back(jet_sub.perp()); // here we put only the pt of the jet - note: this can be negative
734  } else {
735  // subtract the background using scalars
736  // fj::PseudoJet jet_sub = fInclusiveJets[i] - area * fMedUsedForBgSub_;
737  Double_t area = fClustSeq->area(fInclusiveJets[i]);
738  // standard subtraction
739  Double_t pt_sub = fInclusiveJets[i].perp() - fMedUsedForBgSub * area;
740  fSubtractedJetsPt.push_back(pt_sub); // here we put only the pt of the jet - note: this can be negative
741  }
742  }
743 }
744 
745 //_________________________________________________________________________________________________
747  //Do generic subtraction for jet mass
748 #ifdef FASTJET_VERSION
749  CreateGenSub();
750 
751  // Define jet shape
752  AliJetShapeMass shapeMass;
753 
754  // clear the generic subtractor info vector
755  fGenSubtractorInfoJetMass.clear();
756  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
757  fj::contrib::GenericSubtractorInfo info;
758  if(fInclusiveJets[i].perp()>1.e-4)
759  double subtracted_shape = (*fGenSubtractor)(shapeMass, fInclusiveJets[i], info);
760  fGenSubtractorInfoJetMass.push_back(info);
761  }
762 #endif
763  return 0;
764 }
765 
766 //_________________________________________________________________________________________________
767 Int_t AliFJWrapper::DoGenericSubtractionGR(Int_t ijet) {
768  //Do generic subtraction for jet mass
769 #ifdef FASTJET_VERSION
770  CreateGenSub();
771 
772  if(ijet>fInclusiveJets.size()) return 0;
773 
774  fGRNumerator.clear();
775  fGRDenominator.clear();
776  fGRNumeratorSub.clear();
777  fGRDenominatorSub.clear();
778 
779  // Define jet shape
780  for(Double_t r = 0.; r<fRMax; r+=fDRStep) {
781  AliJetShapeGRNum shapeGRNum(r,fDRStep);
782  AliJetShapeGRDen shapeGRDen(r,fDRStep);
783 
784  // clear the generic subtractor info vector
785  fGenSubtractorInfoGRNum.clear();
786  fGenSubtractorInfoGRDen.clear();
787  fj::contrib::GenericSubtractorInfo infoNum;
788  fj::contrib::GenericSubtractorInfo infoDen;
789  if(fInclusiveJets[ijet].perp()>1.e-4) {
790  double sub_num = (*fGenSubtractor)(shapeGRNum, fInclusiveJets[ijet], infoNum);
791  double sub_den = (*fGenSubtractor)(shapeGRDen, fInclusiveJets[ijet], infoDen);
792  }
793  fGenSubtractorInfoGRNum.push_back(infoNum);
794  fGenSubtractorInfoGRDen.push_back(infoDen);
795  fGRNumerator.push_back(infoNum.unsubtracted());
796  fGRDenominator.push_back(infoDen.unsubtracted());
797  fGRNumeratorSub.push_back(infoNum.second_order_subtracted());
798  fGRDenominatorSub.push_back(infoDen.second_order_subtracted());
799  }
800 #endif
801  return 0;
802 }
803 //_________________________________________________________________________________________________
805  //Do generic subtraction for jet mass
806 #ifdef FASTJET_VERSION
807  CreateGenSub();
808 
809  // Define jet shape
810  AliJetShapeAngularity shapeAngularity;
811 
812  // clear the generic subtractor info vector
813  fGenSubtractorInfoJetAngularity.clear();
814  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
815  fj::contrib::GenericSubtractorInfo infoAng;
816  if(fInclusiveJets[i].perp()>1.e-4)
817  double subtracted_shape = (*fGenSubtractor)(shapeAngularity, fInclusiveJets[i], infoAng);
818  fGenSubtractorInfoJetAngularity.push_back(infoAng);
819  }
820 #endif
821  return 0;
822 }
823 //_________________________________________________________________________________________________
825  //Do generic subtraction for jet mass
826 #ifdef FASTJET_VERSION
827  CreateGenSub();
828 
829  // Define jet shape
830  AliJetShapepTD shapepTD;
831 
832  // clear the generic subtractor info vector
833  fGenSubtractorInfoJetpTD.clear();
834  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
835  fj::contrib::GenericSubtractorInfo infopTD;
836  if(fInclusiveJets[i].perp()>1.e-4)
837  double subtracted_shape = (*fGenSubtractor)(shapepTD, fInclusiveJets[i], infopTD);
838  fGenSubtractorInfoJetpTD.push_back(infopTD);
839  }
840 #endif
841  return 0;
842 }
843 //_________________________________________________________________________________________________
845  //Do generic subtraction for jet mass
846 #ifdef FASTJET_VERSION
847  CreateGenSub();
848 
849  // Define jet shape
850  AliJetShapeCircularity shapecircularity;
851 
852  // clear the generic subtractor info vector
853  fGenSubtractorInfoJetCircularity.clear();
854  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
855  fj::contrib::GenericSubtractorInfo infoCirc;
856  if(fInclusiveJets[i].perp()>1.e-4)
857  double subtracted_shape = (*fGenSubtractor)(shapecircularity, fInclusiveJets[i], infoCirc);
858  fGenSubtractorInfoJetCircularity.push_back(infoCirc);
859  }
860 #endif
861  return 0;
862 }
863 //_________________________________________________________________________________________________
865  //Do generic subtraction for jet mass
866 #ifdef FASTJET_VERSION
867  CreateGenSub();
868 
869  // Define jet shape
870  AliJetShapeSigma2 shapesigma2;
871 
872  // clear the generic subtractor info vector
873  fGenSubtractorInfoJetSigma2.clear();
874  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
875  fj::contrib::GenericSubtractorInfo infoSigma;
876  if(fInclusiveJets[i].perp()>1.e-4)
877  double subtracted_shape = (*fGenSubtractor)(shapesigma2, fInclusiveJets[i], infoSigma);
878  fGenSubtractorInfoJetSigma2.push_back(infoSigma);
879  }
880 #endif
881  return 0;
882 }
883 //_________________________________________________________________________________________________
885  //Do generic subtraction for jet mass
886 #ifdef FASTJET_VERSION
887  CreateGenSub();
888 
889  // Define jet shape
890  AliJetShapeConstituent shapeconst;
891 
892  // clear the generic subtractor info vector
893  fGenSubtractorInfoJetConstituent.clear();
894  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
895  fj::contrib::GenericSubtractorInfo infoConst;
896  if(fInclusiveJets[i].perp()>1.e-4)
897  double subtracted_shape = (*fGenSubtractor)(shapeconst, fInclusiveJets[i], infoConst);
898  fGenSubtractorInfoJetConstituent.push_back(infoConst);
899  }
900 #endif
901  return 0;
902 }
903 
904 //_________________________________________________________________________________________________
906  //Do generic subtraction for jet mass
907 #ifdef FASTJET_VERSION
908  CreateGenSub();
909 
910  // Define jet shape
911  AliJetShapeLeSub shapeLeSub;
912 
913  // clear the generic subtractor info vector
914  fGenSubtractorInfoJetLeSub.clear();
915  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
916  fj::contrib::GenericSubtractorInfo infoLeSub;
917  if(fInclusiveJets[i].perp()>1.e-4)
918  double subtracted_shape = (*fGenSubtractor)(shapeLeSub, fInclusiveJets[i], infoLeSub);
919  fGenSubtractorInfoJetLeSub.push_back(infoLeSub);
920  }
921 #endif
922  return 0;
923 }
924 
925 //_________________________________________________________________________________________________
927  //Do generic subtraction for 1subjettiness
928 #ifdef FASTJET_VERSION
929  CreateGenSub();
930 
931  // Define jet shape
932  AliJetShape1subjettiness_kt shape1subjettiness_kt;
933 
934  // clear the generic subtractor info vector
935  fGenSubtractorInfoJet1subjettiness_kt.clear();
936  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
937  fj::contrib::GenericSubtractorInfo info1subjettiness_kt;
938  if(fInclusiveJets[i].perp()>1.e-4)
939  double subtracted_shape = (*fGenSubtractor)(shape1subjettiness_kt, fInclusiveJets[i], info1subjettiness_kt);
940  fGenSubtractorInfoJet1subjettiness_kt.push_back(info1subjettiness_kt);
941  }
942 #endif
943  return 0;
944 }
945 
946 //_________________________________________________________________________________________________
948  //Do generic subtraction for 2subjettiness
949 #ifdef FASTJET_VERSION
950  CreateGenSub();
951 
952  // Define jet shape
953  AliJetShape2subjettiness_kt shape2subjettiness_kt;
954 
955  // clear the generic subtractor info vector
956  fGenSubtractorInfoJet2subjettiness_kt.clear();
957  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
958  fj::contrib::GenericSubtractorInfo info2subjettiness_kt;
959  if(fInclusiveJets[i].perp()>1.e-4)
960  double subtracted_shape = (*fGenSubtractor)(shape2subjettiness_kt, fInclusiveJets[i], info2subjettiness_kt);
961  fGenSubtractorInfoJet2subjettiness_kt.push_back(info2subjettiness_kt);
962  }
963 #endif
964  return 0;
965 }
966 
967 //_________________________________________________________________________________________________
969  //Do generic subtraction for 3subjettiness
970 #ifdef FASTJET_VERSION
971  CreateGenSub();
972 
973  // Define jet shape
974  AliJetShape3subjettiness_kt shape3subjettiness_kt;
975 
976  // clear the generic subtractor info vector
977  fGenSubtractorInfoJet3subjettiness_kt.clear();
978  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
979  fj::contrib::GenericSubtractorInfo info3subjettiness_kt;
980  if(fInclusiveJets[i].perp()>1.e-4)
981  double subtracted_shape = (*fGenSubtractor)(shape3subjettiness_kt, fInclusiveJets[i], info3subjettiness_kt);
982  fGenSubtractorInfoJet3subjettiness_kt.push_back(info3subjettiness_kt);
983  }
984 #endif
985  return 0;
986 }
987 
988 //_________________________________________________________________________________________________
990  //Do generic subtraction for 2subjettiness axes opening angle
991 #ifdef FASTJET_VERSION
992  CreateGenSub();
993 
994  // Define jet shape
995  AliJetShapeOpeningAngle_kt shapeOpeningAngle_kt;
996 
997  // clear the generic subtractor info vector
998  fGenSubtractorInfoJetOpeningAngle_kt.clear();
999  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
1000  fj::contrib::GenericSubtractorInfo infoOpeningAngle_kt;
1001  if(fInclusiveJets[i].perp()>1.e-4)
1002  double subtracted_shape = (*fGenSubtractor)(shapeOpeningAngle_kt, fInclusiveJets[i], infoOpeningAngle_kt);
1003  fGenSubtractorInfoJetOpeningAngle_kt.push_back(infoOpeningAngle_kt);
1004  }
1005 #endif
1006  return 0;
1007 }
1008 
1009 //_________________________________________________________________________________________________
1011  //Do constituent subtraction
1012 #ifdef FASTJET_VERSION
1013  CreateConstituentSub();
1014  // fConstituentSubtractor->set_alpha(/* double alpha */);
1015  // fConstituentSubtractor->set_max_deltaR(/* double max_deltaR */);
1016 
1017  //clear constituent subtracted jets
1018  fConstituentSubtrJets.clear();
1019  for (unsigned i = 0; i < fInclusiveJets.size(); i++) {
1020  fj::PseudoJet subtracted_jet(0.,0.,0.,0.);
1021  if(fInclusiveJets[i].perp()>0.)
1022  subtracted_jet = (*fConstituentSubtractor)(fInclusiveJets[i]);
1023  fConstituentSubtrJets.push_back(subtracted_jet);
1024  }
1025  if(fConstituentSubtractor) { delete fConstituentSubtractor; fConstituentSubtractor = NULL; }
1026 
1027 #endif
1028  return 0;
1029 }
1030 
1031 //_________________________________________________________________________________________________
1032 Int_t AliFJWrapper::CreateGenSub() {
1033  //Do generic subtraction for jet mass
1034  #ifdef FASTJET_VERSION
1035  if (fGenSubtractor) { delete fGenSubtractor; } // protect against memory leaks
1036 
1037  if (fUseExternalBkg)
1038  { fGenSubtractor = new fj::contrib::GenericSubtractor(fRho,fRhom); }
1039  else
1040  {
1041  fGenSubtractor = new fj::contrib::GenericSubtractor(fBkrdEstimator);
1042  #if FASTJET_VERSION_NUMBER >= 30100
1043  fGenSubtractor->set_common_bge_for_rho_and_rhom(); // see contrib 1.020 GenericSubtractor.hh line 62
1044  #endif
1045  }
1046 
1047  #endif
1048  return 0;
1049 }
1050 
1051 //_________________________________________________________________________________________________
1052 Int_t AliFJWrapper::CreateConstituentSub() {
1053  //Do generic subtraction for jet mass
1054  #ifdef FASTJET_VERSION
1055  if (fConstituentSubtractor) { delete fConstituentSubtractor; } // protect against memory leaks
1056 
1057  // see ConstituentSubtractor.hh signatures
1058  // ConstituentSubtractor(double rho, double rhom=0, double alpha=0, double maxDeltaR=-1)
1059  if (fUseExternalBkg) { fConstituentSubtractor = new fj::contrib::ConstituentSubtractor(fRho,fRhom); }
1060  else { fConstituentSubtractor = new fj::contrib::ConstituentSubtractor(fBkrdEstimator); }
1061 
1062  #endif
1063  return 0;
1064 }
1065 
1066 //_________________________________________________________________________________________________
1067 void AliFJWrapper::SetupAlgorithmfromOpt(const char *option)
1068 {
1069  // Setup algorithm from char.
1070 
1071  std::string opt(option);
1072 
1073  if (!opt.compare("kt")) fAlgor = fj::kt_algorithm;
1074  if (!opt.compare("antikt")) fAlgor = fj::antikt_algorithm;
1075  if (!opt.compare("cambridge")) fAlgor = fj::cambridge_algorithm;
1076  if (!opt.compare("genkt")) fAlgor = fj::genkt_algorithm;
1077  if (!opt.compare("cambridge_passive")) fAlgor = fj::cambridge_for_passive_algorithm;
1078  if (!opt.compare("genkt_passive")) fAlgor = fj::genkt_for_passive_algorithm;
1079  if (!opt.compare("ee_kt")) fAlgor = fj::ee_kt_algorithm;
1080  if (!opt.compare("ee_genkt")) fAlgor = fj::ee_genkt_algorithm;
1081  if (!opt.compare("plugin")) fAlgor = fj::plugin_algorithm;
1082 }
1083 
1084 //_________________________________________________________________________________________________
1085 void AliFJWrapper::SetupAreaTypefromOpt(const char *option)
1086 {
1087  // Setup area type from char.
1088 
1089  std::string opt(option);
1090 
1091  if (!opt.compare("active")) fAreaType = fj::active_area;
1092  if (!opt.compare("invalid")) fAreaType = fj::invalid_area;
1093  if (!opt.compare("active_area_explicit_ghosts")) fAreaType = fj::active_area_explicit_ghosts;
1094  if (!opt.compare("one_ghost_passive")) fAreaType = fj::one_ghost_passive_area;
1095  if (!opt.compare("passive")) fAreaType = fj::passive_area;
1096  if (!opt.compare("voronoi")) fAreaType = fj::voronoi_area;
1097 }
1098 
1099 //_________________________________________________________________________________________________
1100 void AliFJWrapper::SetupSchemefromOpt(const char *option)
1101 {
1102  //
1103  // setup scheme from char
1104  //
1105 
1106  std::string opt(option);
1107 
1108  if (!opt.compare("BIpt")) fScheme = fj::BIpt_scheme;
1109  if (!opt.compare("BIpt2")) fScheme = fj::BIpt2_scheme;
1110  if (!opt.compare("E")) fScheme = fj::E_scheme;
1111  if (!opt.compare("pt")) fScheme = fj::pt_scheme;
1112  if (!opt.compare("pt2")) fScheme = fj::pt2_scheme;
1113  if (!opt.compare("Et")) fScheme = fj::Et_scheme;
1114  if (!opt.compare("Et2")) fScheme = fj::Et2_scheme;
1115 }
1116 
1117 //_________________________________________________________________________________________________
1118 void AliFJWrapper::SetupStrategyfromOpt(const char *option)
1119 {
1120  // Setup strategy from char.
1121 
1122  std::string opt(option);
1123 
1124  if (!opt.compare("Best")) fStrategy = fj::Best;
1125  if (!opt.compare("N2MinHeapTiled")) fStrategy = fj::N2MinHeapTiled;
1126  if (!opt.compare("N2Tiled")) fStrategy = fj::N2Tiled;
1127  if (!opt.compare("N2PoorTiled")) fStrategy = fj::N2PoorTiled;
1128  if (!opt.compare("N2Plain")) fStrategy = fj::N2Plain;
1129  if (!opt.compare("N3Dumb")) fStrategy = fj::N3Dumb;
1130  if (!opt.compare("NlnN")) fStrategy = fj::NlnN;
1131  if (!opt.compare("NlnN3pi")) fStrategy = fj::NlnN3pi;
1132  if (!opt.compare("NlnN4pi")) fStrategy = fj::NlnN4pi;
1133  if (!opt.compare("NlnNCam4pi")) fStrategy = fj::NlnNCam4pi;
1134  if (!opt.compare("NlnNCam2pi2R")) fStrategy = fj::NlnNCam2pi2R;
1135  if (!opt.compare("NlnNCam")) fStrategy = fj::NlnNCam;
1136  if (!opt.compare("plugin")) fStrategy = fj::plugin_strategy;
1137 }
1138 
1139 Double_t AliFJWrapper::NSubjettiness(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Int_t Option){
1140 
1141 
1142  //Option 0=Nsubjettiness result, 1=opening angle between axes in Eta-Phi plane, 2=Distance between axes in Eta-Phi plane
1143 
1144  fJetDef = new fj::JetDefinition(fAlgor, fR*100, fScheme, fStrategy ); //the *2 is becasue of a handful of jets that end up missing a track for some reason.
1145 
1146  try {
1147  fClustSeqSA = new fastjet::ClusterSequence(fInputVectors, *fJetDef);
1148  // ClustSeqSA = new fastjet::ClusterSequenceArea(fInputVectors, *fJetDef, *fAreaDef);
1149  } catch (fj::Error) {
1150  AliError(" [w] FJ Exception caught.");
1151  return -1;
1152  }
1153  fFilteredJets.clear();
1154  fFilteredJets = fClustSeqSA->inclusive_jets(fMinJetPt-0.1); //becasue this is < not <=
1155  Double_t Result=-1;
1156  std::vector<fastjet::PseudoJet> SubJet_Axes;
1157  fj::PseudoJet SubJet1_Axis;
1158  fj::PseudoJet SubJet2_Axis;
1159  if (Algorithm==0){
1160  fj::contrib::Nsubjettiness nSub(N, fj::contrib::KT_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1161  Result= nSub.result(fFilteredJets[0]);
1162  SubJet_Axes=nSub.currentAxes();
1163  }
1164  else if (Algorithm==1) {
1165  fj::contrib::Nsubjettiness nSub(N, fj::contrib::CA_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1166  Result= nSub.result(fFilteredJets[0]);
1167  SubJet_Axes=nSub.currentAxes();
1168  }
1169  else if (Algorithm==2){
1170  fj::contrib::Nsubjettiness nSub(N, fj::contrib::AntiKT_Axes(Radius), fj::contrib::NormalizedMeasure(Beta,fR));
1171  Result= nSub.result(fFilteredJets[0]);
1172  SubJet_Axes=nSub.currentAxes();
1173  }
1174  else if (Algorithm==3) {
1175  fj::contrib::Nsubjettiness nSub(N, fj::contrib::WTA_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1176  Result= nSub.result(fFilteredJets[0]);
1177  SubJet_Axes=nSub.currentAxes();
1178  }
1179  else if (Algorithm==4) {
1180  fj::contrib::Nsubjettiness nSub(N, fj::contrib::WTA_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1181  Result= nSub.result(fFilteredJets[0]);
1182  SubJet_Axes=nSub.currentAxes();
1183  }
1184  else if (Algorithm==5){
1185  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1186  Result= nSub.result(fFilteredJets[0]);
1187  SubJet_Axes=nSub.currentAxes();
1188  }
1189  else if (Algorithm==6){
1190  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1191  Result= nSub.result(fFilteredJets[0]);
1192  SubJet_Axes=nSub.currentAxes();
1193  }
1194  else if (Algorithm==7){
1195  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_AntiKT_Axes(Radius), fj::contrib::NormalizedMeasure(Beta,fR));
1196  Result= nSub.result(fFilteredJets[0]);
1197  SubJet_Axes=nSub.currentAxes();
1198  }
1199  else if (Algorithm==8){
1200  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_WTA_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1201  Result= nSub.result(fFilteredJets[0]);
1202  SubJet_Axes=nSub.currentAxes();
1203  }
1204  else if (Algorithm==9){
1205  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_WTA_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,fR));
1206  Result= nSub.result(fFilteredJets[0]);
1207  SubJet_Axes=nSub.currentAxes();
1208  }
1209  else if (Algorithm==10){
1210  fj::contrib::Nsubjettiness nSub(N, fj::contrib::MultiPass_Axes(100), fj::contrib::NormalizedMeasure(Beta,fR));
1211  Result= nSub.result(fFilteredJets[0]);
1212  SubJet_Axes=nSub.currentAxes();
1213  }
1214 
1215  SubJet1_Axis=SubJet_Axes[0];
1216  Double_t SubJet1_Eta=SubJet1_Axis.pseudorapidity();
1217  Double_t SubJet2_Eta;
1218  Double_t SubJet1_Phi=SubJet1_Axis.phi();
1219  if(SubJet1_Phi < -1*TMath::Pi()) SubJet1_Phi += (2*TMath::Pi());
1220  else if (SubJet1_Phi > TMath::Pi()) SubJet1_Phi -= (2*TMath::Pi());
1221  Double_t SubJet2_Phi;
1222  Double_t DeltaPhi=-5;
1223  if (SubJet_Axes.size()>1){
1224  SubJet2_Axis=SubJet_Axes[1];
1225  SubJet2_Eta=SubJet2_Axis.pseudorapidity();
1226  SubJet2_Phi=SubJet2_Axis.phi();
1227  if(SubJet2_Phi < -1*TMath::Pi()) SubJet2_Phi += (2*TMath::Pi());
1228  else if (SubJet2_Phi > TMath::Pi()) SubJet2_Phi -= (2*TMath::Pi());
1229  DeltaPhi=SubJet1_Phi-SubJet2_Phi;
1230  if(DeltaPhi < -1*TMath::Pi()) DeltaPhi += (2*TMath::Pi());
1231  else if (DeltaPhi > TMath::Pi()) DeltaPhi -= (2*TMath::Pi());
1232  }
1233 
1234  if (Option==0) return Result;
1235  else if (Option==1 && SubJet_Axes.size()>1 && N==2) return TMath::Sqrt(TMath::Power(SubJet1_Eta-SubJet2_Eta,2)+TMath::Power(DeltaPhi,2));
1236  else if (Option==2 && SubJet_Axes.size()>1 && N==2) return TMath::Sqrt(TMath::Power(SubJet1_Eta-SubJet2_Eta,2)+TMath::Power(DeltaPhi,2));
1237  else return -2;
1238 }
1239 
1240 
1241 
1242 Double32_t AliFJWrapper::NSubjettinessDerivativeSub(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Double_t JetR, fastjet::PseudoJet jet, Int_t Option){ //For derivative subtraction
1243 
1244  Double_t Result=-1;
1245  std::vector<fastjet::PseudoJet> SubJet_Axes;
1246  fj::PseudoJet SubJet1_Axis;
1247  fj::PseudoJet SubJet2_Axis;
1248  if (Algorithm==0){
1249  fj::contrib::Nsubjettiness nSub(N, fj::contrib::KT_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1250  Result= nSub.result(jet);
1251  SubJet_Axes=nSub.currentAxes();
1252  }
1253  else if (Algorithm==1) {
1254  fj::contrib::Nsubjettiness nSub(N, fj::contrib::CA_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1255  Result= nSub.result(jet);
1256  SubJet_Axes=nSub.currentAxes();
1257  }
1258  else if (Algorithm==2){
1259  fj::contrib::Nsubjettiness nSub(N, fj::contrib::AntiKT_Axes(Radius), fj::contrib::NormalizedMeasure(Beta,JetR));
1260  Result= nSub.result(jet);
1261  SubJet_Axes=nSub.currentAxes();
1262  }
1263  else if (Algorithm==3) {
1264  fj::contrib::Nsubjettiness nSub(N, fj::contrib::WTA_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1265  Result= nSub.result(jet);
1266  SubJet_Axes=nSub.currentAxes();
1267  }
1268  else if (Algorithm==4) {
1269  fj::contrib::Nsubjettiness nSub(N, fj::contrib::WTA_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1270  Result= nSub.result(jet);
1271  SubJet_Axes=nSub.currentAxes();
1272  }
1273  else if (Algorithm==5){
1274  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1275  Result= nSub.result(jet);
1276  SubJet_Axes=nSub.currentAxes();
1277  }
1278  else if (Algorithm==6){
1279  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1280  Result= nSub.result(jet);
1281  SubJet_Axes=nSub.currentAxes();
1282  }
1283  else if (Algorithm==7){
1284  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_AntiKT_Axes(Radius), fj::contrib::NormalizedMeasure(Beta,JetR));
1285  Result= nSub.result(jet);
1286  SubJet_Axes=nSub.currentAxes();
1287  }
1288  else if (Algorithm==8){
1289  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_WTA_KT_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1290  Result= nSub.result(jet);
1291  SubJet_Axes=nSub.currentAxes();
1292  }
1293  else if (Algorithm==9){
1294  fj::contrib::Nsubjettiness nSub(N, fj::contrib::OnePass_WTA_CA_Axes(), fj::contrib::NormalizedMeasure(Beta,JetR));
1295  Result= nSub.result(jet);
1296  SubJet_Axes=nSub.currentAxes();
1297  }
1298  else if (Algorithm==10){
1299  fj::contrib::Nsubjettiness nSub(N, fj::contrib::MultiPass_Axes(100), fj::contrib::NormalizedMeasure(Beta,JetR));
1300  Result= nSub.result(jet);
1301  SubJet_Axes=nSub.currentAxes();
1302  }
1303 
1304  SubJet1_Axis=SubJet_Axes[0];
1305  Double_t SubJet1_Eta=SubJet1_Axis.pseudorapidity();
1306  Double_t SubJet2_Eta;
1307  Double_t SubJet1_Phi=SubJet1_Axis.phi();
1308  if(SubJet1_Phi < -1*TMath::Pi()) SubJet1_Phi += (2*TMath::Pi());
1309  else if (SubJet1_Phi > TMath::Pi()) SubJet1_Phi -= (2*TMath::Pi());
1310  Double_t SubJet2_Phi;
1311  Double_t DeltaPhi=-5;
1312  if (SubJet_Axes.size()>1){
1313  SubJet2_Axis=SubJet_Axes[1];
1314  SubJet2_Eta=SubJet2_Axis.pseudorapidity();
1315  SubJet2_Phi=SubJet2_Axis.phi();
1316  if(SubJet2_Phi < -1*TMath::Pi()) SubJet2_Phi += (2*TMath::Pi());
1317  else if (SubJet2_Phi > TMath::Pi()) SubJet2_Phi -= (2*TMath::Pi());
1318  DeltaPhi=SubJet1_Phi-SubJet2_Phi;
1319  if(DeltaPhi < -1*TMath::Pi()) DeltaPhi += (2*TMath::Pi());
1320  else if (DeltaPhi > TMath::Pi()) DeltaPhi -= (2*TMath::Pi());
1321  }
1322 
1323  if (Option==0) return Result;
1324  else if (Option==1 && SubJet_Axes.size()>1 && N==2) return TMath::Sqrt(TMath::Power(SubJet1_Eta-SubJet2_Eta,2)+TMath::Power(DeltaPhi,2));
1325  else if (Option==2 && SubJet_Axes.size()>1 && N==2) return TMath::Sqrt(TMath::Power(SubJet1_Eta-SubJet2_Eta,2)+TMath::Power(DeltaPhi,2));
1326  else return -2;
1327 
1328 }
1329 
1330 
1331 
1332 
1333 
1334 
1335 
1336 
1337 
1338 #endif
fastjet::PseudoJet GetJetAreaVector(UInt_t idx) const
std::vector< double > fGRNumeratorSub
Definition: AliFJWrapper.h:182
void SetRMaxAndStep(Double_t rmax, Double_t dr)
Definition: AliFJWrapper.h:111
void SetUseExternalBkg(Bool_t b, Double_t rho, Double_t rhom)
Definition: AliFJWrapper.h:110
Double_t GetFilteredJetArea(UInt_t idx) const
std::vector< fastjet::PseudoJet > fInputVectors
Definition: AliFJWrapper.h:118
Double_t fR
Definition: AliFJWrapper.h:144
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:135
Int_t fPluginAlgor
Definition: AliFJWrapper.h:150
fastjet::ClusterSequenceActiveAreaExplicitGhosts * fClustSeqActGhosts
Definition: AliFJWrapper.h:136
virtual Int_t DoGenericSubtractionJetOpeningAngle_kt()
Double_t fKtScatter
Definition: AliFJWrapper.h:148
Int_t fNGhostRepeats
Definition: AliFJWrapper.h:141
virtual Int_t DoGenericSubtractionJetpTD()
Double_t GetMedianUsedForBgSubtraction() const
Definition: AliFJWrapper.h:37
virtual void AddInputVectors(const std::vector< fastjet::PseudoJet > &vecs, Int_t offsetIndex=-99999)
std::vector< fastjet::PseudoJet > fInclusiveJets
Definition: AliFJWrapper.h:120
fastjet::JetAlgorithm fAlgor
Definition: AliFJWrapper.h:138
const char * GetTitle() const
Definition: AliFJWrapper.h:39
fastjet::RecombinationScheme fScheme
Definition: AliFJWrapper.h:139
Double_t fRMax
Definition: AliFJWrapper.h:178
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 SetMinJetPt(Double_t MinPt)
Definition: AliFJWrapper.h:113
void SetMaxRap(Double_t maxrap)
Definition: AliFJWrapper.h:97
Double_t GetJetSubtractedPt(UInt_t idx) const
Double_t fMinJetPt
Definition: AliFJWrapper.h:145
fastjet::PseudoJet GetFilteredJetAreaVector(UInt_t idx) const
void SetRecombScheme(const fastjet::RecombinationScheme &scheme)
Definition: AliFJWrapper.h:93
fastjet::RangeDefinition * fRange
Definition: AliFJWrapper.h:130
fastjet::ClusterSequenceActiveAreaExplicitGhosts * GetClusterSequenceGhosts() const
Definition: AliFJWrapper.h:30
fastjet::ClusterSequence * GetClusterSequenceSA() const
Definition: AliFJWrapper.h:29
void SetStrategy(const fastjet::Strategy &strat)
Definition: AliFJWrapper.h:91
std::vector< fastjet::PseudoJet > fFilteredJets
Definition: AliFJWrapper.h:121
void SetupAreaTypefromOpt(const char *option)
virtual void ClearMemory()
fastjet::AreaType fAreaType
Definition: AliFJWrapper.h:140
void SetR(Double_t r)
Definition: AliFJWrapper.h:98
Double_t fMaxRap
Definition: AliFJWrapper.h:143
fastjet::ClusterSequenceArea * GetClusterSequence() const
Definition: AliFJWrapper.h:28
Double_t * sigma
const std::vector< fastjet::PseudoJet > & GetInclusiveJets() const
Definition: AliFJWrapper.h:33
AliFJWrapper & operator=(const AliFJWrapper &wrapper)
Double_t GetJetArea(UInt_t idx) const
std::vector< double > fSubtractedJetsPt
Definition: AliFJWrapper.h:122
fastjet::JetDefinition::Plugin * fPlugin
Definition: AliFJWrapper.h:128
void SetLegacyMode(Bool_t mode)
Definition: AliFJWrapper.h:108
Bool_t fDoFilterArea
Definition: AliFJWrapper.h:173
void SetAlgorithm(const fastjet::JetAlgorithm &algor)
Definition: AliFJWrapper.h:92
virtual Int_t DoGenericSubtractionJetCircularity()
void SetupSchemefromOpt(const char *option)
const std::vector< fastjet::PseudoJet > & GetInputGhosts() const
Definition: AliFJWrapper.h:32
void SetUseArea4Vector(Bool_t useA4v)
Definition: AliFJWrapper.h:103
fastjet::ClusterSequenceArea * fClustSeq
Definition: AliFJWrapper.h:134
void SetName(const char *name)
Definition: AliFJWrapper.h:89
virtual void SubtractBackground(const Double_t median_pt=-1)
std::vector< fastjet::PseudoJet > fInputGhosts
Definition: AliFJWrapper.h:119
Double_t fGridScatter
Definition: AliFJWrapper.h:147
virtual void Clear(const Option_t *="")
Bool_t fLegacyMode
Definition: AliFJWrapper.h:174
fastjet::VoronoiAreaSpec * fVorAreaSpec
Definition: AliFJWrapper.h:125
Int_t mode
Definition: anaM.C:40
Double_t fDRStep
Definition: AliFJWrapper.h:179
std::vector< double > fGRDenominatorSub
Definition: AliFJWrapper.h:183
Bool_t GetLegacyMode()
Definition: AliFJWrapper.h:46
fastjet::GhostedAreaSpec * fGhostedAreaSpec
Definition: AliFJWrapper.h:126
Bool_t fUseArea4Vector
Definition: AliFJWrapper.h:153
void SetNRepeats(Int_t nrepeat)
Definition: AliFJWrapper.h:95
virtual std::vector< double > GetGRNumeratorSub() const
Definition: AliFJWrapper.h:68
virtual ~AliFJWrapper()
std::vector< double > fGRDenominator
Definition: AliFJWrapper.h:181
virtual Int_t DoGenericSubtractionJet1subjettiness_kt()
Double_t fGhostArea
Definition: AliFJWrapper.h:142
virtual void GetMedianAndSigma(Double_t &median, Double_t &sigma, Int_t remove=0) const
TString fName
Definition: AliFJWrapper.h:116
void SetGhostArea(Double_t gharea)
Definition: AliFJWrapper.h:96
Double_t fMedUsedForBgSub
Definition: AliFJWrapper.h:152
std::vector< fastjet::PseudoJet > GetJetConstituents(UInt_t idx) const
fastjet::JetDefinition * fJetDef
Definition: AliFJWrapper.h:127
Double_t NSubjettiness(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Int_t Option=0)
virtual Int_t DoGenericSubtractionJetSigma2()
void SetupStrategyfromOpt(const char *option)
const std::vector< fastjet::PseudoJet > & GetInputVectors() const
Definition: AliFJWrapper.h:31
fastjet::AreaDefinition * fAreaDef
Definition: AliFJWrapper.h:124
virtual std::vector< double > GetGRDenominatorSub() const
Definition: AliFJWrapper.h:69
Bool_t fUseExternalBkg
Definition: AliFJWrapper.h:175
void SetTitle(const char *title)
Definition: AliFJWrapper.h:90
void SetPluginAlgor(Int_t plugin)
Definition: AliFJWrapper.h:102
virtual Int_t DoConstituentSubtraction()
Double_t fRhom
Definition: AliFJWrapper.h:177
virtual std::vector< double > GetGRDenominator() const
Definition: AliFJWrapper.h:67
void SetRhoRhom(Double_t rho, Double_t rhom)
Definition: AliFJWrapper.h:112
void SetGridScatter(Double_t gridSc)
Definition: AliFJWrapper.h:99
std::vector< fastjet::PseudoJet > fConstituentSubtrJets
Definition: AliFJWrapper.h:123
void SetMeanGhostKt(Double_t meankt)
Definition: AliFJWrapper.h:101
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)
virtual Int_t DoGenericSubtractionJet3subjettiness_kt()
void SetLegacyFJ()
std::vector< fastjet::PseudoJet > GetFilteredJetConstituents(UInt_t idx) const
Double32_t NSubjettinessDerivativeSub(Int_t N, Int_t Algorithm, Double_t Radius, Double_t Beta, Double_t JetR, fastjet::PseudoJet jet, Int_t Option=0)
void SetKtScatter(Double_t ktSc)
Definition: AliFJWrapper.h:100
virtual std::vector< double > GetGRNumerator() const
Definition: AliFJWrapper.h:66
virtual const char * ClassName() const
Definition: AliFJWrapper.h:23
virtual void RemoveLastInputVector()
void SetAreaType(const fastjet::AreaType &atype)
Definition: AliFJWrapper.h:94
Double_t fRho
Definition: AliFJWrapper.h:176
Bool_t GetDoFilterArea()
Definition: AliFJWrapper.h:47
virtual Int_t DoGenericSubtractionJetLeSub()
virtual Int_t DoGenericSubtractionJetConstituent()
const char * GetName() const
Definition: AliFJWrapper.h:38
const std::vector< fastjet::PseudoJet > & GetFilteredJets() const
Definition: AliFJWrapper.h:34
virtual Int_t DoGenericSubtractionJet2subjettiness_kt()
std::vector< double > fGRNumerator
Definition: AliFJWrapper.h:180
Double_t fMeanGhostKt
Definition: AliFJWrapper.h:149
fastjet::Strategy fStrategy
Definition: AliFJWrapper.h:137
TString fTitle
Definition: AliFJWrapper.h:117