FastJet 3.0alpha3
Protected Member Functions
fastjet::BackgroundEstimator Class Reference

Class to estimate the density of the background per unit area. More...

#include <BackgroundEstimator.hh>

Collaboration diagram for fastjet::BackgroundEstimator:
Collaboration graph
[legend]

List of all members.

Public Member Functions

constructors and destructors
 BackgroundEstimator (const ClusterSequenceAreaBase &csa, const Selector &rho_range)
 ctor from a ClusterSequenceAreaBase with area
 BackgroundEstimator (const std::vector< PseudoJet > &jets, const Selector &rho_range)
 ctor from a list of jets
 BackgroundEstimator (const Selector &rho_range)
 Constructor that just sets the rho range.
 BackgroundEstimator ()
 default ctor.
 ~BackgroundEstimator ()
 default dtor
double rho () const
 @ name retrieving fundamental information
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.
double mean_area () const
 @ name retrieving additional useful information
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 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.
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 $p_t/A$.
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 ()
 return the pointer to the jet density class
void set_rescaling_class (const FunctionOfPseudoJet< double > *rescaling_class)
 Set a pointer to a class that calculates the rescaling factor as a function of the jet (position).
const FunctionOfPseudoJet
< double > * 
rescaling_class ()
 return the pointer to the jet density class

Protected Member Functions

