32 #include "fastjet/Selector.hh"
34 #include "fastjet/GhostedAreaSpec.hh"
39 FASTJET_BEGIN_NAMESPACE
46 std::vector<PseudoJet> Selector::operator()(
const std::vector<PseudoJet> & jets)
const {
47 std::vector<PseudoJet> result;
52 for (std::vector<PseudoJet>::const_iterator jet = jets.begin();
53 jet != jets.end(); jet++) {
54 if (worker_local->
pass(*jet)) result.push_back(*jet);
59 std::vector<const PseudoJet *> jetptrs(jets.size());
60 for (
unsigned i = 0; i < jets.size(); i++) {
61 jetptrs[i] = & jets[i];
64 for (
unsigned i = 0; i < jetptrs.size(); i++) {
65 if (jetptrs[i]) result.push_back(jets[i]);
74 unsigned int Selector::count(
const std::vector<PseudoJet> & jets)
const {
80 for (
unsigned i = 0; i < jets.size(); i++) {
81 if (worker_local->
pass(jets[i])) n++;
84 std::vector<const PseudoJet *> jetptrs(jets.size());
85 for (
unsigned i = 0; i < jets.size(); i++) {
86 jetptrs[i] = & jets[i];
89 for (
unsigned i = 0; i < jetptrs.size(); i++) {
101 void Selector::sift(
const std::vector<PseudoJet> & jets,
102 std::vector<PseudoJet> & jets_that_pass,
103 std::vector<PseudoJet> & jets_that_fail
107 jets_that_pass.clear();
108 jets_that_fail.clear();
112 for (
unsigned i = 0; i < jets.size(); i++) {
113 if (worker_local->
pass(jets[i])) {
114 jets_that_pass.push_back(jets[i]);
116 jets_that_fail.push_back(jets[i]);
120 std::vector<const PseudoJet *> jetptrs(jets.size());
121 for (
unsigned i = 0; i < jets.size(); i++) {
122 jetptrs[i] = & jets[i];
125 for (
unsigned i = 0; i < jetptrs.size(); i++) {
127 jets_that_pass.push_back(jets[i]);
129 jets_that_fail.push_back(jets[i]);
137 double Selector::area()
const{
138 return area(gas::def_ghost_area);
142 double Selector::area(
double ghost_area)
const{
146 if (_worker->has_known_area())
return _worker->known_area();
149 double rapmin, rapmax;
150 get_rapidity_extent(rapmin, rapmax);
152 std::vector<PseudoJet> ghosts;
156 return ghost_spec.ghost_area() * ((*this)(ghosts)).size();
165 bool SelectorWorker::has_finite_area()
const {
166 if (! is_geometric())
return false;
167 double rapmin, rapmax;
168 get_rapidity_extent(rapmin, rapmax);
169 return (rapmax != std::numeric_limits<double>::infinity())
170 && (-rapmin != std::numeric_limits<double>::infinity());
187 virtual bool pass(
const PseudoJet &)
const {
193 virtual void terminator(vector<const PseudoJet *> &)
const {
199 virtual string description()
const {
return "Identity";}
202 virtual bool is_geometric()
const {
return true;}
207 Selector SelectorIdentity() {
208 return Selector(
new SW_Identity);
218 class SW_Not :
public SelectorWorker {
221 SW_Not(
const Selector & s) : _s(s) {}
224 virtual SelectorWorker* copy(){
return new SW_Not(*
this);}
228 virtual bool pass(
const PseudoJet & jet)
const {
230 if (!applies_jet_by_jet())
231 throw Error(
"Cannot apply this selector worker to an individual jet");
233 return ! _s.pass(jet);
237 virtual bool applies_jet_by_jet()
const {
return _s.applies_jet_by_jet();}
240 virtual void terminator(vector<const PseudoJet *> & jets)
const {
243 if (applies_jet_by_jet()){
244 SelectorWorker::terminator(jets);
249 vector<const PseudoJet *> s_jets = jets;
250 _s.worker()->terminator(s_jets);
254 for (
unsigned int i=0; i<s_jets.size(); i++){
255 if (s_jets[i]) jets[i] = NULL;
260 virtual string description()
const {
262 ostr <<
"!(" << _s.description() <<
")";
267 virtual bool is_geometric()
const {
return _s.is_geometric();}
270 virtual bool takes_reference()
const {
return _s.takes_reference();}
273 virtual void set_reference(
const PseudoJet &ref) { _s.set_reference(ref);}
288 class SW_BinaryOperator:
public SelectorWorker {
291 SW_BinaryOperator(
const Selector & s1,
const Selector & s2) : _s1(s1), _s2(s2) {
295 _applies_jet_by_jet = _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
298 _takes_reference = _s1.takes_reference() || _s2.takes_reference();
301 _is_geometric = _s1.is_geometric() && _s2.is_geometric();
305 virtual bool applies_jet_by_jet()
const {
return _applies_jet_by_jet;}
308 virtual bool takes_reference()
const{
309 return _takes_reference;
313 virtual void set_reference(
const PseudoJet ¢re){
314 _s1.set_reference(centre);
315 _s2.set_reference(centre);
319 virtual bool is_geometric()
const {
return _is_geometric;}
323 bool _applies_jet_by_jet;
324 bool _takes_reference;
332 class SW_And:
public SW_BinaryOperator {
335 SW_And(
const Selector & s1,
const Selector & s2) : SW_BinaryOperator(s1,s2){}
338 virtual SelectorWorker* copy(){
return new SW_And(*
this);}
342 virtual bool pass(
const PseudoJet & jet)
const {
344 if (!applies_jet_by_jet())
345 throw Error(
"Cannot apply this selector worker to an individual jet");
347 return _s1.pass(jet) && _s2.pass(jet);
351 virtual void terminator(vector<const PseudoJet *> & jets)
const {
354 if (applies_jet_by_jet()){
355 SelectorWorker::terminator(jets);
360 vector<const PseudoJet *> s1_jets = jets;
361 _s1.worker()->terminator(s1_jets);
364 _s2.worker()->terminator(jets);
367 for (
unsigned int i=0; i<jets.size(); i++){
368 if (! s1_jets[i]) jets[i] = NULL;
373 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const {
374 double s1min, s1max, s2min, s2max;
375 _s1.get_rapidity_extent(s1min, s1max);
376 _s2.get_rapidity_extent(s2min, s2max);
377 rapmax = min(s1max, s2max);
378 rapmin = max(s1min, s2min);
382 virtual string description()
const {
384 ostr <<
"(" << _s1.description() <<
" && " << _s2.description() <<
")";
399 class SW_Or:
public SW_BinaryOperator {
402 SW_Or(
const Selector & s1,
const Selector & s2) : SW_BinaryOperator(s1,s2) {}
405 virtual SelectorWorker* copy(){
return new SW_Or(*
this);}
409 virtual bool pass(
const PseudoJet & jet)
const {
411 if (!applies_jet_by_jet())
412 throw Error(
"Cannot apply this selector worker to an individual jet");
414 return _s1.pass(jet) || _s2.pass(jet);
418 virtual bool applies_jet_by_jet()
const {
422 return _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
426 virtual void terminator(vector<const PseudoJet *> & jets)
const {
429 if (applies_jet_by_jet()){
430 SelectorWorker::terminator(jets);
435 vector<const PseudoJet *> s1_jets = jets;
436 _s1.worker()->terminator(s1_jets);
439 _s2.worker()->terminator(jets);
443 for (
unsigned int i=0; i<jets.size(); i++){
444 if (s1_jets[i]) jets[i] = s1_jets[i];
449 virtual string description()
const {
451 ostr <<
"(" << _s1.description() <<
" || " << _s2.description() <<
")";
456 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const {
457 double s1min, s1max, s2min, s2max;
458 _s1.get_rapidity_extent(s1min, s1max);
459 _s2.get_rapidity_extent(s2min, s2max);
460 rapmax = max(s1max, s2max);
461 rapmin = min(s1min, s2min);
473 class SW_Mult:
public SW_And {
476 SW_Mult(
const Selector & s1,
const Selector & s2) : SW_And(s1,s2) {}
479 virtual SelectorWorker* copy(){
return new SW_Mult(*
this);}
482 virtual void terminator(vector<const PseudoJet *> & jets)
const {
485 if (applies_jet_by_jet()){
486 SelectorWorker::terminator(jets);
491 _s2.worker()->terminator(jets);
494 _s1.worker()->terminator(jets);
498 virtual string description()
const {
500 ostr <<
"(" << _s1.description() <<
" * " << _s2.description() <<
")";
508 return Selector(
new SW_Mult(s1,s2));
523 QuantityBase(
double q) : _q(q){}
524 virtual ~QuantityBase(){}
525 virtual double operator()(
const PseudoJet & jet )
const =0;
526 virtual string description()
const =0;
527 virtual bool is_geometric()
const {
return false;}
528 virtual double comparison_value()
const {
return _q;}
529 virtual double description_value()
const {
return comparison_value();}
535 class QuantitySquareBase :
public QuantityBase{
537 QuantitySquareBase(
double sqrtq) : QuantityBase(sqrtq*sqrtq), _sqrtq(sqrtq){}
538 virtual double description_value()
const {
return _sqrtq;}
544 template<
typename QuantityType>
545 class SW_QuantityMin :
public SelectorWorker{
548 SW_QuantityMin(
double qmin) : _qmin(qmin) {}
551 virtual bool pass(
const PseudoJet & jet)
const {
return _qmin(jet) >= _qmin.comparison_value();}
554 virtual string description()
const {
556 ostr << _qmin.description() <<
" >= " << _qmin.description_value();
560 virtual bool is_geometric()
const {
return _qmin.is_geometric();}
568 template<
typename QuantityType>
569 class SW_QuantityMax :
public SelectorWorker {
572 SW_QuantityMax(
double qmax) : _qmax(qmax) {}
575 virtual bool pass(
const PseudoJet & jet)
const {
return _qmax(jet) <= _qmax.comparison_value();}
578 virtual string description()
const {
580 ostr << _qmax.description() <<
" <= " << _qmax.description_value();
584 virtual bool is_geometric()
const {
return _qmax.is_geometric();}
592 template<
typename QuantityType>
593 class SW_QuantityRange :
public SelectorWorker {
596 SW_QuantityRange(
double qmin,
double qmax) : _qmin(qmin), _qmax(qmax) {}
599 virtual bool pass(
const PseudoJet & jet)
const {
600 double q = _qmin(jet);
601 return (q >= _qmin.comparison_value()) && (q <= _qmax.comparison_value());
605 virtual string description()
const {
607 ostr << _qmin.description_value() <<
" <= " << _qmin.description() <<
" <= " << _qmax.description_value();
611 virtual bool is_geometric()
const {
return _qmin.is_geometric();}
621 class QuantityPt2 :
public QuantitySquareBase{
623 QuantityPt2(
double pt) : QuantitySquareBase(pt){}
624 virtual double operator()(
const PseudoJet & jet )
const {
return jet.perp2();}
625 virtual string description()
const {
return "pt";}
630 return Selector(
new SW_QuantityMin<QuantityPt2>(ptmin));
635 return Selector(
new SW_QuantityMax<QuantityPt2>(ptmax));
640 return Selector(
new SW_QuantityRange<QuantityPt2>(ptmin, ptmax));
646 class QuantityEt2 :
public QuantitySquareBase{
648 QuantityEt2(
double Et) : QuantitySquareBase(Et){}
649 virtual double operator()(
const PseudoJet & jet )
const {
return jet.Et2();}
650 virtual string description()
const {
return "Et";}
655 return Selector(
new SW_QuantityMin<QuantityEt2>(Etmin));
660 return Selector(
new SW_QuantityMax<QuantityEt2>(Etmax));
665 return Selector(
new SW_QuantityRange<QuantityEt2>(Etmin, Etmax));
671 class QuantityE :
public QuantityBase{
673 QuantityE(
double E) : QuantityBase(E){}
674 virtual double operator()(
const PseudoJet & jet )
const {
return jet.E();}
675 virtual string description()
const {
return "E";}
680 return Selector(
new SW_QuantityMin<QuantityE>(Emin));
685 return Selector(
new SW_QuantityMax<QuantityE>(Emax));
690 return Selector(
new SW_QuantityRange<QuantityE>(Emin, Emax));
696 class QuantityM2 :
public QuantitySquareBase{
698 QuantityM2(
double m) : QuantitySquareBase(m){}
699 virtual double operator()(
const PseudoJet & jet )
const {
return jet.m2();}
700 virtual string description()
const {
return "mass";}
705 return Selector(
new SW_QuantityMin<QuantityM2>(mmin));
710 return Selector(
new SW_QuantityMax<QuantityM2>(mmax));
715 return Selector(
new SW_QuantityRange<QuantityM2>(mmin, mmax));
722 class QuantityRap :
public QuantityBase{
724 QuantityRap(
double rap) : QuantityBase(rap){}
725 virtual double operator()(
const PseudoJet & jet )
const {
return jet.rap();}
726 virtual string description()
const {
return "rap";}
727 virtual bool is_geometric()
const {
return true;}
732 class SW_RapMin :
public SW_QuantityMin<QuantityRap>{
734 SW_RapMin(
double rapmin) : SW_QuantityMin<QuantityRap>(rapmin){}
735 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
736 rapmax = std::numeric_limits<double>::max();
737 rapmin = _qmin.comparison_value();
742 class SW_RapMax :
public SW_QuantityMax<QuantityRap>{
744 SW_RapMax(
double rapmax) : SW_QuantityMax<QuantityRap>(rapmax){}
745 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
746 rapmax = _qmax.comparison_value();
747 rapmin = -std::numeric_limits<double>::max();
752 class SW_RapRange :
public SW_QuantityRange<QuantityRap>{
754 SW_RapRange(
double rapmin,
double rapmax) : SW_QuantityRange<QuantityRap>(rapmin, rapmax){
755 assert(rapmin<=rapmax);
757 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
758 rapmax = _qmax.comparison_value();
759 rapmin = _qmin.comparison_value();
761 virtual bool has_known_area()
const {
return true;}
762 virtual double known_area()
const {
763 return twopi * (_qmax.comparison_value()-_qmin.comparison_value());
769 return Selector(
new SW_RapMin(rapmin));
774 return Selector(
new SW_RapMax(rapmax));
779 return Selector(
new SW_RapRange(rapmin, rapmax));
785 class QuantityAbsRap :
public QuantityBase{
787 QuantityAbsRap(
double absrap) : QuantityBase(absrap){}
788 virtual double operator()(
const PseudoJet & jet )
const {
return abs(jet.rap());}
789 virtual string description()
const {
return "|rap|";}
790 virtual bool is_geometric()
const {
return true;}
795 class SW_AbsRapMax :
public SW_QuantityMax<QuantityAbsRap>{
797 SW_AbsRapMax(
double absrapmax) : SW_QuantityMax<QuantityAbsRap>(absrapmax){}
798 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
799 rapmax = _qmax.comparison_value();
800 rapmin = -_qmax.comparison_value();
802 virtual bool has_known_area()
const {
return true;}
803 virtual double known_area()
const {
804 return twopi * 2 * _qmax.comparison_value();
809 class SW_AbsRapRange :
public SW_QuantityRange<QuantityAbsRap>{
811 SW_AbsRapRange(
double absrapmin,
double absrapmax) : SW_QuantityRange<QuantityAbsRap>(absrapmin, absrapmax){}
812 virtual void get_rapidity_extent(
double &rapmin,
double & rapmax)
const{
813 rapmax = _qmax.comparison_value();
814 rapmin = -_qmax.comparison_value();
816 virtual bool has_known_area()
const {
return true;}
817 virtual double known_area()
const {
818 return twopi * 2 * (_qmax.comparison_value()-max(_qmin.comparison_value(),0.0));
824 return Selector(
new SW_QuantityMin<QuantityAbsRap>(absrapmin));
829 return Selector(
new SW_AbsRapMax(absrapmax));
834 return Selector(
new SW_AbsRapRange(rapmin, rapmax));
840 class QuantityEta :
public QuantityBase{
842 QuantityEta(
double eta) : QuantityBase(eta){}
843 virtual double operator()(
const PseudoJet & jet )
const {
return jet.eta();}
844 virtual string description()
const {
return "eta";}
850 return Selector(
new SW_QuantityMin<QuantityEta>(etamin));
855 return Selector(
new SW_QuantityMax<QuantityEta>(etamax));
860 return Selector(
new SW_QuantityRange<QuantityEta>(etamin, etamax));
866 class QuantityAbsEta :
public QuantityBase{
868 QuantityAbsEta(
double abseta) : QuantityBase(abseta){}
869 virtual double operator()(
const PseudoJet & jet )
const {
return abs(jet.eta());}
870 virtual string description()
const {
return "|eta|";}
871 virtual bool is_geometric()
const {
return true;}
876 return Selector(
new SW_QuantityMin<QuantityAbsEta>(absetamin));
881 return Selector(
new SW_QuantityMax<QuantityAbsEta>(absetamax));
886 return Selector(
new SW_QuantityRange<QuantityAbsEta>(absetamin, absetamax));
896 class SW_PhiRange :
public SelectorWorker {
899 SW_PhiRange(
double phimin,
double phimax) : _phimin(phimin), _phimax(phimax){
900 assert(_phimin<_phimax);
901 assert(_phimin>-twopi);
902 assert(_phimax<2*twopi);
904 _phispan = _phimax - _phimin;
908 virtual bool pass(
const PseudoJet & jet)
const {
909 double dphi=jet.phi()-_phimin;
910 if (dphi >= twopi) dphi -= twopi;
911 if (dphi < 0) dphi += twopi;
912 return (dphi <= _phispan);
916 virtual string description()
const {
918 ostr << _phimin <<
" <= phi <= " << _phimax;
922 virtual bool is_geometric()
const {
return true;}
933 return Selector(
new SW_PhiRange(phimin, phimax));
938 class SW_RapPhiRange :
public SW_And{
940 SW_RapPhiRange(
double rapmin,
double rapmax,
double phimin,
double phimax)
942 _known_area = ((phimax-phimin > twopi) ? twopi : phimax-phimin) * (rapmax-rapmin);
946 virtual double known_area()
const{
955 return Selector(
new SW_RapPhiRange(rapmin, rapmax, phimin, phimax));
961 class SW_NHardest :
public SelectorWorker {
964 SW_NHardest(
unsigned int n) : _n(n) {};
969 virtual bool pass(
const PseudoJet &)
const {
970 if (!applies_jet_by_jet())
971 throw Error(
"Cannot apply this selector worker to an individual jet");
977 virtual void terminator(vector<const PseudoJet *> & jets)
const {
979 if (jets.size() < _n)
return;
986 vector<double> minus_pt2(jets.size());
987 vector<unsigned int> indices(jets.size());
989 for (
unsigned int i=0; i<jets.size(); i++){
995 minus_pt2[i] = jets[i] ? -jets[i]->perp2() : 0.0;
998 IndexedSortHelper sort_helper(& minus_pt2);
1000 partial_sort(indices.begin(), indices.begin()+_n, indices.end(), sort_helper);
1002 for (
unsigned int i=_n; i<jets.size(); i++)
1003 jets[indices[i]] = NULL;
1007 virtual bool applies_jet_by_jet()
const {
return false;}
1010 virtual string description()
const {
1012 ostr << _n <<
" hardest";
1023 return Selector(
new SW_NHardest(n));
1034 class SW_WithReference :
public SelectorWorker{
1037 SW_WithReference() : _is_initialised(false){};
1040 virtual bool takes_reference()
const {
return true;}
1043 virtual void set_reference(
const PseudoJet ¢re){
1044 _is_initialised =
true;
1045 _reference = centre;
1049 PseudoJet _reference;
1050 bool _is_initialised;
1055 class SW_Circle :
public SW_WithReference {
1057 SW_Circle(
const double &radius) : _radius2(radius*radius) {}
1060 virtual SelectorWorker* copy(){
return new SW_Circle(*
this);}
1064 virtual bool pass(
const PseudoJet & jet)
const {
1066 if (! _is_initialised)
1067 throw Error(
"To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1069 return jet.squared_distance(_reference) <= _radius2;
1073 virtual string description()
const {
1075 ostr <<
"distance from the centre <= " << sqrt(_radius2);
1080 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1082 if (! _is_initialised)
1083 throw Error(
"To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1085 rapmax = _reference.rap()+sqrt(_radius2);
1086 rapmin = _reference.rap()-sqrt(_radius2);
1089 virtual bool is_geometric()
const {
return true;}
1090 virtual bool has_finite_area()
const {
return true;}
1091 virtual bool has_known_area()
const {
return true;}
1092 virtual double known_area()
const {
1093 return pi * _radius2;
1103 return Selector(
new SW_Circle(radius));
1110 class SW_Doughnut :
public SW_WithReference {
1112 SW_Doughnut(
const double &radius_in,
const double &radius_out)
1113 : _radius_in2(radius_in*radius_in), _radius_out2(radius_out*radius_out) {}
1116 virtual SelectorWorker* copy(){
return new SW_Doughnut(*
this);}
1120 virtual bool pass(
const PseudoJet & jet)
const {
1122 if (! _is_initialised)
1123 throw Error(
"To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1125 double distance2 = jet.squared_distance(_reference);
1127 return (distance2 <= _radius_out2) && (distance2 >= _radius_in2);
1131 virtual string description()
const {
1133 ostr << sqrt(_radius_in2) <<
" <= distance from the centre <= " << sqrt(_radius_out2);
1138 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1140 if (! _is_initialised)
1141 throw Error(
"To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1143 rapmax = _reference.rap()+sqrt(_radius_out2);
1144 rapmin = _reference.rap()-sqrt(_radius_out2);
1147 virtual bool is_geometric()
const {
return true;}
1148 virtual bool has_finite_area()
const {
return true;}
1149 virtual bool has_known_area()
const {
return true;}
1150 virtual double known_area()
const {
1151 return pi * (_radius_out2-_radius_in2);
1155 double _radius_in2, _radius_out2;
1162 return Selector(
new SW_Doughnut(radius_in, radius_out));
1168 class SW_Strip :
public SW_WithReference {
1170 SW_Strip(
const double &delta) : _delta(delta) {}
1173 virtual SelectorWorker* copy(){
return new SW_Strip(*
this);}
1177 virtual bool pass(
const PseudoJet & jet)
const {
1179 if (! _is_initialised)
1180 throw Error(
"To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1182 return abs(jet.rap()-_reference.rap()) <= _delta;
1186 virtual string description()
const {
1188 ostr <<
"|rap - rap_reference| <= " << _delta;
1193 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1195 if (! _is_initialised)
1196 throw Error(
"To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1198 rapmax = _reference.rap()+_delta;
1199 rapmin = _reference.rap()-_delta;
1202 virtual bool is_geometric()
const {
return true;}
1203 virtual bool has_finite_area()
const {
return true;}
1204 virtual bool has_known_area()
const {
return true;}
1205 virtual double known_area()
const {
1206 return twopi * 2 * _delta;
1216 return Selector(
new SW_Strip(half_width));
1224 class SW_Rectangle :
public SW_WithReference {
1226 SW_Rectangle(
const double &delta_rap,
const double &delta_phi)
1227 : _delta_rap(delta_rap), _delta_phi(delta_phi) {}
1230 virtual SelectorWorker* copy(){
return new SW_Rectangle(*
this);}
1234 virtual bool pass(
const PseudoJet & jet)
const {
1236 if (! _is_initialised)
1237 throw Error(
"To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1239 return (abs(jet.rap()-_reference.rap()) <= _delta_rap) && (abs(jet.delta_phi_to(_reference)) <= _delta_phi);
1243 virtual string description()
const {
1245 ostr <<
"|rap - rap_reference| <= " << _delta_rap <<
" && |phi - phi_reference| <= " << _delta_phi ;
1250 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1252 if (! _is_initialised)
1253 throw Error(
"To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1255 rapmax = _reference.rap()+_delta_rap;
1256 rapmin = _reference.rap()-_delta_rap;
1259 virtual bool is_geometric()
const {
return true;}
1260 virtual bool has_finite_area()
const {
return true;}
1261 virtual bool has_known_area()
const {
return true;}
1262 virtual double known_area()
const {
1263 return 4 * _delta_rap * _delta_phi;
1267 double _delta_rap, _delta_phi;
1273 return Selector(
new SW_Rectangle(half_rap_width, half_phi_width));
1280 class SW_PtFractionMin :
public SW_WithReference {
1283 SW_PtFractionMin(
double fraction) : _fraction2(fraction*fraction){}
1286 virtual SelectorWorker* copy(){
return new SW_PtFractionMin(*
this);}
1290 virtual bool pass(
const PseudoJet & jet)
const {
1292 if (! _is_initialised)
1293 throw Error(
"To use a SelectorPtFractionMin (or any selector that requires a reference), you first have to call set_reference(...)");
1296 return (jet.perp2() >= _fraction2*_reference.perp2());
1300 virtual string description()
const {
1302 ostr <<
"pt >= " << sqrt(_fraction2) <<
"* pt_ref";
1314 return Selector(
new SW_PtFractionMin(fraction));
1324 class SW_IsZero :
public SelectorWorker {
1330 virtual bool pass(
const PseudoJet & jet)
const {
1335 virtual string description()
const {
return "zero";}
1348 class SW_IsPureGhost :
public SelectorWorker {
1354 virtual bool pass(
const PseudoJet & jet)
const {
1356 if (!jet.has_area())
return false;
1359 return jet.is_pure_ghost();
1363 virtual string description()
const {
return "pure ghost";}
1369 return Selector(
new SW_IsPureGhost());
1382 class SW_RangeDefinition :
public SelectorWorker{
1385 SW_RangeDefinition(
const RangeDefinition &range) : _range(&range){}
1388 virtual bool pass(
const PseudoJet & jet)
const {
1389 return _range->is_in_range(jet);
1393 virtual string description()
const {
1394 return _range->description();
1398 virtual void get_rapidity_extent(
double & rapmin,
double & rapmax)
const{
1399 _range->get_rap_limits(rapmin, rapmax);
1403 virtual bool is_geometric()
const {
return true;}
1406 virtual bool has_known_area()
const {
return true;}
1409 virtual double known_area()
const{
1410 return _range->area();
1414 const RangeDefinition *_range;
1424 _worker.reset(
new SW_RangeDefinition(range));
1426 #endif // __FJCORE__
1435 _worker.reset(
new SW_And(*
this, b));
1442 _worker.reset(
new SW_Or(*
this, b));
1446 FASTJET_END_NAMESPACE