FastJet 3.0.2
Public Member Functions | Protected Member Functions
fastjet::ClusterSequenceAreaBase Class Reference

base class that sets interface for extensions of ClusterSequence that provide information about the area of each jet More...

#include <fastjet/ClusterSequenceAreaBase.hh>

Inheritance diagram for fastjet::ClusterSequenceAreaBase:
Inheritance graph
[legend]
Collaboration diagram for fastjet::ClusterSequenceAreaBase:
Collaboration graph
[legend]

List of all members.

Public Member Functions

template<class L >
 ClusterSequenceAreaBase (const std::vector< L > &pseudojets, const JetDefinition &jet_def_in, const bool &writeout_combinations=false)
 a constructor which just carries out the construction of the parent class
 ClusterSequenceAreaBase ()
 default constructor
virtual ~ClusterSequenceAreaBase ()
 destructor
virtual double area (const PseudoJet &) const
 return the area associated with the given jet; this base class returns 0.
virtual double area_error (const PseudoJet &) const
 return the error (uncertainty) associated with the determination of the area of this jet; this base class returns 0.
virtual PseudoJet area_4vector (const PseudoJet &) const
 return a PseudoJet whose 4-vector is defined by the following integral
virtual bool is_pure_ghost (const PseudoJet &) const
 true if a jet is made exclusively of ghosts
virtual bool has_explicit_ghosts () const
 returns true if ghosts are explicitly included within jets for this ClusterSequence;
virtual double empty_area (const Selector &selector) const
 return the total area, corresponding to the given Selector, that is free of jets, in general based on the inclusive jets.
double empty_area_from_jets (const std::vector< PseudoJet > &all_jets, const Selector &selector) const
 return the total area, corresponding to the given Selector, that is free of jets, based on the supplied all_jets
virtual double n_empty_jets (const Selector &selector) const
 return something similar to the number of pure ghost jets in the given selector's range in an active area case.
double median_pt_per_unit_area (const Selector &selector) const
 the median of (pt/area) for jets contained within the selector range, making use also of the info on n_empty_jets
double median_pt_per_unit_area_4vector (const Selector &selector) const
 the median of (pt/area_4vector) for jets contained within the selector range, making use also of the info on n_empty_jets
double median_pt_per_unit_something (const Selector &selector, bool use_area_4vector) const
 the function that does the work for median_pt_per_unit_area and median_pt_per_unit_area_4vector:
virtual void get_median_rho_and_sigma (const Selector &selector, bool use_area_4vector, double &median, double &sigma, double &mean_area) const
 using jets withing the selector range (and with 4-vector areas if use_area_4vector), calculate the median pt/area, as well as an "error" (uncertainty), which is defined as the 1-sigma half-width of the distribution of pt/A, obtained by looking for the point below which we have (1-0.6827)/2 of the jets (including empty jets).
