29 #include "fastjet/tools/Pruner.hh" 30 #include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh" 31 #include "fastjet/Selector.hh" 40 FASTJET_BEGIN_NAMESPACE
55 : _jet_def(jet_def), _zcut(0), _Rcut_factor(0),
56 _zcut_dyn(zcut_dyn), _Rcut_dyn(Rcut_dyn), _get_recombiner_from_jet(false) {
57 assert(_zcut_dyn != 0 && _Rcut_dyn != 0);
65 throw Error(
"Pruner: trying to apply the Pruner transformer to a jet that has no constituents");
70 bool do_areas = jet.
has_area() && _check_explicit_ghosts(jet);
73 double Rcut = (_Rcut_dyn) ? (*_Rcut_dyn)(jet) : _Rcut_factor * 2.0*jet.
m()/jet.
perp();
74 double zcut = (_zcut_dyn) ? (*_zcut_dyn)(jet) : _zcut;
75 PruningPlugin * pruning_plugin;
79 if (_get_recombiner_from_jet) {
81 _get_common_recombiner(jet);
82 if (common_recombiner) {
91 pruning_plugin =
new PruningPlugin(jet_def, zcut, Rcut);
95 pruning_plugin =
new PruningPlugin(_jet_def, zcut, Rcut);
98 pruning_plugin =
new PruningPlugin(_jet_def, zcut, Rcut);
110 vector<PseudoJet> particles, ghosts;
115 double ghost_area = (ghosts.size()) ? ghosts[0].area() : 0.01;
138 bool Pruner::_check_explicit_ghosts(
const PseudoJet &jet)
const{
146 vector<PseudoJet> pieces = jet.
pieces();
147 for (
unsigned int i=0;i<pieces.size(); i++)
148 if (!_check_explicit_ghosts(pieces[i]))
return false;
173 vector<PseudoJet> pieces = jet.
pieces();
174 if (pieces.size() == 0)
return 0;
176 for (
unsigned int i=1;i<pieces.size(); i++)
177 if (_get_common_recombiner(pieces[i]) != reco)
return 0;
190 oss <<
"Pruner with jet_definition = (" << _jet_def.
description() <<
")";
192 oss <<
", dynamic zcut (" << _zcut_dyn->
description() <<
")" 193 <<
", dynamic Rcut (" << _Rcut_dyn->
description() <<
")";
195 oss <<
", zcut = " << _zcut
196 <<
", Rcut_factor = " << _Rcut_factor;
220 void PruningRecombiner::recombine(
const PseudoJet &pa,
224 _recombiner->recombine(pa, pb, p);
231 double pt2a = pa.
perp2();
232 double pt2b = pb.
perp2();
236 if (pt2a<_zcut2*p.
perp2()){
242 if (pt2b<_zcut2*p.
perp2()) {
251 string PruningRecombiner::description()
const{
253 oss <<
"Pruning recombiner with zcut = " << sqrt(_zcut2)
254 <<
", Rcut = " << sqrt(_Rcut2)
255 <<
", and underlying recombiner = " << _recombiner->
description();
268 PruningRecombiner pruning_recombiner(_zcut, _Rcut, _jet_def.
recombiner());
274 const vector<ClusterSequence::history_element> & internal_hist = internal_cs.
history();
277 vector<bool> kept(internal_hist.size(),
true);
278 const vector<unsigned int> &pr_rej = pruning_recombiner.rejected();
279 for (
unsigned int i=0;i<pr_rej.size(); i++) kept[pr_rej[i]]=
false;
285 vector<unsigned int> internal2input(internal_hist.size());
286 for (
unsigned int i=0; i<input_cs.
jets().size(); i++)
287 internal2input[i] = i;
289 for (
unsigned int i=input_cs.
jets().size(); i<internal_hist.size(); i++){
293 if (he.
parent2 == ClusterSequence::BeamJet){
294 int internal_jetp_index = internal_hist[he.parent1].
jetp_index;
295 int internal_hist_index = internal_cs.jets()[internal_jetp_index].cluster_hist_index();
297 int input_jetp_index = input_cs.
history()[internal2input[internal_hist_index]].jetp_index;
307 if (!kept[he.parent1]){
308 internal2input[i]=internal2input[he.
parent2];
314 internal2input[i]=internal2input[he.parent1];
324 internal2input[i]=input_cs.
jets()[new_index].cluster_hist_index();
333 string PruningPlugin::description()
const{
335 oss <<
"Pruning plugin with jet_definition = (" << _jet_def.
description()
336 <<
"), zcut = " << _zcut
337 <<
", Rcut = " << _Rcut;
342 FASTJET_END_NAMESPACE
const JetDefinition & jet_def() const
return a reference to the jet definition
void set_recombination_scheme(RecombinationScheme)
set the recombination scheme to the one provided
const ClusterSequenceAreaBase * validated_csab() const
shorthand for validated_cluster_sequence_area_base()
virtual PseudoJet result(const PseudoJet &jet) const
action on a single jet
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
virtual bool has_area() const
check if it has a defined area
Selector SelectorNHardest(unsigned int n)
select the n hardest objects
const std::vector< history_element > & history() const
allow the user to access the raw internal history.
bool has_associated_cluster_sequence() const
returns true if this PseudoJet has an associated ClusterSequence.
std::vector< PseudoJet > extra_jets() const
return the other jets that may have been found along with the result of the pruning The resulting vec...
void delete_plugin_when_unused()
allows to let the JetDefinition handle the deletion of the plugin when it is no longer used ...
const ClusterSequence * validated_cs() const
shorthand for validated_cluster_sequence()
double squared_distance(const PseudoJet &other) const
returns squared cylinder (rap-phi) distance between this jet and another
int jetp_index
index in _history where the current jet is recombined with another jet to form its child...
void plugin_record_ij_recombination(int jet_i, int jet_j, double dij, int &newjet_k)
record the fact that there has been a recombination between jets()[jet_i] and jets()[jet_k], with the specified dij, and return the index (newjet_k) allocated to the new jet, whose momentum is assumed to be the 4-vector sum of that of jet_i and jet_j
virtual bool has_pieces() const
returns true if a jet has pieces
int parent2
index in _history where first parent of this jet was created (InexistentParent if this jet is an orig...
double dij
index in the _jets vector where we will find the
A class that will provide the recombination scheme facilities and/or allow a user to extend these fac...
void plugin_record_iB_recombination(int jet_i, double diB)
record the fact that there has been a recombination between jets()[jet_i] and the beam...
virtual std::vector< PseudoJet > constituents() const
retrieve the constituents.
double perp() const
returns the scalar transverse momentum
const Recombiner * recombiner() const
return a pointer to the currently defined recombiner.
The structure associated with a PseudoJet thas has gone through a Pruner transformer.
virtual std::string description() const
description
virtual std::vector< PseudoJet > pieces() const
retrieve the pieces that make up the jet.
base class corresponding to errors that can be thrown by FastJet
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP) ...
Like ClusterSequence with computation of the active jet area with the addition of explicit ghosts...
int cluster_hist_index() const
return the cluster_hist_index, intended to be used by clustering routines.
RecombinationScheme
the various recombination schemes
double perp2() const
returns the squared transverse momentum
an implementation of C++0x shared pointers (or boost's)
virtual bool has_constituents() const
returns true if the PseudoJet has constituents
virtual bool has_explicit_ghosts() const
returns true if ghosts are explicitly included within jets for this ClusterSequence; ...
void sift(const std::vector< PseudoJet > &jets, std::vector< PseudoJet > &jets_that_pass, std::vector< PseudoJet > &jets_that_fail) const
sift the input jets into two vectors – those that pass the selector and those that do not ...
a single element in the clustering history
void set_recombiner(const Recombiner *recomb)
set the recombiner class to the one provided
std::vector< PseudoJet > inclusive_jets(const double &ptmin=0.0) const
return a vector of all jets (in the sense of the inclusive algorithm) with pt >= ptmin.
std::string description() const
return a string describing what kind of PseudoJet we are dealing with
Selector SelectorIsPureGhost()
select objects that are (or are only made of) ghosts.
std::string description() const
return a textual description of the current jet definition
An abstract base class that will provide the recombination scheme facilities and/or allow a user to e...
void delete_self_when_unused()
by calling this routine you tell the ClusterSequence to delete itself when all the Pseudojets associa...
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
const std::vector< PseudoJet > & jets() const
allow the user to access the internally stored _jets() array, which contains both the initial particl...
virtual std::string description() const
returns a description of the function (an empty string by default)
class that is intended to hold a full definition of the jet clusterer