void _median_and_stddev (const std::vector< double > &quantity_vector, double n_empty_jets, double &median, double &stand_dev_if_gaussian, bool do_fj2_calculation=false) const
 given a quantity in a vector (e.g.

Detailed Description

Class to estimate the density of the background per unit area.

For a given event, this class calculated the median of the distribution of pt/Area for all jets in the event that pass some selection criterion.

Events are passed either in the form of a ClusterSequenceArea (in which case the jets used as 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).

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 innacurate (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 you supply a ClusterSequenceArea class with explicit ghosts (ActiveAreaExplicitGhosts). This is _recommended_!

Definition at line 154 of file BackgroundEstimator.hh.


Constructor & Destructor Documentation

fastjet::BackgroundEstimator::BackgroundEstimator ( const ClusterSequenceAreaBase csa,
const Selector rho_range 
)

ctor from a ClusterSequenceAreaBase with area

Parameters:
csathe ClusterSequenceArea to use
rho_rangea Selector that specifies the range over which jets will be considered

Pre-conditions:

  • one should be able to estimate the "empty area" (i.e. the area not occupied by jets). This is feasible if at least one of the following conditions is satisfied: ( i) the ClusterSequence has explicit ghosts (ii) the range has a computable area.
  • the jet algorithm must be suited for median computation (otherwise a warning will be issues)

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 92 of file BackgroundEstimator.cc.

References reset(), and set_cluster_sequence().

  : _rho_range(rho_range){

  // initialise things properly
  reset();

  // tell the BGE about the cluster sequence
  set_cluster_sequence(csa);
}
fastjet::BackgroundEstimator::BackgroundEstimator ( const std::vector< PseudoJet > &  jets,
const Selector rho_range 
)

ctor from a list of jets

Parameters:
jetsthe jets to use
rho_rangethe range over which jets will be considered

Pre-conditions:

  • all the jets must come from the same cluster sequence
  • As for the above ctor, one needs to be able to estimate the empty area. The conditions that the ClusterSequence must have explicit ghosts is checked on the ClusterSequence shared by the jets.
  • As for the above ctor, the jet alg must be adequate
fastjet::BackgroundEstimator::BackgroundEstimator ( const Selector rho_range) [inline]

Constructor that just sets the rho range.

The CSA or jets actually used for estimating the background must be passed later via set_cluster_sequence(...) or set_jets(...)

Parameters:
rho_rangethe range over which jets will be considered

Definition at line 199 of file BackgroundEstimator.hh.

                                                 : _rho_range(rho_range) {
    reset();
  }
fastjet::BackgroundEstimator::BackgroundEstimator ( ) [inline]

default ctor.

In order to use the BG estimator you then need to call set_selector(...) and also set_cluster_sequence(...) or set_jets(...)

Definition at line 206 of file BackgroundEstimator.hh.

{reset();}

Member Function Documentation

double fastjet::BackgroundEstimator::rho ( ) const [inline]

@ name retrieving fundamental information

get rho, the median background density oer unit area

Definition at line 219 of file BackgroundEstimator.hh.

                     {
    if (_rho_range.takes_reference())
      throw Error("The background estimation is obtained from a selector that takes a reference jet. rho(PseudoJet) should be used in that case");
    _recompute_if_needed();
    return _rho;
  }
double fastjet::BackgroundEstimator::rho ( const PseudoJet  jet) [inline]

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.

Definition at line 240 of file BackgroundEstimator.hh.

                                  {
    _recompute_if_needed(jet);
    double our_rho = _rho;
    if (_rescaling_class != 0) { 
      our_rho *= (*_rescaling_class)(jet);
    }
    return our_rho;
  }
double fastjet::BackgroundEstimator::sigma ( const PseudoJet jet) [inline]

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.

Definition at line 255 of file BackgroundEstimator.hh.

                                     {
    _recompute_if_needed(jet);
    double our_sigma = _sigma;
    if (_rescaling_class != 0) { 
      our_sigma *= (*_rescaling_class)(jet);
    }
    return our_sigma;
  }
double fastjet::BackgroundEstimator::mean_area ( ) const [inline]

@ name retrieving additional useful information

Returns the mean area of the jets used to actually compute the background properties in the last call of rho() or sigma()

Definition at line 271 of file BackgroundEstimator.hh.

                          {
    _recompute_if_needed();
    return _mean_area;
  }
double fastjet::BackgroundEstimator::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 296 of file BackgroundEstimator.hh.

                           {
    _recompute_if_needed();
    return _empty_area;
  }
double fastjet::BackgroundEstimator::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 $ 0.55 \pi R^2 $ (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 312 of file BackgroundEstimator.hh.

                             {
    _recompute_if_needed();
    return _n_empty_jets;
  }
void fastjet::BackgroundEstimator::set_cluster_sequence ( const ClusterSequenceAreaBase csa)

(re)set the cluster sequence (with area support) to be used by future calls to rho() etc.

Parameters:
csathe cluster sequence area

Pre-conditions:

  • one should be able to estimate the "empty area" (i.e. the area not occupied by jets). This is feasible if at least one of the following conditions is satisfied: ( i) the ClusterSequence has explicit ghosts (ii) the range has a computable area.
  • the jet algorithm must be suited for median computation (otherwise a warning will be issues)

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 123 of file BackgroundEstimator.cc.

References fastjet::ClusterSequenceAreaBase::has_explicit_ghosts(), fastjet::Selector::has_finite_area(), fastjet::ClusterSequence::inclusive_jets(), and fastjet::ClusterSequence::structure_shared_ptr().

Referenced by BackgroundEstimator().

                                                                                  {
  _csi = csa.structure_shared_ptr();

  // sanity checks
  //---------------
  //  (i) check the alg is appropriate
  _check_jet_alg_good_for_median();

  //  (ii) check that, if there are no explicit ghosts, the selector has a finite area
  if ((!csa.has_explicit_ghosts()) && (!_rho_range.has_finite_area())){
    throw Error("BackgroundEstimator: either an area with explicit ghosts (recommended) or a Selector with finite area is needed (to allow for the computation of the empty area)");
  }

  // get the initial list of jets
  _included_jets = csa.inclusive_jets();

  _uptodate = false;
}
void fastjet::BackgroundEstimator::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 144 of file BackgroundEstimator.cc.

References fastjet::SharedPtr< T >::get(), fastjet::ClusterSequenceAreaBase::has_explicit_ghosts(), fastjet::Selector::has_finite_area(), and fastjet::ClusterSequenceStructure::validated_csab().

                                                                {
  
  if (! jets.size())
    throw Error("BackgroundEstimator::BackgroundEstimator: At least one jet is needed to compute the background properties");

  // sanity checks
  //---------------
  //  (o) check that there is an underlying CS shared by all the jets
  if (! (jets[0].has_associated_cluster_sequence()) && (jets[0].has_area()))
    throw Error("BackgroundEstimator::BackgroundEstimator: the jets used to estimate the background properties must be associated with a valid ClusterSequenceAreaBase");

  _csi = jets[0].structure_shared_ptr();
  ClusterSequenceStructure * csi = dynamic_cast<ClusterSequenceStructure*>(_csi());
  const ClusterSequenceAreaBase * csab = csi->validated_csab();

  for (unsigned int i=1;i<jets.size(); i++){
    if (! jets[i].has_associated_cluster_sequence()) // area automatic if the next test succeeds
      throw Error("BackgroundEstimator::set_jets(...): the jets used to estimate the background properties must be associated with a valid ClusterSequenceAreaBase");

    if (jets[i].structure_shared_ptr().get() != _csi.get())
      throw Error("BackgroundEstimator::set_jets(...): all the jets used to estimate the background properties must share the same ClusterSequence");
  }

  //  (i) check the alg is appropriate
  _check_jet_alg_good_for_median();

  //  (ii) check that, if there are no explicit ghosts, the selector has a finite area
  if ((!csab->has_explicit_ghosts()) && (!_rho_range.has_finite_area())){
    throw Error("BackgroundEstimator: either an area with explicit ghosts (recommended) or a Selector with finite area is needed (to allow for the computation of the empty area)");
  }


  // get the initial list of jets
  _included_jets = jets;

  // ensure recalculation of quantities that need it
  _uptodate = false;
}
void fastjet::BackgroundEstimator::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 $p_t/A$.

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.

Parameters:
use_itwhether one uses the 4-vector area or not (true by default)

Definition at line 369 of file BackgroundEstimator.hh.

Referenced by reset().

                                               {
    _use_area_4vector = use_it;
    _uptodate = false;
  }  
void fastjet::BackgroundEstimator::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 381 of file BackgroundEstimator.hh.

Referenced by reset().

                                                            {
    _provide_fj2_sigma = provide_fj2_sigma;
    _uptodate = false;
  }
void fastjet::BackgroundEstimator::set_jet_density_class ( const FunctionOfPseudoJet< double > *  jet_density_class) [inline]

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).

Definition at line 389 of file BackgroundEstimator.hh.

                                                                                    {
    _jet_density_class = jet_density_class;
    _uptodate = false;
  }
void fastjet::BackgroundEstimator::set_rescaling_class ( const FunctionOfPseudoJet< double > *  rescaling_class) [inline]

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).

Definition at line 404 of file BackgroundEstimator.hh.

                                                                                {
    _rescaling_class = rescaling_class;
    _uptodate = false;
  }
void fastjet::BackgroundEstimator::_median_and_stddev ( const std::vector< double > &  quantity_vector,
double  n_empty_jets,
double &  median,
double &  stand_dev_if_gaussian,
bool  do_fj2_calculation = false 
) const [protected]

given a quantity in a vector (e.g.

pt_over_area) and knowledge about the number of empty jets, calculate the median and stand_dev_if_gaussian (roughly from the 16th percentile)

If do_fj2_calculation is set to true then this performs FastJet 2.X estimation of the standard deviation, which has a spurious offset in the limit of a small number of jets.

Definition at line 297 of file BackgroundEstimator.cc.

References n_jets_used().

                                                                            {

  // this check is redundant (the code below behaves sensibly even
  // with a zero size), but serves as a reminder of what happens if
  // the quantity vector is zero-sized
  if (quantity_vector.size() == 0) {
    median = 0;
    stand_dev_if_gaussian = 0;
    return;
  }

  vector<double> sorted_quantity_vector = quantity_vector;
  sort(sorted_quantity_vector.begin(), sorted_quantity_vector.end());

  // now get the median & error, accounting for empty jets
  // define the fractions of distribution at median, median-1sigma
  double posn[2] = {0.5, (1.0-0.6827)/2.0};
  double res[2];

  int n_jets_used = sorted_quantity_vector.size();
  double total_njets = n_jets_used + _n_empty_jets;

  for (int i = 0; i < 2; i++) {
    double nj_median_pos;
    if (do_fj2_calculation) {
      nj_median_pos = (total_njets-1)*posn[i] - n_empty_jets;
    } else {
      nj_median_pos = (total_njets)*posn[i] - n_empty_jets - 0.5;
    }

    double nj_median_ratio;
    if (nj_median_pos >= 0 && sorted_quantity_vector.size() > 1) {
      int int_nj_median = int(nj_median_pos);
      nj_median_ratio =
        sorted_quantity_vector[int_nj_median] * (int_nj_median+1-nj_median_pos)
        + sorted_quantity_vector[int_nj_median+1] * (nj_median_pos - int_nj_median);
    } else if (nj_median_pos > -0.5 && sorted_quantity_vector.size() >= 1 && !do_fj2_calculation) {
      // in the LHS of this "bin", just keep a constant value (we could have
      // interpolated to zero, but this might misbehave in cases where all jets
      // are active, because it would go to zero too fast)
      nj_median_ratio = sorted_quantity_vector[0];
    } else {
      nj_median_ratio = 0.0;
    }
    res[i] = nj_median_ratio;
  }
  
  median = res[0];
  stand_dev_if_gaussian = res[0] - res[1];
}

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