virtual void get_median_rho_and_sigma (const std::vector< PseudoJet > &all_jets, const Selector &selector, bool use_area_4vector, double &median, double &sigma, double &mean_area, bool all_are_inclusive=false) const
 a more advanced version of get_median_rho_and_sigma, which allows one to use any "view" of the event containing all jets (so that, e.g.
virtual void get_median_rho_and_sigma (const Selector &selector, bool use_area_4vector, double &median, double &sigma) const
 same as the full version of get_median_rho_and_error, but without access to the mean_area
virtual void parabolic_pt_per_unit_area (double &a, double &b, const Selector &selector, double exclude_above=-1.0, bool use_area_4vector=false) const
 fits a form pt_per_unit_area(y) = a + b*y^2 in the selector range.
std::vector< PseudoJetsubtracted_jets (const double rho, const double ptmin=0.0) const
 return a vector of all subtracted jets, using area_4vector, given rho.
std::vector< PseudoJetsubtracted_jets (const Selector &selector, const double ptmin=0.0) const
 return a vector of subtracted jets, using area_4vector.
PseudoJet subtracted_jet (const PseudoJet &jet, const double rho) const
 return a subtracted jet, using area_4vector, given rho
PseudoJet subtracted_jet (const PseudoJet &jet, const Selector &selector) const
 return a subtracted jet, using area_4vector; note that this is potentially inefficient if repeatedly used for many different jets, because rho will be recalculated each time around.
double subtracted_pt (const PseudoJet &jet, const double rho, bool use_area_4vector=false) const
 return the subtracted pt, given rho
double subtracted_pt (const PseudoJet &jet, const Selector &selector, bool use_area_4vector=false) const
 return the subtracted pt; note that this is potentially inefficient if repeatedly used for many different jets, because rho will be recalculated each time around.

Protected Member Functions

void _check_selector_good_for_median (const Selector &selector) const
 check the selector is suited for the computations i.e. applies jet by jet and has a finite area

Detailed Description

base class that sets interface for extensions of ClusterSequence that provide information about the area of each jet

the virtual functions here all return 0, since no area determination is implemented.

Definition at line 45 of file ClusterSequenceAreaBase.hh.


Member Function Documentation

virtual double fastjet::ClusterSequenceAreaBase::area ( const PseudoJet ) const [inline, virtual]

return the area associated with the given jet; this base class returns 0.

Reimplemented in fastjet::ClusterSequenceActiveArea, fastjet::ClusterSequenceActiveAreaExplicitGhosts, fastjet::ClusterSequenceArea, and fastjet::ClusterSequenceVoronoiArea.

Definition at line 67 of file ClusterSequenceAreaBase.hh.

virtual double fastjet::ClusterSequenceAreaBase::area_error ( const PseudoJet ) const [inline, virtual]

return the error (uncertainty) associated with the determination of the area of this jet; this base class returns 0.

Reimplemented in fastjet::ClusterSequenceActiveArea, fastjet::ClusterSequenceArea, and fastjet::ClusterSequenceVoronoiArea.

Definition at line 71 of file ClusterSequenceAreaBase.hh.

virtual PseudoJet fastjet::ClusterSequenceAreaBase::area_4vector ( const PseudoJet ) const [inline, virtual]

return a PseudoJet whose 4-vector is defined by the following integral

drap d PseudoJet("rap,phi,pt=one") * Theta("rap,phi inside jet boundary")

where PseudoJet("rap,phi,pt=one") is a 4-vector with the given rapidity (rap), azimuth (phi) and pt=1, while Theta("rap,phi inside jet boundary") is a function that is 1 when rap,phi define a direction inside the jet boundary and 0 otherwise.

This base class returns a null 4-vector.

Reimplemented in fastjet::ClusterSequenceActiveArea, fastjet::ClusterSequenceActiveAreaExplicitGhosts, fastjet::ClusterSequenceArea, and fastjet::ClusterSequenceVoronoiArea.

Definition at line 84 of file ClusterSequenceAreaBase.hh.

virtual bool fastjet::ClusterSequenceAreaBase::is_pure_ghost ( const PseudoJet ) const [inline, virtual]

true if a jet is made exclusively of ghosts

NB: most area classes do not give any explicit ghost jets, but some do, and they should replace this function with their own version.

Reimplemented in fastjet::ClusterSequenceActiveAreaExplicitGhosts, and fastjet::ClusterSequenceArea.

Definition at line 92 of file ClusterSequenceAreaBase.hh.

virtual bool fastjet::ClusterSequenceAreaBase::has_explicit_ghosts ( ) const [inline, virtual]

returns true if ghosts are explicitly included within jets for this ClusterSequence;

Derived classes that do include explicit ghosts should provide an alternative version of this routine and set it properly.

Reimplemented in fastjet::ClusterSequenceActiveAreaExplicitGhosts, and fastjet::ClusterSequenceArea.

Definition at line 101 of file ClusterSequenceAreaBase.hh.

double fastjet::ClusterSequenceAreaBase::empty_area ( const Selector selector) const [virtual]

return the total area, corresponding to the given Selector, that is free of jets, in general based on the inclusive jets.

return the total area, within the selector's range, that is free of jets.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Calculate this as (range area) - {i in range} A_i

for ClusterSequences with explicit ghosts, assume that there will never be any empty area, i.e. it is always filled in by pure ghosts jets. This holds for seq.rec. algorithms

Reimplemented in fastjet::ClusterSequenceActiveArea, fastjet::ClusterSequenceActiveAreaExplicitGhosts, fastjet::ClusterSequenceArea, and fastjet::ClusterSequencePassiveArea.

Definition at line 55 of file ClusterSequenceAreaBase.cc.

double fastjet::ClusterSequenceAreaBase::empty_area_from_jets ( const std::vector< PseudoJet > &  all_jets,
const Selector selector 
) const

return the total area, corresponding to the given Selector, that is free of jets, based on the supplied all_jets

return the total area, within range, that is free of jets.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Calculate this as (range area) - {i in range} A_i

Definition at line 67 of file ClusterSequenceAreaBase.cc.

virtual double fastjet::ClusterSequenceAreaBase::n_empty_jets ( const Selector selector) const [inline, virtual]

return something similar to the number of pure ghost jets in the given selector's range in an active area case.

For the local implementation we return empty_area/(0.55 pi R^2), based on measured properties of ghost jets with kt and cam (cf arXiv:0802.1188).

Note that the number returned is a double.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Reimplemented in fastjet::ClusterSequence1GhostPassiveArea, fastjet::ClusterSequenceActiveArea, and fastjet::ClusterSequenceArea.

Definition at line 133 of file ClusterSequenceAreaBase.hh.

double fastjet::ClusterSequenceAreaBase::median_pt_per_unit_area ( const Selector selector) const

the median of (pt/area) for jets contained within the selector range, making use also of the info on n_empty_jets

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Definition at line 79 of file ClusterSequenceAreaBase.cc.

double fastjet::ClusterSequenceAreaBase::median_pt_per_unit_area_4vector ( const Selector selector) const

the median of (pt/area_4vector) for jets contained within the selector range, making use also of the info on n_empty_jets

The selector passed as an argument has to have a finite area and apply jet-by-jet

Definition at line 83 of file ClusterSequenceAreaBase.cc.

double fastjet::ClusterSequenceAreaBase::median_pt_per_unit_something ( const Selector selector,
bool  use_area_4vector 
) const

the function that does the work for median_pt_per_unit_area and median_pt_per_unit_area_4vector:

the median of (pt/area) for jets contained within range, counting the empty area as if it were made up of a collection of empty jets each of area (0.55 * pi R^2).

  • something_is_area_4vect = false -> use plain area
  • something_is_area_4vect = true -> use 4-vector area

Definition at line 92 of file ClusterSequenceAreaBase.cc.

void fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma ( const Selector selector,
bool  use_area_4vector,
double &  median,
double &  sigma,
double &  mean_area 
) const [virtual]

using jets withing the selector range (and with 4-vector areas if use_area_4vector), calculate the median pt/area, as well as an "error" (uncertainty), which is defined as the 1-sigma half-width of the distribution of pt/A, obtained by looking for the point below which we have (1-0.6827)/2 of the jets (including empty jets).

The subtraction for a jet with uncorrected pt pt^U and area A is

pt^S = pt^U - median*A +- sigma*sqrt(A)

where the error is only that associated with the fluctuations in the noise and not that associated with the noise having caused changes in the hard-particle content of the jet.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

NB: subtraction may also be done with 4-vector area of course, and this is recommended for jets with larger values of R, as long as rho has also been determined with a 4-vector area; using a scalar area causes one to neglect terms of relative order $R^2/8$ in the jet $p_t$.

Reimplemented in fastjet::ClusterSequenceArea.

Definition at line 162 of file ClusterSequenceAreaBase.cc.

virtual void fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma ( const std::vector< PseudoJet > &  all_jets,
const Selector selector,
bool  use_area_4vector,
double &  median,
double &  sigma,
double &  mean_area,
bool  all_are_inclusive = false 
) const [virtual]

a more advanced version of get_median_rho_and_sigma, which allows one to use any "view" of the event containing all jets (so that, e.g.

one might use Cam on a different resolution scale without have to rerun the algorithm).

By default it will assume that "all" are not inclusive jets, so that in dealing with empty area it has to calculate the number of empty jets based on the empty area and the the observed <area> of jets rather than a surmised area

Note that for small effective radii, this can cause problems because the harder jets get an area >> <ghost-jet-area> and so the estimate comes out all wrong. In these situations it is highly advisable to use an area with explicit ghosts, since then the "empty" jets are actually visible.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Reimplemented in fastjet::ClusterSequenceArea.

virtual void fastjet::ClusterSequenceAreaBase::get_median_rho_and_sigma ( const Selector selector,
bool  use_area_4vector,
double &  median,
double &  sigma 
) const [inline, virtual]

same as the full version of get_median_rho_and_error, but without access to the mean_area

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Reimplemented in fastjet::ClusterSequenceArea.

Definition at line 221 of file ClusterSequenceAreaBase.hh.

void fastjet::ClusterSequenceAreaBase::parabolic_pt_per_unit_area ( double &  a,
double &  b,
const Selector selector,
double  exclude_above = -1.0,
bool  use_area_4vector = false 
) const [virtual]

fits a form pt_per_unit_area(y) = a + b*y^2 in the selector range.

fits a form pt_per_unit_area(y) = a + b*y^2 for jets in range.

exclude_above allows one to exclude large values of pt/area from fit. (if negative, the cut is discarded) use_area_4vector = true uses the 4vector areas.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

exclude_above allows one to exclude large values of pt/area from fit. use_area_4vector = true uses the 4vector areas.

Reimplemented in fastjet::ClusterSequenceArea.

Definition at line 106 of file ClusterSequenceAreaBase.cc.

vector< PseudoJet > fastjet::ClusterSequenceAreaBase::subtracted_jets ( const double  rho,
const double  ptmin = 0.0 
) const

return a vector of all subtracted jets, using area_4vector, given rho.

Only inclusive_jets above ptmin are subtracted and returned. the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()), i.e. not necessarily ordered in pt once subtracted

