FastJet 3.0.2
|
Class to estimate the pt density of the background per unit area, using the median of the distribution of pt/area from jets that pass some selection criterion. More...
#include <fastjet/tools/JetMedianBackgroundEstimator.hh>
Public Member Functions | |
constructors and destructors | |
JetMedianBackgroundEstimator (const Selector &rho_range, const JetDefinition &jet_def, const AreaDefinition &area_def) | |
Constructor that sets the rho range as well as the jet definition and area definition to be used to cluster the particles. | |
JetMedianBackgroundEstimator (const Selector &rho_range, const ClusterSequenceAreaBase &csa) | |
ctor from a ClusterSequenceAreaBase with area | |
JetMedianBackgroundEstimator (const Selector &rho_range=SelectorIdentity()) | |
Default constructor that optionally sets the rho range. | |
~JetMedianBackgroundEstimator () | |
default dtor | |
setting a new event | |
virtual void | set_particles (const std::vector< PseudoJet > &particles) |
tell the background estimator that it has a new event, composed of the specified particles. | |
void | set_cluster_sequence (const ClusterSequenceAreaBase &csa) |
(re)set the cluster sequence (with area support) to be used by future calls to rho() etc. | |
void | set_jets (const std::vector< PseudoJet > &jets) |
(re)set the jets (which must have area support) to be used by future calls to rho() etc. | |
void | set_selector (const Selector &rho_range_selector) |
(re)set the selector to be used for future calls to rho() etc. | |
retrieving fundamental information | |
double | rho () const |
get rho, the median background density per unit area | |
double | sigma () const |
get sigma, the background fluctuations per unit area | |
double | rho (const PseudoJet &jet) |
get rho, the median background density per unit area, locally at the position of a given jet. | |
double | sigma (const PseudoJet &jet) |
get sigma, the background fluctuations per unit area, locally at the position of a given jet. | |
virtual bool | has_sigma () |
returns true if this background estimator has support for determination of sigma | |
retrieving additional useful information | |
double | mean_area () const |
Returns the mean area of the jets used to actually compute the background properties in the last call of rho() or sigma() | |
unsigned int | n_jets_used () const |
returns the number of jets used to actually compute the background properties in the last call of rho() or sigma() | |
double | empty_area () const |
Returns the estimate of the area (within the range defined by the selector) that is not occupied by jets. | |
double | n_empty_jets () const |
Returns the number of empty jets used when computing the background properties. | |
configuring behaviour | |
void | reset () |
Resets the class to its default state, including the choice to use 4-vector areas. | |
void | set_use_area_4vector (bool use_it=true) |
By default when calculating pt/Area for a jet, it is the transverse component of the 4-vector area that is used in the ratiof . | |
bool | use_area_4vector () const |
check if the estimator uses the 4-vector area or the scalar area | |
void | set_provide_fj2_sigma (bool provide_fj2_sigma=true) |
The FastJet v2.X sigma calculation had a small spurious offset in the limit of a small number of jets. | |
void | set_jet_density_class (const FunctionOfPseudoJet< double > *jet_density_class) |
Set a pointer to a class that calculates the quantity whose median will be calculated; if the pointer is null then pt/area is used (as occurs also if this function is not called). | |
const FunctionOfPseudoJet < double > * | jet_density_class () const |
return the pointer to the jet density class | |
virtual void | set_rescaling_class (const FunctionOfPseudoJet< double > *rescaling_class_in) |
Set a pointer to a class that calculates the rescaling factor as a function of the jet (position). | |
description | |
std::string | description () const |
returns a textual description of the background estimator |
Class to estimate the pt density of the background per unit area, using the median of the distribution of pt/area from jets that pass some selection criterion.
Events are passed either in the form of the event particles (in which they're clustered by the class), a ClusterSequenceArea (in which case the jets used are those returned by "inclusive_jets()") or directly as a set of jets.
The selection criterion is typically a geometrical one (e.g. all jets with |y|<2) sometimes supplemented with some kinematical restriction (e.g. exclusion of the two hardest jets). It is passed to the class through a Selector.
Beware: by default, to correctly handle partially empty events, the class attempts to calculate an "empty area", based (schematically) on
range.total_area() - sum_{jets_in_range} jets.area()
For ranges with small areas, this can be inaccurate (particularly relevant in dense events where empty_area should be zero and ends up not being zero).
This calculation of empty area can be avoided if a ClusterSequenceArea class with explicit ghosts (ActiveAreaExplicitGhosts) is used. This is _recommended_ unless speed requirements cause you to use Voronoi areas. For speedy background estimation you could also consider using GridMedianBackgroundEstimator.
Definition at line 78 of file JetMedianBackgroundEstimator.hh.
fastjet::JetMedianBackgroundEstimator::JetMedianBackgroundEstimator | ( | const Selector & | rho_range, |
const JetDefinition & | jet_def, | ||
const AreaDefinition & | area_def | ||
) |
Constructor that sets the rho range as well as the jet definition and area definition to be used to cluster the particles.
Prior to the estimation of rho, one has to provide the particles to cluster using set_particles(...)
rho_range | the Selector specifying which jets will be considered |
jet_def | the jet definition to use for the clustering |
area_def | the area definition to use for the clustering |
Definition at line 71 of file JetMedianBackgroundEstimator.cc.
fastjet::JetMedianBackgroundEstimator::JetMedianBackgroundEstimator | ( | const Selector & | rho_range, |
const ClusterSequenceAreaBase & | csa | ||
) |
ctor from a ClusterSequenceAreaBase with area
rho_range | the Selector specifying which jets will be considered |
csa | the ClusterSequenceArea to use |
Pre-conditions:
Note that selectors with e.g. hardest-jets exclusion do not have a well-defined area. For this reasons, it is STRONGLY advised to use an area with explicit ghosts.
Definition at line 88 of file JetMedianBackgroundEstimator.cc.
fastjet::JetMedianBackgroundEstimator::JetMedianBackgroundEstimator | ( | const Selector & | rho_range = SelectorIdentity() | ) | [inline] |
Default constructor that optionally sets the rho range.
The configuration must be done later calling set_cluster_sequence(...) or set_jets(...).
rho_range | the Selector specifying which jets will be considered |
Definition at line 122 of file JetMedianBackgroundEstimator.hh.
void fastjet::JetMedianBackgroundEstimator::set_particles | ( | const std::vector< PseudoJet > & | particles | ) | [virtual] |
tell the background estimator that it has a new event, composed of the specified particles.
Implements fastjet::BackgroundEstimatorBase.
Definition at line 40 of file GridMedianBackgroundEstimator.cc.
void fastjet::JetMedianBackgroundEstimator::set_cluster_sequence | ( | const ClusterSequenceAreaBase & | csa | ) |
(re)set the cluster sequence (with area support) to be used by future calls to rho() etc.
csa | the cluster sequence area |
Pre-conditions:
Note that selectors with e.g. hardest-jets exclusion do not have a well-defined area. For this reasons, it is STRONGLY advised to use an area with explicit ghosts.
Definition at line 147 of file JetMedianBackgroundEstimator.cc.
void fastjet::JetMedianBackgroundEstimator::set_jets | ( | const std::vector< PseudoJet > & | jets | ) |
(re)set the jets (which must have area support) to be used by future calls to rho() etc.
; for the conditions that must be satisfied by the jets, see the Constructor that takes jets.
Definition at line 171 of file JetMedianBackgroundEstimator.cc.
double fastjet::JetMedianBackgroundEstimator::rho | ( | const PseudoJet & | jet | ) | [virtual] |
get rho, the median background density per unit area, locally at the position of a given jet.
If the Selector associated with the range takes a reference jet (i.e. is relocatable), then for subsequent operations the Selector has that jet set as its reference.
Implements fastjet::BackgroundEstimatorBase.
Definition at line 237 of file JetMedianBackgroundEstimator.cc.
double fastjet::JetMedianBackgroundEstimator::sigma | ( | const PseudoJet & | jet | ) | [virtual] |
get sigma, the background fluctuations per unit area, locally at the position of a given jet.
If the Selector associated with the range takes a reference jet (i.e. is relocatable), then for subsequent operations the Selector has that jet set as its reference.
Reimplemented from fastjet::BackgroundEstimatorBase.
Definition at line 252 of file JetMedianBackgroundEstimator.cc.
double fastjet::JetMedianBackgroundEstimator::empty_area | ( | ) | const [inline] |
Returns the estimate of the area (within the range defined by the selector) that is not occupied by jets.
The value is that for the last call of rho() or sigma()
The answer is defined to be zero if the area calculation involved explicit ghosts; if the area calculation was an active area, then use is made of the active area's internal list of pure ghost jets (taking those that pass the selector); otherwise it is based on the difference between the selector's total area and the area of the jets that pass the selector.
The result here is just the cached result of the corresponding call to the ClusterSequenceAreaBase function.
Definition at line 235 of file JetMedianBackgroundEstimator.hh.
double fastjet::JetMedianBackgroundEstimator::n_empty_jets | ( | ) | const [inline] |
Returns the number of empty jets used when computing the background properties.
The value is that for the last call of rho() or sigma().
If the area has explicit ghosts the result is zero; for active areas it is the number of internal pure ghost jets that pass the selector; otherwise it is deduced from the empty area, divided by (the average pure-ghost-jet area).
The result here is just the cached result of the corresponding call to the ClusterSequenceAreaBase function.
Definition at line 251 of file JetMedianBackgroundEstimator.hh.
void fastjet::JetMedianBackgroundEstimator::set_use_area_4vector | ( | bool | use_it = true | ) | [inline] |
By default when calculating pt/Area for a jet, it is the transverse component of the 4-vector area that is used in the ratiof .
Calling this function with a "false" argument causes the scalar area to be used instead.
While the difference between the two choices is usually small, for high-precision work it is usually the 4-vector area that is to be preferred.
use_it | whether one uses the 4-vector area or not (true by default) |
Definition at line 278 of file JetMedianBackgroundEstimator.hh.
void fastjet::JetMedianBackgroundEstimator::set_provide_fj2_sigma | ( | bool | provide_fj2_sigma = true | ) | [inline] |
The FastJet v2.X sigma calculation had a small spurious offset in the limit of a small number of jets.
This is fixed by default in versions 3 upwards. The old behaviour can be obtained with a call to this function.
Definition at line 290 of file JetMedianBackgroundEstimator.hh.
void fastjet::JetMedianBackgroundEstimator::set_jet_density_class | ( | const FunctionOfPseudoJet< double > * | jet_density_class | ) |
Set a pointer to a class that calculates the quantity whose median will be calculated; if the pointer is null then pt/area is used (as occurs also if this function is not called).
Note that this is still preliminary in FastJet 3.0 and that backward compatibility is not guaranteed in future releases of FastJet
Definition at line 290 of file JetMedianBackgroundEstimator.cc.
virtual void fastjet::JetMedianBackgroundEstimator::set_rescaling_class | ( | const FunctionOfPseudoJet< double > * | rescaling_class_in | ) | [inline, virtual] |
Set a pointer to a class that calculates the rescaling factor as a function of the jet (position).
Note that the rescaling factor is used both in the determination of the "global" rho (the pt/A of each jet is divided by this factor) and when asking for a local rho (the result is multiplied by this factor).
The BackgroundRescalingYPolynomial class can be used to get a rescaling that depends just on rapidity.
Reimplemented from fastjet::BackgroundEstimatorBase.
Definition at line 317 of file JetMedianBackgroundEstimator.hh.