1 #ifndef __FASTJET_TOOLS_FILTER_HH__ 2 #define __FASTJET_TOOLS_FILTER_HH__ 32 #include <fastjet/ClusterSequence.hh> 33 #include <fastjet/Selector.hh> 34 #include <fastjet/CompositeJetStructure.hh> 35 #include <fastjet/tools/Transformer.hh> 39 FASTJET_BEGIN_NAMESPACE
43 class FilterStructure;
115 _subjet_def(subjet_def), _Rfiltfunc(0), _Rfilt(-1), _selector(selector), _rho(rho), _subtractor(0) {}
123 _Rfiltfunc(0), _Rfilt(Rfilt), _selector(selector), _rho(rho), _subtractor(0) {
125 throw Error(
"Attempt to create a Filter with a negative filtering radius");
135 _Rfiltfunc(Rfilt_func), _Rfilt(-1), _selector(selector), _rho(rho), _subtractor(0) {}
152 virtual std::string description()
const;
161 void _set_filtered_elements(
const PseudoJet & jet,
162 std::vector<PseudoJet> & filtered_elements,
164 bool & discard_area)
const;
167 void _set_filtered_elements_cafilt(
const PseudoJet & jet,
168 std::vector<PseudoJet> & filtered_elements,
172 void _set_filtered_elements_generic(
const PseudoJet & jet,
173 std::vector<PseudoJet> & filtered_elements,
175 bool do_areas)
const;
180 std::vector<PseudoJet> & kept,
181 std::vector<PseudoJet> & rejected,
183 const bool discard_area)
const;
188 bool _get_all_pieces(
const PseudoJet &jet, std::vector<PseudoJet> &all_pieces)
const;
194 bool _check_ca(
const std::vector<PseudoJet> &all_pieces)
const;
201 bool _check_explicit_ghosts(
const std::vector<PseudoJet> &all_pieces)
const;
203 bool _uses_subtraction()
const {
return (_subtractor || _rho != 0);}
233 virtual std::string
description()
const {
return "Filtered PseudoJet"; }
245 const std::vector<PseudoJet> &
rejected()
const {
return _rejected;}
255 FASTJET_END_NAMESPACE
257 #endif // __FASTJET_TOOLS_FILTER_HH__ Filter(JetDefinition subjet_def, Selector selector, double rho=0.0)
define a filter that decomposes a jet into subjets using a generic JetDefinition and then keeps only ...
Class that helps perform filtering (Butterworth, Davison, Rubin and Salam, arXiv:0802.2470) and trimming (Krohn, Thaler and Wang, arXiv:0912.1342) on jets, optionally in conjunction with subtraction (Cacciari and Salam, arXiv:0707.1378).
The structure for a jet made of pieces.
FilterStructure(const std::vector< PseudoJet > &pieces_in, const JetDefinition::Recombiner *rec=0)
constructor from an original ClusterSequenceInfo We just share the original ClusterSequenceWrapper an...
Filter(double Rfilt, Selector selector, double rho=0.0)
Same as the full constructor (see above) but just specifying the radius By default, Cambridge-Aachen is used If the jet (or all its pieces) is obtained with a non-default recombiner, that one will be used.
Filter()
trivial ctor Note: this is just for derived classes a Filter initialised through this constructor wil...
base class corresponding to errors that can be thrown by FastJet
Class to contain structure information for a filtered jet.
virtual std::string description() const
description
const std::vector< PseudoJet > & rejected() const
returns the subjets that were not kept during the filtering procedure (subtracted if the filter reque...
virtual ~FilterStructure()
virtual dtor to allow further overloading
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Filter(FunctionOfPseudoJet< double > *Rfilt_func, Selector selector, double rho=0.0)
Same as the full constructor (see above) but just specifying a filtering radius that will depend on t...
virtual ~Filter()
default dtor
std::vector< PseudoJet > _rejected
the subjets rejected by the filter
An abstract base class that will provide the recombination scheme facilities and/or allow a user to e...
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
class that is intended to hold a full definition of the jet clusterer
void set_subtractor(const Transformer *subtractor)
Set a subtractor that is applied to all individual subjets before deciding which ones to keep...