FastJet 3.0.2
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes
fastjet::ClusterSequenceActiveArea Class Reference

Like ClusterSequence with computation of the active jet area. More...

#include <fastjet/ClusterSequenceActiveArea.hh>

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

List of all members.

Public Types

enum  mean_pt_strategies {
  median = 0, non_ghost_median, pttot_over_areatot, pttot_over_areatot_cut,
  mean_ratio_cut, play, median_4vector
}
 enum providing a variety of tentative strategies for estimating the background (e.g. More...

Public Member Functions

 ClusterSequenceActiveArea ()
 default constructor
template<class L >
 ClusterSequenceActiveArea (const std::vector< L > &pseudojets, const JetDefinition &jet_def_in, const GhostedAreaSpec &ghost_spec, const bool &writeout_combinations=false)
 constructor based on JetDefinition and GhostedAreaSpec
virtual double area (const PseudoJet &jet) const
 return the area associated with the given jet; this base class returns 0.
virtual double area_error (const PseudoJet &jet) 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 &jet) const
 return a PseudoJet whose 4-vector is defined by the following integral
double pt_per_unit_area (mean_pt_strategies strat=median, double range=2.0) const
 return the transverse momentum per unit area according to one of the above strategies; for some strategies (those with "cut" in their name) the parameter "range" allows one to exclude a subset of the jets for the background estimation, those that have pt/area > median(pt/area)*range.
virtual double empty_area (const Selector &selector) const
 rewrite the empty area from the parent class, so as to use all info at our disposal return the total area, corresponding to a given Selector, that consists of ghost jets or unclustered ghosts
virtual double n_empty_jets (const Selector &selector) const
 return the true number of empty jets (replaces ClusterSequenceAreaBase::n_empty_jets(...))

Protected Member Functions

void _resize_and_zero_AA ()
void _initialise_AA (const JetDefinition &jet_def, const GhostedAreaSpec &ghost_spec, const bool &writeout_combinations, bool &continue_running)
void _run_AA (const GhostedAreaSpec &ghost_spec)
void _postprocess_AA (const GhostedAreaSpec &ghost_spec)
 run the postprocessing for the active area (and derived classes)
void _initialise_and_run_AA (const JetDefinition &jet_def, const GhostedAreaSpec &ghost_spec, const bool &writeout_combinations=false)
 does the initialisation and running specific to the active areas class
void _transfer_ghost_free_history (const ClusterSequenceActiveAreaExplicitGhosts &clust_seq)
 transfer the history (and jet-momenta) from clust_seq to our own internal structure while removing ghosts
void _transfer_areas (const std::vector< int > &unique_hist_order, const ClusterSequenceActiveAreaExplicitGhosts &)
 transfer areas from the ClusterSequenceActiveAreaExplicitGhosts object into our internal area bookkeeping...
bool has_dangerous_particles () const
 returns true if there are any particles whose transverse momentum if so low that there's a risk of the ghosts having modified the clustering sequence

Protected Attributes

std::valarray< double > _average_area
 child classes benefit from having these at their disposal
std::valarray< double > _average_area2
std::valarray< PseudoJet_average_area_4vector

Detailed Description

Like ClusterSequence with computation of the active jet area.

Class that behaves essentially like ClusterSequence except that it also provides access to the area of a jet (which will be a random quantity... Figure out what to do about seeds later...)

This class should not be used directly. Rather use ClusterSequenceArea with the appropriate AreaDefinition

Definition at line 61 of file ClusterSequenceActiveArea.hh.


Member Enumeration Documentation

enum providing a variety of tentative strategies for estimating the background (e.g.

non-jet) activity in a highly populated event; the one that has been most extensively tested is median.

Definition at line 85 of file ClusterSequenceActiveArea.hh.


Member Function Documentation

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

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

Reimplemented from fastjet::ClusterSequenceAreaBase.

Definition at line 74 of file ClusterSequenceActiveArea.hh.

virtual double fastjet::ClusterSequenceActiveArea::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 from fastjet::ClusterSequenceAreaBase.

Definition at line 76 of file ClusterSequenceActiveArea.hh.

virtual PseudoJet fastjet::ClusterSequenceActiveArea::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 from fastjet::ClusterSequenceAreaBase.

Definition at line 79 of file ClusterSequenceActiveArea.hh.

double fastjet::ClusterSequenceActiveArea::pt_per_unit_area ( mean_pt_strategies  strat = median,
double  range = 2.0 
) const

return the transverse momentum per unit area according to one of the above strategies; for some strategies (those with "cut" in their name) the parameter "range" allows one to exclude a subset of the jets for the background estimation, those that have pt/area > median(pt/area)*range.

NB: This call is OBSOLETE; use media_pt_per_unit_area from the

Definition at line 274 of file ClusterSequenceActiveArea.cc.

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

rewrite the empty area from the parent class, so as to use all info at our disposal return the total area, corresponding to a given Selector, that consists of ghost jets or unclustered ghosts

The selector passed as an argument needs to apply jet by jet.

Reimplemented from fastjet::ClusterSequenceAreaBase.

Reimplemented in fastjet::ClusterSequencePassiveArea.

Definition at line 443 of file ClusterSequenceActiveArea.cc.

void fastjet::ClusterSequenceActiveArea::_initialise_and_run_AA ( const JetDefinition jet_def,
const GhostedAreaSpec ghost_spec,
const bool &  writeout_combinations = false 
) [protected]

does the initialisation and running specific to the active areas class

global routine for running active area

Definition at line 51 of file ClusterSequenceActiveArea.cc.

void fastjet::ClusterSequenceActiveArea::_transfer_areas ( const std::vector< int > &  unique_hist_order,
const ClusterSequenceActiveAreaExplicitGhosts  
) [protected]

transfer areas from the ClusterSequenceActiveAreaExplicitGhosts object into our internal area bookkeeping...

Definition at line 564 of file ClusterSequenceActiveArea.cc.


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