32 #include "fastjet/Selector.hh"
33 #include "fastjet/GhostedAreaSpec.hh"
37 FASTJET_BEGIN_NAMESPACE
44 std::vector<PseudoJet> Selector::operator()(
const std::vector<PseudoJet> & jets)
const {
45 std::vector<PseudoJet> result;
50 for (std::vector<PseudoJet>::const_iterator jet = jets.begin();
51 jet != jets.end(); jet++) {
52 if (worker_local->
pass(*jet)) result.push_back(*jet);
57 std::vector<const PseudoJet *> jetptrs(jets.size());
58 for (
unsigned i = 0; i < jets.size(); i++) {
59 jetptrs[i] = & jets[i];
62 for (
unsigned i = 0; i < jetptrs.size(); i++) {
63 if (jetptrs[i]) result.push_back(jets[i]);
72 unsigned int Selector::count(
const std::vector<PseudoJet> & jets)
const {
78 for (
unsigned i = 0; i < jets.size(); i++) {
79 if (worker_local->
pass(jets[i])) n++;
82 std::vector<const PseudoJet *> jetptrs(jets.size());
83 for (
unsigned i = 0; i < jets.size(); i++) {
84 jetptrs[i] = & jets[i];
87 for (
unsigned i = 0; i < jetptrs.size(); i++) {
99 void Selector::sift(
const std::vector<PseudoJet> & jets,
100 std::vector<PseudoJet> & jets_that_pass,
101 std::vector<PseudoJet> & jets_that_fail
105 jets_that_pass.clear();
106 jets_that_fail.clear();
110 for (
unsigned i = 0; i < jets.size(); i++) {
111 if (worker_local->
pass(jets[i])) {
112 jets_that_pass.push_back(jets[i]);
114 jets_that_fail.push_back(jets[i]);
118 std::vector<const PseudoJet *> jetptrs(jets.size());
119 for (
unsigned i = 0; i < jets.size(); i++) {
120 jetptrs[i] = & jets[i];
123 for (
unsigned i = 0; i < jetptrs.size(); i++) {
125 jets_that_pass.push_back(jets[i]);
127 jets_that_fail.push_back(jets[i]);
134 double Selector::area()
const{
135 return area(gas::def_ghost_area);
139 double Selector::area(
double ghost_area)
const{
143 if (_worker->has_known_area())
return _worker->known_area();
146 double rapmin, rapmax;
147 get_rapidity_extent(rapmin, rapmax);
149 std::vector<PseudoJet> ghosts;
153 return ghost_spec.ghost_area() * ((*this)(ghosts)).size();
161 bool SelectorWorker::has_finite_area()
const {
162 if (! is_geometric())
return false;
163 double rapmin, rapmax;
164 get_rapidity_extent(rapmin, rapmax);
165 return (rapmax != std::numeric_limits<double>::infinity())
166 && (-rapmin != std::numeric_limits<double>::infinity());
183 virtual bool pass(
const PseudoJet &)
const {
189 virtual void terminator(vector<const PseudoJet *> &)
const {
195 virtual string description()
const {
return "Identity";}
198 virtual bool is_geometric()
const {
return true;}
203 Selector SelectorIdentity() {
204 return Selector(
new SW_Identity);
214 class SW_Not :
public SelectorWorker {
217 SW_Not(
const Selector & s) : _s(s) {}
220 virtual SelectorWorker* copy(){
return new SW_Not(*
this);}
224 virtual bool pass(
const PseudoJet & jet)
const {
226 if (!applies_jet_by_jet())
227 throw Error(
"Cannot apply this selector worker to an individual jet");
229 return ! _s.pass(jet);
233 virtual bool applies_jet_by_jet()
const {
return _s.applies_jet_by_jet();}
236 virtual void terminator(vector<const PseudoJet *> & jets)
const {
239 if (applies_jet_by_jet()){
240 SelectorWorker::terminator(jets);
245 vector<const PseudoJet *> s_jets = jets;
246 _s.worker()->terminator(s_jets);
250 for (
unsigned int i=0; i<s_jets.size(); i++){
251 if (s_jets[i]) jets[i] = NULL;
256 virtual string description()
const {
258 ostr <<
"!(" << _s.description() <<
")";
263 virtual bool is_geometric()
const {
return _s.is_geometric();}
266 virtual bool takes_reference()
const {
return _s.takes_reference();}
269 virtual void set_reference(
const PseudoJet &ref) { _s.set_reference(ref);}
284 class SW_BinaryOperator:
public SelectorWorker {
287 SW_BinaryOperator(
const Selector & s1,
const Selector & s2) : _s1(s1), _s2(s2) {
291 _applies_jet_by_jet = _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
294 _takes_reference = _s1.takes_reference() || _s2.takes_reference();
297 _is_geometric = _s1.is_geometric() && _s2.is_geometric();
301 virtual bool applies_jet_by_jet()
const {
return _applies_jet_by_jet;}
304 virtual bool takes_reference()
const{
305 return _takes_reference;
309 virtual void set_reference(
const PseudoJet ¢re){
310 _s1.set_reference(centre);
311 _s2.set_reference(centre);
315 virtual bool is_geometric()
const {
return _is_geometric;}
319 bool _applies_jet_by_jet;
320 bool _takes_reference;
328 class SW_And:
public SW_BinaryOperator {
331 SW_And(
const Selector & s1,
const Selector & s2) : SW_BinaryOperator(s1,s2){}
334 virtual SelectorWorker* copy(){
return new SW_And(*
this);}
338 virtual bool pass(
const PseudoJet & jet)
const {
340 if (!applies_jet_by_jet())
341 throw Error(
"Cannot apply this selector worker to an individual jet");
343 return _s1.pass(jet) && _s2.pass(jet);
347 virtual void terminator(vector<const PseudoJet *> & jets)
const {
350 if (applies_jet_by_jet()){
351 SelectorWorker::terminator(jets);
356 vector<const PseudoJet *> s1_jets = jets;
357 _s1.worker()->terminator(s1_jets);
360 _s2.worker()->terminator(jets);
363 for (
unsigned int i=0; i<jets.size(); i++){
364 if (! s1_jets[i]) jets[i] = NULL;
369 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const {
370 double s1min, s1max, s2min, s2max;
371 _s1.get_rapidity_extent(s1min, s1max);
372 _s2.get_rapidity_extent(s2min, s2max);
373 rapmax = min(s1max, s2max);
374 rapmin = max(s1min, s2min);
378 virtual string description()
const {
380 ostr <<
"(" << _s1.description() <<
" && " << _s2.description() <<
")";
395 class SW_Or:
public SW_BinaryOperator {
398 SW_Or(
const Selector & s1,
const Selector & s2) : SW_BinaryOperator(s1,s2) {}
401 virtual SelectorWorker* copy(){
return new SW_Or(*
this);}
405 virtual bool pass(
const PseudoJet & jet)
const {
407 if (!applies_jet_by_jet())
408 throw Error(
"Cannot apply this selector worker to an individual jet");
410 return _s1.pass(jet) || _s2.pass(jet);
414 virtual bool applies_jet_by_jet()
const {
418 return _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
422 virtual void terminator(vector<const PseudoJet *> & jets)
const {
425 if (applies_jet_by_jet()){
426 SelectorWorker::terminator(jets);
431 vector<const PseudoJet *> s1_jets = jets;
432 _s1.worker()->terminator(s1_jets);
435 _s2.worker()->terminator(jets);
439 for (
unsigned int i=0; i<jets.size(); i++){
440 if (s1_jets[i]) jets[i] = s1_jets[i];
445 virtual string description()
const {
447 ostr <<
"(" << _s1.description() <<
" || " << _s2.description() <<
")";
452 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const {
453 double s1min, s1max, s2min, s2max;
454 _s1.get_rapidity_extent(s1min, s1max);
455 _s2.get_rapidity_extent(s2min, s2max);
456 rapmax = max(s1max, s2max);
457 rapmin = min(s1min, s2min);
469 class SW_Mult:
public SW_And {
472 SW_Mult(
const Selector & s1,
const Selector & s2) : SW_And(s1,s2) {}
475 virtual SelectorWorker* copy(){
return new SW_Mult(*
this);}
478 virtual void terminator(vector<const PseudoJet *> & jets)
const {
481 if (applies_jet_by_jet()){
482 SelectorWorker::terminator(jets);
487 _s2.worker()->terminator(jets);
490 _s1.worker()->terminator(jets);
494 virtual string description()
const {
496 ostr <<
"(" << _s1.description() <<
" * " << _s2.description() <<
")";
504 return Selector(
new SW_Mult(s1,s2));
519 QuantityBase(
double q) : _q(q){}
520 virtual ~QuantityBase(){}
521 virtual double operator()(
const PseudoJet & jet )
const =0;
522 virtual string description()
const =0;
523 virtual bool is_geometric()
const {
return false;}
524 virtual double comparison_value()
const {
return _q;}
525 virtual double description_value()
const {
return comparison_value();}
531 class QuantitySquareBase :
public QuantityBase{
533 QuantitySquareBase(
double sqrtq) : QuantityBase(sqrtq*sqrtq), _sqrtq(sqrtq){}
534 virtual double description_value()
const {
return _sqrtq;}
540 template<
typename QuantityType>
541 class SW_QuantityMin :
public SelectorWorker{
544 SW_QuantityMin(
double qmin) : _qmin(qmin) {}
547 virtual bool pass(
const PseudoJet & jet)
const {
return _qmin(jet) >= _qmin.comparison_value();}
550 virtual string description()
const {
552 ostr << _qmin.description() <<
" >= " << _qmin.description_value();
556 virtual bool is_geometric()
const {
return _qmin.is_geometric();}
564 template<
typename QuantityType>
565 class SW_QuantityMax :
public SelectorWorker {
568 SW_QuantityMax(
double qmax) : _qmax(qmax) {}
571 virtual bool pass(
const PseudoJet & jet)
const {
return _qmax(jet) <= _qmax.comparison_value();}
574 virtual string description()
const {
576 ostr << _qmax.description() <<
" <= " << _qmax.description_value();
580 virtual bool is_geometric()
const {
return _qmax.is_geometric();}
588 template<
typename QuantityType>
589 class SW_QuantityRange :
public SelectorWorker {
592 SW_QuantityRange(
double qmin,
double qmax) : _qmin(qmin), _qmax(qmax) {}
595 virtual bool pass(
const PseudoJet & jet)
const {
596 double q = _qmin(jet);
597 return (q >= _qmin.comparison_value()) && (q <= _qmax.comparison_value());
601 virtual string description()
const {
603 ostr << _qmin.description_value() <<
" <= " << _qmin.description() <<
" <= " << _qmax.description_value();
607 virtual bool is_geometric()
const {
return _qmin.is_geometric();}
617 class QuantityPt2 :
public QuantitySquareBase{
619 QuantityPt2(
double pt) : QuantitySquareBase(pt){}
620 virtual double operator()(
const PseudoJet & jet )
const {
return jet.perp2();}
621 virtual string description()
const {
return "pt";}
626 return Selector(
new SW_QuantityMin<QuantityPt2>(ptmin));
631 return Selector(
new SW_QuantityMax<QuantityPt2>(ptmax));
636 return Selector(
new SW_QuantityRange<QuantityPt2>(ptmin, ptmax));
642 class QuantityEt2 :
public QuantitySquareBase{
644 QuantityEt2(
double Et) : QuantitySquareBase(Et){}
645 virtual double operator()(
const PseudoJet & jet )
const {
return jet.Et2();}
646 virtual string description()
const {
return "Et";}
651 return Selector(
new SW_QuantityMin<QuantityEt2>(Etmin));
656 return Selector(
new SW_QuantityMax<QuantityEt2>(Etmax));
661 return Selector(
new SW_QuantityRange<QuantityEt2>(Etmin, Etmax));
667 class QuantityE :
public QuantityBase{
669 QuantityE(
double E) : QuantityBase(E){}
670 virtual double operator()(
const PseudoJet & jet )
const {
return jet.E();}
671 virtual string description()
const {
return "E";}
676 return Selector(
new SW_QuantityMin<QuantityE>(Emin));
681 return Selector(
new SW_QuantityMax<QuantityE>(Emax));
686 return Selector(
new SW_QuantityRange<QuantityE>(Emin, Emax));
692 class QuantityM2 :
public QuantitySquareBase{
694 QuantityM2(
double m) : QuantitySquareBase(m){}
695 virtual double operator()(
const PseudoJet & jet )
const {
return jet.m2();}
696 virtual string description()
const {
return "mass";}
701 return Selector(
new SW_QuantityMin<QuantityM2>(mmin));
706 return Selector(
new SW_QuantityMax<QuantityM2>(mmax));
711 return Selector(
new SW_QuantityRange<QuantityM2>(mmin, mmax));
718 class QuantityRap :
public QuantityBase{
720 QuantityRap(
double rap) : QuantityBase(rap){}
721 virtual double operator()(
const PseudoJet & jet )
const {
return jet.rap();}
722 virtual string description()
const {
return "rap";}
723 virtual bool is_geometric()
const {
return true;}
728 class SW_RapMin :
public SW_QuantityMin<QuantityRap>{
730 SW_RapMin(
double rapmin) : SW_QuantityMin<QuantityRap>(rapmin){}
731 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
732 rapmax = std::numeric_limits<double>::max();
733 rapmin = _qmin.comparison_value();
738 class SW_RapMax :
public SW_QuantityMax<QuantityRap>{
740 SW_RapMax(
double rapmax) : SW_QuantityMax<QuantityRap>(rapmax){}
741 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
742 rapmax = _qmax.comparison_value();
743 rapmin = -std::numeric_limits<double>::max();
748 class SW_RapRange :
public SW_QuantityRange<QuantityRap>{
750 SW_RapRange(
double rapmin,
double rapmax) : SW_QuantityRange<QuantityRap>(rapmin, rapmax){
751 assert(rapmin<=rapmax);
753 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
754 rapmax = _qmax.comparison_value();
755 rapmin = _qmin.comparison_value();
757 virtual bool has_known_area()
const {
return true;}
758 virtual double known_area()
const {
759 return twopi * (_qmax.comparison_value()-_qmin.comparison_value());
765 return Selector(
new SW_RapMin(rapmin));
770 return Selector(
new SW_RapMax(rapmax));
775 return Selector(
new SW_RapRange(rapmin, rapmax));
781 class QuantityAbsRap :
public QuantityBase{
783 QuantityAbsRap(
double absrap) : QuantityBase(absrap){}
784 virtual double operator()(
const PseudoJet & jet )
const {
return abs(jet.rap());}
785 virtual string description()
const {
return "|rap|";}
786 virtual bool is_geometric()
const {
return true;}
791 class SW_AbsRapMax :
public SW_QuantityMax<QuantityAbsRap>{
793 SW_AbsRapMax(
double absrapmax) : SW_QuantityMax<QuantityAbsRap>(absrapmax){}
794 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
795 rapmax = _qmax.comparison_value();
796 rapmin = -_qmax.comparison_value();
798 virtual bool has_known_area()
const {
return true;}
799 virtual double known_area()
const {
800 return twopi * 2 * _qmax.comparison_value();
805 class SW_AbsRapRange :
public SW_QuantityRange<QuantityAbsRap>{
807 SW_AbsRapRange(
double absrapmin,
double absrapmax) : SW_QuantityRange<QuantityAbsRap>(absrapmin, absrapmax){}
808 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
809 rapmax = _qmax.comparison_value();
810 rapmin = -_qmax.comparison_value();
812 virtual bool has_known_area()
const {
return true;}
813 virtual double known_area()
const {
814 return twopi * 2 * (_qmax.comparison_value()-max(_qmin.comparison_value(),0.0));
820 return Selector(
new SW_QuantityMin<QuantityAbsRap>(absrapmin));
825 return Selector(
new SW_AbsRapMax(absrapmax));
830 return Selector(
new SW_AbsRapRange(rapmin, rapmax));
836 class QuantityEta :
public QuantityBase{
838 QuantityEta(
double eta) : QuantityBase(eta){}
839 virtual double operator()(
const PseudoJet & jet )
const {
return jet.eta();}
840 virtual string description()
const {
return "eta";}
846 return Selector(
new SW_QuantityMin<QuantityEta>(etamin));
851 return Selector(
new SW_QuantityMax<QuantityEta>(etamax));
856 return Selector(
new SW_QuantityRange<QuantityEta>(etamin, etamax));
862 class QuantityAbsEta :
public QuantityBase{
864 QuantityAbsEta(
double abseta) : QuantityBase(abseta){}
865 virtual double operator()(
const PseudoJet & jet )
const {
return abs(jet.eta());}
866 virtual string description()
const {
return "|eta|";}
867 virtual bool is_geometric()
const {
return true;}
872 return Selector(
new SW_QuantityMin<QuantityAbsEta>(absetamin));
877 return Selector(
new SW_QuantityMax<QuantityAbsEta>(absetamax));
882 return Selector(
new SW_QuantityRange<QuantityAbsEta>(absetamin, absetamax));
892 class SW_PhiRange :
public SelectorWorker {
895 SW_PhiRange(
double phimin,
double phimax) : _phimin(phimin), _phimax(phimax){
896 assert(_phimin<_phimax);
897 assert(_phimin>-twopi);
898 assert(_phimax<2*twopi);
900 _phispan = _phimax - _phimin;
904 virtual bool pass(
const PseudoJet & jet)
const {
905 double dphi=jet.phi()-_phimin;
906 if (dphi >= twopi) dphi -= twopi;
907 if (dphi < 0) dphi += twopi;
908 return (dphi <= _phispan);
912 virtual string description()
const {
914 ostr << _phimin <<
" <= phi <= " << _phimax;
918 virtual bool is_geometric()
const {
return true;}
929 return Selector(
new SW_PhiRange(phimin, phimax));
934 class SW_RapPhiRange :
public SW_And{
936 SW_RapPhiRange(
double rapmin,
double rapmax,
double phimin,
double phimax)
938 _known_area = ((phimax-phimin > twopi) ? twopi : phimax-phimin) * (rapmax-rapmin);
942 virtual double known_area()
const{
951 return Selector(
new SW_RapPhiRange(rapmin, rapmax, phimin, phimax));
957 class SW_NHardest :
public SelectorWorker {
960 SW_NHardest(
unsigned int n) : _n(n) {};
965 virtual bool pass(
const PseudoJet &)
const {
966 if (!applies_jet_by_jet())
967 throw Error(
"Cannot apply this selector worker to an individual jet");
973 virtual void terminator(vector<const PseudoJet *> & jets)
const {
975 if (jets.size() < _n)
return;
982 vector<double> minus_pt2(jets.size());
983 vector<unsigned int> indices(jets.size());
985 for (
unsigned int i=0; i<jets.size(); i++){
991 minus_pt2[i] = jets[i] ? -jets[i]->perp2() : 0.0;
994 IndexedSortHelper sort_helper(& minus_pt2);
996 partial_sort(indices.begin(), indices.begin()+_n, indices.end(), sort_helper);
998 for (
unsigned int i=_n; i<jets.size(); i++)
999 jets[indices[i]] = NULL;
1003 virtual bool applies_jet_by_jet()
const {
return false;}
1006 virtual string description()
const {
1008 ostr << _n <<
" hardest";
1019 return Selector(
new SW_NHardest(n));
1030 class SW_WithReference :
public SelectorWorker{
1033 SW_WithReference() : _is_initialised(false){};
1036 virtual bool takes_reference()
const {
return true;}
1039 virtual void set_reference(
const PseudoJet ¢re){
1040 _is_initialised =
true;
1041 _reference = centre;
1045 PseudoJet _reference;
1046 bool _is_initialised;
1051 class SW_Circle :
public SW_WithReference {
1053 SW_Circle(
const double &radius) : _radius2(radius*radius) {}
1056 virtual SelectorWorker* copy(){
return new SW_Circle(*
this);}
1060 virtual bool pass(
const PseudoJet & jet)
const {
1062 if (! _is_initialised)
1063 throw Error(
"To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1065 return jet.squared_distance(_reference) <= _radius2;
1069 virtual string description()
const {
1071 ostr <<
"distance from the centre <= " << sqrt(_radius2);
1076 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1078 if (! _is_initialised)
1079 throw Error(
"To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1081 rapmax = _reference.rap()+sqrt(_radius2);
1082 rapmin = _reference.rap()-sqrt(_radius2);
1085 virtual bool is_geometric()
const {
return true;}
1086 virtual bool has_finite_area()
const {
return true;}
1087 virtual bool has_known_area()
const {
return true;}
1088 virtual double known_area()
const {
1089 return pi * _radius2;
1099 return Selector(
new SW_Circle(radius));
1106 class SW_Doughnut :
public SW_WithReference {
1108 SW_Doughnut(
const double &radius_in,
const double &radius_out)
1109 : _radius_in2(radius_in*radius_in), _radius_out2(radius_out*radius_out) {}
1112 virtual SelectorWorker* copy(){
return new SW_Doughnut(*
this);}
1116 virtual bool pass(
const PseudoJet & jet)
const {
1118 if (! _is_initialised)
1119 throw Error(
"To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1121 double distance2 = jet.squared_distance(_reference);
1123 return (distance2 <= _radius_out2) && (distance2 >= _radius_in2);
1127 virtual string description()
const {
1129 ostr << sqrt(_radius_in2) <<
" <= distance from the centre <= " << sqrt(_radius_out2);
1134 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1136 if (! _is_initialised)
1137 throw Error(
"To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1139 rapmax = _reference.rap()+sqrt(_radius_out2);
1140 rapmin = _reference.rap()-sqrt(_radius_out2);
1143 virtual bool is_geometric()
const {
return true;}
1144 virtual bool has_finite_area()
const {
return true;}
1145 virtual bool has_known_area()
const {
return true;}
1146 virtual double known_area()
const {
1147 return pi * (_radius_out2-_radius_in2);
1151 double _radius_in2, _radius_out2;
1158 return Selector(
new SW_Doughnut(radius_in, radius_out));
1164 class SW_Strip :
public SW_WithReference {
1166 SW_Strip(
const double &delta) : _delta(delta) {}
1169 virtual SelectorWorker* copy(){
return new SW_Strip(*
this);}
1173 virtual bool pass(
const PseudoJet & jet)
const {
1175 if (! _is_initialised)
1176 throw Error(
"To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1178 return abs(jet.rap()-_reference.rap()) <= _delta;
1182 virtual string description()
const {
1184 ostr <<
"|rap - rap_reference| <= " << _delta;
1189 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1191 if (! _is_initialised)
1192 throw Error(
"To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1194 rapmax = _reference.rap()+_delta;
1195 rapmin = _reference.rap()-_delta;
1198 virtual bool is_geometric()
const {
return true;}
1199 virtual bool has_finite_area()
const {
return true;}
1200 virtual bool has_known_area()
const {
return true;}
1201 virtual double known_area()
const {
1202 return twopi * 2 * _delta;
1212 return Selector(
new SW_Strip(half_width));
1220 class SW_Rectangle :
public SW_WithReference {
1222 SW_Rectangle(
const double &delta_rap,
const double &delta_phi)
1223 : _delta_rap(delta_rap), _delta_phi(delta_phi) {}
1226 virtual SelectorWorker* copy(){
return new SW_Rectangle(*
this);}
1230 virtual bool pass(
const PseudoJet & jet)
const {
1232 if (! _is_initialised)
1233 throw Error(
"To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1235 return (abs(jet.rap()-_reference.rap()) <= _delta_rap) && (abs(jet.delta_phi_to(_reference)) <= _delta_phi);
1239 virtual string description()
const {
1241 ostr <<
"|rap - rap_reference| <= " << _delta_rap <<
" && |phi - phi_reference| <= " << _delta_phi ;
1246 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1248 if (! _is_initialised)
1249 throw Error(
"To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1251 rapmax = _reference.rap()+_delta_rap;
1252 rapmin = _reference.rap()-_delta_rap;
1255 virtual bool is_geometric()
const {
return true;}
1256 virtual bool has_finite_area()
const {
return true;}
1257 virtual bool has_known_area()
const {
return true;}
1258 virtual double known_area()
const {
1259 return 4 * _delta_rap * _delta_phi;
1263 double _delta_rap, _delta_phi;
1269 return Selector(
new SW_Rectangle(half_rap_width, half_phi_width));
1276 class SW_PtFractionMin :
public SW_WithReference {
1279 SW_PtFractionMin(
double fraction) : _fraction2(fraction*fraction){}
1282 virtual SelectorWorker* copy(){
return new SW_PtFractionMin(*
this);}
1286 virtual bool pass(
const PseudoJet & jet)
const {
1288 if (! _is_initialised)
1289 throw Error(
"To use a SelectorPtFractionMin (or any selector that requires a reference), you first have to call set_reference(...)");
1292 return (jet.perp2() >= _fraction2*_reference.perp2());
1296 virtual string description()
const {
1298 ostr <<
"pt >= " << sqrt(_fraction2) <<
"* pt_ref";
1310 return Selector(
new SW_PtFractionMin(fraction));
1320 class SW_IsZero :
public SelectorWorker {
1326 virtual bool pass(
const PseudoJet & jet)
const {
1331 virtual string description()
const {
return "zero";}
1343 class SW_IsPureGhost :
public SelectorWorker {
1349 virtual bool pass(
const PseudoJet & jet)
const {
1351 if (!jet.has_area())
return false;
1354 return jet.is_pure_ghost();
1358 virtual string description()
const {
return "pure ghost";}
1364 return Selector(
new SW_IsPureGhost());
1378 class SW_RangeDefinition :
public SelectorWorker{
1381 SW_RangeDefinition(
const RangeDefinition &range) : _range(&range){}
1384 virtual bool pass(
const PseudoJet & jet)
const {
1385 return _range->is_in_range(jet);
1389 virtual string description()
const {
1390 return _range->description();
1394 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1395 _range->get_rap_limits(rapmin, rapmax);
1399 virtual bool is_geometric()
const {
return true;}
1402 virtual bool has_known_area()
const {
return true;}
1405 virtual double known_area()
const{
1406 return _range->area();
1410 const RangeDefinition *_range;
1420 _worker.reset(
new SW_RangeDefinition(range));
1430 _worker.reset(
new SW_And(*
this, b));
1437 _worker.reset(
new SW_Or(*
this, b));
1441 FASTJET_END_NAMESPACE