Definition at line 292 of file ClusterSequenceAreaBase.cc.

vector< PseudoJet > fastjet::ClusterSequenceAreaBase::subtracted_jets ( const Selector selector,
const double  ptmin = 0.0 
) const

return a vector of subtracted jets, using area_4vector.

Only inclusive_jets above ptmin are subtracted and returned. the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()), i.e. not necessarily ordered in pt once subtracted

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Only inclusive_jets above ptmin are subtracted and returned. the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()), i.e. not necessarily ordered in pt once subtracted

Definition at line 308 of file ClusterSequenceAreaBase.cc.

PseudoJet fastjet::ClusterSequenceAreaBase::subtracted_jet ( const PseudoJet jet,
const Selector selector 
) const

return a subtracted jet, using area_4vector; note that this is potentially inefficient if repeatedly used for many different jets, because rho will be recalculated each time around.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Definition at line 341 of file ClusterSequenceAreaBase.cc.

double fastjet::ClusterSequenceAreaBase::subtracted_pt ( const PseudoJet jet,
const Selector selector,
bool  use_area_4vector = false 
) const

return the subtracted pt; note that this is potentially inefficient if repeatedly used for many different jets, because rho will be recalculated each time around.

The selector passed as an argument has to have a finite area and apply jet-by-jet (see the BackgroundEstimator and Subtractor tools for more generic usages)

Definition at line 365 of file ClusterSequenceAreaBase.cc.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends