Fork me on GitHub

source: git/external/fastjet/Selector.cc@ 942a705

Last change on this file since 942a705 was b7b836a, checked in by Pavel Demin <pavel-demin@…>, 6 years ago

update FastJet library to 3.3.1 and FastJet Contrib library to 1.036

  • Property mode set to 100644
File size: 48.8 KB
Line 
1//FJSTARTHEADER
2// $Id: Selector.cc 4354 2018-04-22 07:12:37Z salam $
3//
4// Copyright (c) 2005-2018, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5//
6//----------------------------------------------------------------------
7// This file is part of FastJet.
8//
9// FastJet is free software; you can redistribute it and/or modify
10// it under the terms of the GNU General Public License as published by
11// the Free Software Foundation; either version 2 of the License, or
12// (at your option) any later version.
13//
14// The algorithms that underlie FastJet have required considerable
15// development. They are described in the original FastJet paper,
16// hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17// FastJet as part of work towards a scientific publication, please
18// quote the version you use and include a citation to the manual and
19// optionally also to hep-ph/0512210.
20//
21// FastJet is distributed in the hope that it will be useful,
22// but WITHOUT ANY WARRANTY; without even the implied warranty of
23// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24// GNU General Public License for more details.
25//
26// You should have received a copy of the GNU General Public License
27// along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28//----------------------------------------------------------------------
29//FJENDHEADER
30
31
32#include <sstream>
33#include <algorithm>
34#include "fastjet/Selector.hh"
35#ifndef __FJCORE__
36#include "fastjet/GhostedAreaSpec.hh" // for area support
37#endif // __FJCORE__
38
39using namespace std;
40
41FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
42
43//----------------------------------------------------------------------
44// implementations of some of the more complex bits of Selector
45//----------------------------------------------------------------------
46
47// implementation of the operator() acting on a vector of jets
48std::vector<PseudoJet> Selector::operator()(const std::vector<PseudoJet> & jets) const {
49 std::vector<PseudoJet> result;
50 const SelectorWorker * worker_local = validated_worker();
51 if (worker_local->applies_jet_by_jet()) {
52 //if (false) {
53 // for workers that apply jet by jet, this is more efficient
54 for (std::vector<PseudoJet>::const_iterator jet = jets.begin();
55 jet != jets.end(); jet++) {
56 if (worker_local->pass(*jet)) result.push_back(*jet);
57 }
58 } else {
59 // for workers that can only be applied to entire vectors,
60 // go through the following
61 std::vector<const PseudoJet *> jetptrs(jets.size());
62 for (unsigned i = 0; i < jets.size(); i++) {
63 jetptrs[i] = & jets[i];
64 }
65 worker_local->terminator(jetptrs);
66 for (unsigned i = 0; i < jetptrs.size(); i++) {
67 if (jetptrs[i]) result.push_back(jets[i]);
68 }
69 }
70 return result;
71}
72
73
74//----------------------------------------------------------------------
75// count the number of jets that pass the cuts
76unsigned int Selector::count(const std::vector<PseudoJet> & jets) const {
77 unsigned n = 0;
78 const SelectorWorker * worker_local = validated_worker();
79
80 // separate strategies according to whether the worker applies jet by jet
81 if (worker_local->applies_jet_by_jet()) {
82 for (unsigned i = 0; i < jets.size(); i++) {
83 if (worker_local->pass(jets[i])) n++;
84 }
85 } else {
86 std::vector<const PseudoJet *> jetptrs(jets.size());
87 for (unsigned i = 0; i < jets.size(); i++) {
88 jetptrs[i] = & jets[i];
89 }
90 worker_local->terminator(jetptrs);
91 for (unsigned i = 0; i < jetptrs.size(); i++) {
92 if (jetptrs[i]) n++;
93 }
94 }
95
96 return n;
97}
98
99//----------------------------------------------------------------------
100// sum the momenta of the jets that pass the cuts
101PseudoJet Selector::sum(const std::vector<PseudoJet> & jets) const {
102 PseudoJet this_sum(0,0,0,0);
103 const SelectorWorker * worker_local = validated_worker();
104
105 // separate strategies according to whether the worker applies jet by jet
106 if (worker_local->applies_jet_by_jet()) {
107 for (unsigned i = 0; i < jets.size(); i++) {
108 if (worker_local->pass(jets[i])) this_sum += jets[i];
109 }
110 } else {
111 std::vector<const PseudoJet *> jetptrs(jets.size());
112 for (unsigned i = 0; i < jets.size(); i++) {
113 jetptrs[i] = & jets[i];
114 }
115 worker_local->terminator(jetptrs);
116 for (unsigned i = 0; i < jetptrs.size(); i++) {
117 if (jetptrs[i]) this_sum += jets[i];
118 }
119 }
120
121 return this_sum;
122}
123
124//----------------------------------------------------------------------
125// sum the (scalar) pt of the jets that pass the cuts
126double Selector::scalar_pt_sum(const std::vector<PseudoJet> & jets) const {
127 double this_sum = 0.0;
128 const SelectorWorker * worker_local = validated_worker();
129
130 // separate strategies according to whether the worker applies jet by jet
131 if (worker_local->applies_jet_by_jet()) {
132 for (unsigned i = 0; i < jets.size(); i++) {
133 if (worker_local->pass(jets[i])) this_sum += jets[i].pt();
134 }
135 } else {
136 std::vector<const PseudoJet *> jetptrs(jets.size());
137 for (unsigned i = 0; i < jets.size(); i++) {
138 jetptrs[i] = & jets[i];
139 }
140 worker_local->terminator(jetptrs);
141 for (unsigned i = 0; i < jetptrs.size(); i++) {
142 if (jetptrs[i]) this_sum += jets[i].pt();
143 }
144 }
145
146 return this_sum;
147}
148
149
150//----------------------------------------------------------------------
151// sift the input jets into two vectors -- those that pass the selector
152// and those that do not
153void Selector::sift(const std::vector<PseudoJet> & jets,
154 std::vector<PseudoJet> & jets_that_pass,
155 std::vector<PseudoJet> & jets_that_fail
156 ) const {
157 const SelectorWorker * worker_local = validated_worker();
158
159 jets_that_pass.clear();
160 jets_that_fail.clear();
161
162 // separate strategies according to whether the worker applies jet by jet
163 if (worker_local->applies_jet_by_jet()) {
164 for (unsigned i = 0; i < jets.size(); i++) {
165 if (worker_local->pass(jets[i])) {
166 jets_that_pass.push_back(jets[i]);
167 } else {
168 jets_that_fail.push_back(jets[i]);
169 }
170 }
171 } else {
172 std::vector<const PseudoJet *> jetptrs(jets.size());
173 for (unsigned i = 0; i < jets.size(); i++) {
174 jetptrs[i] = & jets[i];
175 }
176 worker_local->terminator(jetptrs);
177 for (unsigned i = 0; i < jetptrs.size(); i++) {
178 if (jetptrs[i]) {
179 jets_that_pass.push_back(jets[i]);
180 } else {
181 jets_that_fail.push_back(jets[i]);
182 }
183 }
184 }
185}
186
187#ifndef __FJCORE__
188// area using default ghost area
189double Selector::area() const{
190 return area(gas::def_ghost_area);
191}
192
193// implementation of the Selector's area function
194double Selector::area(double ghost_area) const{
195 if (! is_geometric()) throw InvalidArea();
196
197 // has area will already check we've got a valid worker
198 if (_worker->has_known_area()) return _worker->known_area();
199
200 // generate a set of "ghosts"
201 double rapmin, rapmax;
202 get_rapidity_extent(rapmin, rapmax);
203 GhostedAreaSpec ghost_spec(rapmin, rapmax, 1, ghost_area);
204 std::vector<PseudoJet> ghosts;
205 ghost_spec.add_ghosts(ghosts);
206
207 // check what passes the selection
208 return ghost_spec.ghost_area() * ((*this)(ghosts)).size();
209}
210#endif // __FJCORE__
211
212
213//----------------------------------------------------------------------
214// implementations of some of the more complex bits of SelectorWorker
215//----------------------------------------------------------------------
216// check if it has a finite area
217bool SelectorWorker::has_finite_area() const {
218 if (! is_geometric()) return false;
219 double rapmin, rapmax;
220 get_rapidity_extent(rapmin, rapmax);
221 return (rapmax != std::numeric_limits<double>::infinity())
222 && (-rapmin != std::numeric_limits<double>::infinity());
223}
224
225
226
227//----------------------------------------------------------------------
228// very basic set of selectors (at the moment just the identity!)
229//----------------------------------------------------------------------
230
231//----------------------------------------------------------------------
232/// helper for selecting the n hardest jets
233class SW_Identity : public SelectorWorker {
234public:
235 /// ctor with specification of the number of objects to keep
236 SW_Identity(){}
237
238 /// just let everything pass
239 virtual bool pass(const PseudoJet &) const {
240 return true;
241 }
242
243 /// For each jet that does not pass the cuts, this routine sets the
244 /// pointer to 0.
245 virtual void terminator(vector<const PseudoJet *> &) const {
246 // everything passes, hence nothing to nullify
247 return;
248 }
249
250 /// returns a description of the worker
251 virtual string description() const { return "Identity";}
252
253 /// strictly speaking, this is geometric
254 virtual bool is_geometric() const { return true;}
255};
256
257
258// returns an "identity" selector that lets everything pass
259Selector SelectorIdentity() {
260 return Selector(new SW_Identity);
261}
262
263
264//----------------------------------------------------------------------
265// selector and workers for operators
266//----------------------------------------------------------------------
267
268//----------------------------------------------------------------------
269/// helper for combining selectors with a logical not
270class SW_Not : public SelectorWorker {
271public:
272 /// ctor
273 SW_Not(const Selector & s) : _s(s) {}
274
275 /// return a copy of the current object
276 virtual SelectorWorker* copy(){ return new SW_Not(*this);}
277
278 /// returns true if a given object passes the selection criterium
279 /// this has to be overloaded by derived workers
280 virtual bool pass(const PseudoJet & jet) const {
281 // make sure that the "pass" can be applied on both selectors
282 if (!applies_jet_by_jet())
283 throw Error("Cannot apply this selector worker to an individual jet");
284
285 return ! _s.pass(jet);
286 }
287
288 /// returns true if this can be applied jet by jet
289 virtual bool applies_jet_by_jet() const {return _s.applies_jet_by_jet();}
290
291 /// select the jets in the list that pass both selectors
292 virtual void terminator(vector<const PseudoJet *> & jets) const {
293 // if we can apply the selector jet-by-jet, call the base selector
294 // that does exactly that
295 if (applies_jet_by_jet()){
296 SelectorWorker::terminator(jets);
297 return;
298 }
299
300 // check the effect of the selector we want to negate
301 vector<const PseudoJet *> s_jets = jets;
302 _s.worker()->terminator(s_jets);
303
304 // now apply the negation: all the jets that pass the base
305 // selector (i.e. are not NULL) have to be set to NULL
306 for (unsigned int i=0; i<s_jets.size(); i++){
307 if (s_jets[i]) jets[i] = NULL;
308 }
309 }
310
311 /// returns a description of the worker
312 virtual string description() const {
313 ostringstream ostr;
314 ostr << "!(" << _s.description() << ")";
315 return ostr.str();
316 }
317
318 /// is geometric if the underlying selector is
319 virtual bool is_geometric() const { return _s.is_geometric();}
320
321 /// returns true if the worker can be set_referenced
322 virtual bool takes_reference() const { return _s.takes_reference();}
323
324 /// set the reference jet for this selector
325 virtual void set_reference(const PseudoJet &ref) { _s.set_reference(ref);}
326
327protected:
328 Selector _s;
329};
330
331
332// logical not applied on a selector
333Selector operator!(const Selector & s) {
334 return Selector(new SW_Not(s));
335}
336
337
338//----------------------------------------------------------------------
339/// Base class for binary operators
340class SW_BinaryOperator: public SelectorWorker {
341public:
342 /// ctor
343 SW_BinaryOperator(const Selector & s1, const Selector & s2) : _s1(s1), _s2(s2) {
344 // stores info for more efficient access to the selector's properties
345
346 // we can apply jet by jet only if this is the case for both sub-selectors
347 _applies_jet_by_jet = _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
348
349 // the selector takes a reference if either of the sub-selectors does
350 _takes_reference = _s1.takes_reference() || _s2.takes_reference();
351
352 // we have a well-defined area provided the two objects have one
353 _is_geometric = _s1.is_geometric() && _s2.is_geometric();
354 }
355
356 /// returns true if this can be applied jet by jet
357 virtual bool applies_jet_by_jet() const {return _applies_jet_by_jet;}
358
359 /// returns true if this takes a reference jet
360 virtual bool takes_reference() const{
361 return _takes_reference;
362 }
363
364 /// sets the reference jet
365 virtual void set_reference(const PseudoJet &centre){
366 _s1.set_reference(centre);
367 _s2.set_reference(centre);
368 }
369
370 /// check if it has a finite area
371 virtual bool is_geometric() const { return _is_geometric;}
372
373protected:
374 Selector _s1, _s2;
375 bool _applies_jet_by_jet;
376 bool _takes_reference;
377 bool _is_geometric;
378};
379
380
381
382//----------------------------------------------------------------------
383/// helper for combining selectors with a logical and
384class SW_And: public SW_BinaryOperator {
385public:
386 /// ctor
387 SW_And(const Selector & s1, const Selector & s2) : SW_BinaryOperator(s1,s2){}
388
389 /// return a copy of this
390 virtual SelectorWorker* copy(){ return new SW_And(*this);}
391
392 /// returns true if a given object passes the selection criterium
393 /// this has to be overloaded by derived workers
394 virtual bool pass(const PseudoJet & jet) const {
395 // make sure that the "pass" can be applied on both selectors
396 if (!applies_jet_by_jet())
397 throw Error("Cannot apply this selector worker to an individual jet");
398
399 return _s1.pass(jet) && _s2.pass(jet);
400 }
401
402 /// select the jets in the list that pass both selectors
403 virtual void terminator(vector<const PseudoJet *> & jets) const {
404 // if we can apply the selector jet-by-jet, call the base selector
405 // that does exactly that
406 if (applies_jet_by_jet()){
407 SelectorWorker::terminator(jets);
408 return;
409 }
410
411 // check the effect of the first selector
412 vector<const PseudoJet *> s1_jets = jets;
413 _s1.worker()->terminator(s1_jets);
414
415 // apply the second
416 _s2.worker()->terminator(jets);
417
418 // terminate the jets that wiuld be terminated by _s1
419 for (unsigned int i=0; i<jets.size(); i++){
420 if (! s1_jets[i]) jets[i] = NULL;
421 }
422 }
423
424 /// returns the rapidity range for which it may return "true"
425 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const {
426 double s1min, s1max, s2min, s2max;
427 _s1.get_rapidity_extent(s1min, s1max);
428 _s2.get_rapidity_extent(s2min, s2max);
429 rapmax = min(s1max, s2max);
430 rapmin = max(s1min, s2min);
431 }
432
433 /// returns a description of the worker
434 virtual string description() const {
435 ostringstream ostr;
436 ostr << "(" << _s1.description() << " && " << _s2.description() << ")";
437 return ostr.str();
438 }
439};
440
441
442// logical and between two selectors
443Selector operator&&(const Selector & s1, const Selector & s2) {
444 return Selector(new SW_And(s1,s2));
445}
446
447
448
449//----------------------------------------------------------------------
450/// helper for combining selectors with a logical or
451class SW_Or: public SW_BinaryOperator {
452public:
453 /// ctor
454 SW_Or(const Selector & s1, const Selector & s2) : SW_BinaryOperator(s1,s2) {}
455
456 /// return a copy of this
457 virtual SelectorWorker* copy(){ return new SW_Or(*this);}
458
459 /// returns true if a given object passes the selection criterium
460 /// this has to be overloaded by derived workers
461 virtual bool pass(const PseudoJet & jet) const {
462 // make sure that the "pass" can be applied on both selectors
463 if (!applies_jet_by_jet())
464 throw Error("Cannot apply this selector worker to an individual jet");
465
466 return _s1.pass(jet) || _s2.pass(jet);
467 }
468
469 /// returns true if this can be applied jet by jet
470 virtual bool applies_jet_by_jet() const {
471 // watch out, even though it's the "OR" selector, to be applied jet
472 // by jet, both the base selectors need to be jet-by-jet-applicable,
473 // so the use of a && in the line below
474 return _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet();
475 }
476
477 /// select the jets in the list that pass both selectors
478 virtual void terminator(vector<const PseudoJet *> & jets) const {
479 // if we can apply the selector jet-by-jet, call the base selector
480 // that does exactly that
481 if (applies_jet_by_jet()){
482 SelectorWorker::terminator(jets);
483 return;
484 }
485
486 // check the effect of the first selector
487 vector<const PseudoJet *> s1_jets = jets;
488 _s1.worker()->terminator(s1_jets);
489
490 // apply the second
491 _s2.worker()->terminator(jets);
492
493 // resurrect any jet that has been terminated by the second one
494 // and not by the first one
495 for (unsigned int i=0; i<jets.size(); i++){
496 if (s1_jets[i]) jets[i] = s1_jets[i];
497 }
498 }
499
500 /// returns a description of the worker
501 virtual string description() const {
502 ostringstream ostr;
503 ostr << "(" << _s1.description() << " || " << _s2.description() << ")";
504 return ostr.str();
505 }
506
507 /// returns the rapidity range for which it may return "true"
508 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const {
509 double s1min, s1max, s2min, s2max;
510 _s1.get_rapidity_extent(s1min, s1max);
511 _s2.get_rapidity_extent(s2min, s2max);
512 rapmax = max(s1max, s2max);
513 rapmin = min(s1min, s2min);
514 }
515};
516
517
518// logical or between two selectors
519Selector operator ||(const Selector & s1, const Selector & s2) {
520 return Selector(new SW_Or(s1,s2));
521}
522
523//----------------------------------------------------------------------
524/// helper for multiplying two selectors (in an operator-like way)
525class SW_Mult: public SW_And {
526public:
527 /// ctor
528 SW_Mult(const Selector & s1, const Selector & s2) : SW_And(s1,s2) {}
529
530 /// return a copy of this
531 virtual SelectorWorker* copy(){ return new SW_Mult(*this);}
532
533 /// select the jets in the list that pass both selectors
534 virtual void terminator(vector<const PseudoJet *> & jets) const {
535 // if we can apply the selector jet-by-jet, call the base selector
536 // that does exactly that
537 if (applies_jet_by_jet()){
538 SelectorWorker::terminator(jets);
539 return;
540 }
541
542 // first apply _s2
543 _s2.worker()->terminator(jets);
544
545 // then apply _s1
546 _s1.worker()->terminator(jets);
547 }
548
549 /// returns a description of the worker
550 virtual string description() const {
551 ostringstream ostr;
552 ostr << "(" << _s1.description() << " * " << _s2.description() << ")";
553 return ostr.str();
554 }
555};
556
557
558// logical and between two selectors
559Selector operator*(const Selector & s1, const Selector & s2) {
560 return Selector(new SW_Mult(s1,s2));
561}
562
563
564//----------------------------------------------------------------------
565// selector and workers for kinematic cuts
566//----------------------------------------------------------------------
567
568//----------------------------------------------------------------------
569// a series of basic classes that allow easy implementations of
570// min, max and ranges on a quantity-to-be-defined
571
572// generic holder for a quantity
573class QuantityBase{
574public:
575 QuantityBase(double q) : _q(q){}
576 virtual ~QuantityBase(){}
577 virtual double operator()(const PseudoJet & jet ) const =0;
578 virtual string description() const =0;
579 virtual bool is_geometric() const { return false;}
580 virtual double comparison_value() const {return _q;}
581 virtual double description_value() const {return comparison_value();}
582protected:
583 double _q;
584};
585
586// generic holder for a squared quantity
587class QuantitySquareBase : public QuantityBase{
588public:
589 QuantitySquareBase(double sqrtq) : QuantityBase(sqrtq*sqrtq), _sqrtq(sqrtq){}
590 virtual double description_value() const {return _sqrtq;}
591protected:
592 double _sqrtq;
593};
594
595// generic_quantity >= minimum
596template<typename QuantityType>
597class SW_QuantityMin : public SelectorWorker{
598public:
599 /// detfault ctor (initialises the pt cut)
600 SW_QuantityMin(double qmin) : _qmin(qmin) {}
601
602 /// returns true is the given object passes the selection pt cut
603 virtual bool pass(const PseudoJet & jet) const {return _qmin(jet) >= _qmin.comparison_value();}
604
605 /// returns a description of the worker
606 virtual string description() const {
607 ostringstream ostr;
608 ostr << _qmin.description() << " >= " << _qmin.description_value();
609 return ostr.str();
610 }
611
612 virtual bool is_geometric() const { return _qmin.is_geometric();}
613
614protected:
615 QuantityType _qmin; ///< the cut
616};
617
618
619// generic_quantity <= maximum
620template<typename QuantityType>
621class SW_QuantityMax : public SelectorWorker {
622public:
623 /// detfault ctor (initialises the pt cut)
624 SW_QuantityMax(double qmax) : _qmax(qmax) {}
625
626 /// returns true is the given object passes the selection pt cut
627 virtual bool pass(const PseudoJet & jet) const {return _qmax(jet) <= _qmax.comparison_value();}
628
629 /// returns a description of the worker
630 virtual string description() const {
631 ostringstream ostr;
632 ostr << _qmax.description() << " <= " << _qmax.description_value();
633 return ostr.str();
634 }
635
636 virtual bool is_geometric() const { return _qmax.is_geometric();}
637
638protected:
639 QuantityType _qmax; ///< the cut
640};
641
642
643// generic quantity in [minimum:maximum]
644template<typename QuantityType>
645class SW_QuantityRange : public SelectorWorker {
646public:
647 /// detfault ctor (initialises the pt cut)
648 SW_QuantityRange(double qmin, double qmax) : _qmin(qmin), _qmax(qmax) {}
649
650 /// returns true is the given object passes the selection pt cut
651 virtual bool pass(const PseudoJet & jet) const {
652 double q = _qmin(jet); // we could identically use _qmax
653 return (q >= _qmin.comparison_value()) && (q <= _qmax.comparison_value());
654 }
655
656 /// returns a description of the worker
657 virtual string description() const {
658 ostringstream ostr;
659 ostr << _qmin.description_value() << " <= " << _qmin.description() << " <= " << _qmax.description_value();
660 return ostr.str();
661 }
662
663 virtual bool is_geometric() const { return _qmin.is_geometric();}
664
665protected:
666 QuantityType _qmin; // the lower cut
667 QuantityType _qmax; // the upper cut
668};
669
670
671//----------------------------------------------------------------------
672/// helper class for selecting on pt
673class QuantityPt2 : public QuantitySquareBase{
674public:
675 QuantityPt2(double pt) : QuantitySquareBase(pt){}
676 virtual double operator()(const PseudoJet & jet ) const { return jet.perp2();}
677 virtual string description() const {return "pt";}
678};
679
680// returns a selector for a minimum pt
681Selector SelectorPtMin(double ptmin) {
682 return Selector(new SW_QuantityMin<QuantityPt2>(ptmin));
683}
684
685// returns a selector for a maximum pt
686Selector SelectorPtMax(double ptmax) {
687 return Selector(new SW_QuantityMax<QuantityPt2>(ptmax));
688}
689
690// returns a selector for a pt range
691Selector SelectorPtRange(double ptmin, double ptmax) {
692 return Selector(new SW_QuantityRange<QuantityPt2>(ptmin, ptmax));
693}
694
695
696//----------------------------------------------------------------------
697/// helper class for selecting on transverse energy
698class QuantityEt2 : public QuantitySquareBase{
699public:
700 QuantityEt2(double Et) : QuantitySquareBase(Et){}
701 virtual double operator()(const PseudoJet & jet ) const { return jet.Et2();}
702 virtual string description() const {return "Et";}
703};
704
705// returns a selector for a minimum Et
706Selector SelectorEtMin(double Etmin) {
707 return Selector(new SW_QuantityMin<QuantityEt2>(Etmin));
708}
709
710// returns a selector for a maximum Et
711Selector SelectorEtMax(double Etmax) {
712 return Selector(new SW_QuantityMax<QuantityEt2>(Etmax));
713}
714
715// returns a selector for a Et range
716Selector SelectorEtRange(double Etmin, double Etmax) {
717 return Selector(new SW_QuantityRange<QuantityEt2>(Etmin, Etmax));
718}
719
720
721//----------------------------------------------------------------------
722/// helper class for selecting on energy
723class QuantityE : public QuantityBase{
724public:
725 QuantityE(double E) : QuantityBase(E){}
726 virtual double operator()(const PseudoJet & jet ) const { return jet.E();}
727 virtual string description() const {return "E";}
728};
729
730// returns a selector for a minimum E
731Selector SelectorEMin(double Emin) {
732 return Selector(new SW_QuantityMin<QuantityE>(Emin));
733}
734
735// returns a selector for a maximum E
736Selector SelectorEMax(double Emax) {
737 return Selector(new SW_QuantityMax<QuantityE>(Emax));
738}
739
740// returns a selector for a E range
741Selector SelectorERange(double Emin, double Emax) {
742 return Selector(new SW_QuantityRange<QuantityE>(Emin, Emax));
743}
744
745
746//----------------------------------------------------------------------
747/// helper class for selecting on mass
748class QuantityM2 : public QuantitySquareBase{
749public:
750 QuantityM2(double m) : QuantitySquareBase(m){}
751 virtual double operator()(const PseudoJet & jet ) const { return jet.m2();}
752 virtual string description() const {return "mass";}
753};
754
755// returns a selector for a minimum mass
756Selector SelectorMassMin(double mmin) {
757 return Selector(new SW_QuantityMin<QuantityM2>(mmin));
758}
759
760// returns a selector for a maximum mass
761Selector SelectorMassMax(double mmax) {
762 return Selector(new SW_QuantityMax<QuantityM2>(mmax));
763}
764
765// returns a selector for a mass range
766Selector SelectorMassRange(double mmin, double mmax) {
767 return Selector(new SW_QuantityRange<QuantityM2>(mmin, mmax));
768}
769
770
771
772//----------------------------------------------------------------------
773/// helper for selecting on rapidities: quantity
774class QuantityRap : public QuantityBase{
775public:
776 QuantityRap(double rap) : QuantityBase(rap){}
777 virtual double operator()(const PseudoJet & jet ) const { return jet.rap();}
778 virtual string description() const {return "rap";}
779 virtual bool is_geometric() const { return true;}
780};
781
782
783/// helper for selecting on rapidities: min
784class SW_RapMin : public SW_QuantityMin<QuantityRap>{
785public:
786 SW_RapMin(double rapmin) : SW_QuantityMin<QuantityRap>(rapmin){}
787 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{
788 rapmax = std::numeric_limits<double>::max();
789 rapmin = _qmin.comparison_value();
790 }
791};
792
793/// helper for selecting on rapidities: max
794class SW_RapMax : public SW_QuantityMax<QuantityRap>{
795public:
796 SW_RapMax(double rapmax) : SW_QuantityMax<QuantityRap>(rapmax){}
797 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{
798 rapmax = _qmax.comparison_value();
799 rapmin = -std::numeric_limits<double>::max();
800 }
801};
802
803/// helper for selecting on rapidities: range
804class SW_RapRange : public SW_QuantityRange<QuantityRap>{
805public:
806 SW_RapRange(double rapmin, double rapmax) : SW_QuantityRange<QuantityRap>(rapmin, rapmax){
807 assert(rapmin<=rapmax);
808 }
809 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{
810 rapmax = _qmax.comparison_value();
811 rapmin = _qmin.comparison_value();
812 }
813 virtual bool has_known_area() const { return true;} ///< the area is analytically known
814 virtual double known_area() const {
815 return twopi * (_qmax.comparison_value()-_qmin.comparison_value());
816 }
817};
818
819// returns a selector for a minimum rapidity
820Selector SelectorRapMin(double rapmin) {
821 return Selector(new SW_RapMin(rapmin));
822}
823
824// returns a selector for a maximum rapidity
825Selector SelectorRapMax(double rapmax) {
826 return Selector(new SW_RapMax(rapmax));
827}
828
829// returns a selector for a rapidity range
830Selector SelectorRapRange(double rapmin, double rapmax) {
831 return Selector(new SW_RapRange(rapmin, rapmax));
832}
833
834
835//----------------------------------------------------------------------
836/// helper for selecting on |rapidities|
837class QuantityAbsRap : public QuantityBase{
838public:
839 QuantityAbsRap(double absrap) : QuantityBase(absrap){}
840 virtual double operator()(const PseudoJet & jet ) const { return abs(jet.rap());}
841 virtual string description() const {return "|rap|";}
842 virtual bool is_geometric() const { return true;}
843};
844
845
846/// helper for selecting on |rapidities|: max
847class SW_AbsRapMax : public SW_QuantityMax<QuantityAbsRap>{
848public:
849 SW_AbsRapMax(double absrapmax) : SW_QuantityMax<QuantityAbsRap>(absrapmax){}
850 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{
851 rapmax = _qmax.comparison_value();
852 rapmin = -_qmax.comparison_value();
853 }
854 virtual bool has_known_area() const { return true;} ///< the area is analytically known
855 virtual double known_area() const {
856 return twopi * 2 * _qmax.comparison_value();
857 }
858};
859
860/// helper for selecting on |rapidities|: max
861class SW_AbsRapRange : public SW_QuantityRange<QuantityAbsRap>{
862public:
863 SW_AbsRapRange(double absrapmin, double absrapmax) : SW_QuantityRange<QuantityAbsRap>(absrapmin, absrapmax){}
864 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{
865 rapmax = _qmax.comparison_value();
866 rapmin = -_qmax.comparison_value();
867 }
868 virtual bool has_known_area() const { return true;} ///< the area is analytically known
869 virtual double known_area() const {
870 return twopi * 2 * (_qmax.comparison_value()-max(_qmin.comparison_value(),0.0)); // this should handle properly absrapmin<0
871 }
872};
873
874// returns a selector for a minimum |rapidity|
875Selector SelectorAbsRapMin(double absrapmin) {
876 return Selector(new SW_QuantityMin<QuantityAbsRap>(absrapmin));
877}
878
879// returns a selector for a maximum |rapidity|
880Selector SelectorAbsRapMax(double absrapmax) {
881 return Selector(new SW_AbsRapMax(absrapmax));
882}
883
884// returns a selector for a |rapidity| range
885Selector SelectorAbsRapRange(double rapmin, double rapmax) {
886 return Selector(new SW_AbsRapRange(rapmin, rapmax));
887}
888
889
890//----------------------------------------------------------------------
891/// helper for selecting on pseudo-rapidities
892class QuantityEta : public QuantityBase{
893public:
894 QuantityEta(double eta) : QuantityBase(eta){}
895 virtual double operator()(const PseudoJet & jet ) const { return jet.eta();}
896 virtual string description() const {return "eta";}
897 // virtual bool is_geometric() const { return true;} // not strictly only y and phi-dependent
898};
899
900// returns a selector for a pseudo-minimum rapidity
901Selector SelectorEtaMin(double etamin) {
902 return Selector(new SW_QuantityMin<QuantityEta>(etamin));
903}
904
905// returns a selector for a pseudo-maximum rapidity
906Selector SelectorEtaMax(double etamax) {
907 return Selector(new SW_QuantityMax<QuantityEta>(etamax));
908}
909
910// returns a selector for a pseudo-rapidity range
911Selector SelectorEtaRange(double etamin, double etamax) {
912 return Selector(new SW_QuantityRange<QuantityEta>(etamin, etamax));
913}
914
915
916//----------------------------------------------------------------------
917/// helper for selecting on |pseudo-rapidities|
918class QuantityAbsEta : public QuantityBase{
919public:
920 QuantityAbsEta(double abseta) : QuantityBase(abseta){}
921 virtual double operator()(const PseudoJet & jet ) const { return abs(jet.eta());}
922 virtual string description() const {return "|eta|";}
923 virtual bool is_geometric() const { return true;}
924};
925
926// returns a selector for a minimum |pseudo-rapidity|
927Selector SelectorAbsEtaMin(double absetamin) {
928 return Selector(new SW_QuantityMin<QuantityAbsEta>(absetamin));
929}
930
931// returns a selector for a maximum |pseudo-rapidity|
932Selector SelectorAbsEtaMax(double absetamax) {
933 return Selector(new SW_QuantityMax<QuantityAbsEta>(absetamax));
934}
935
936// returns a selector for a |pseudo-rapidity| range
937Selector SelectorAbsEtaRange(double absetamin, double absetamax) {
938 return Selector(new SW_QuantityRange<QuantityAbsEta>(absetamin, absetamax));
939}
940
941
942//----------------------------------------------------------------------
943/// helper for selecting on azimuthal angle
944///
945/// Note that the bounds have to be specified as min<max
946/// phimin has to be > -2pi
947/// phimax has to be < 4pi
948class SW_PhiRange : public SelectorWorker {
949public:
950 /// detfault ctor (initialises the pt cut)
951 SW_PhiRange(double phimin, double phimax) : _phimin(phimin), _phimax(phimax){
952 assert(_phimin<_phimax);
953 assert(_phimin>-twopi);
954 assert(_phimax<2*twopi);
955
956 _phispan = _phimax - _phimin;
957 }
958
959 /// returns true is the given object passes the selection pt cut
960 virtual bool pass(const PseudoJet & jet) const {
961 double dphi=jet.phi()-_phimin;
962 if (dphi >= twopi) dphi -= twopi;
963 if (dphi < 0) dphi += twopi;
964 return (dphi <= _phispan);
965 }
966
967 /// returns a description of the worker
968 virtual string description() const {
969 ostringstream ostr;
970 ostr << _phimin << " <= phi <= " << _phimax;
971 return ostr.str();
972 }
973
974 virtual bool is_geometric() const { return true;}
975
976protected:
977 double _phimin; // the lower cut
978 double _phimax; // the upper cut
979 double _phispan; // the span of the range
980};
981
982
983// returns a selector for a phi range
984Selector SelectorPhiRange(double phimin, double phimax) {
985 return Selector(new SW_PhiRange(phimin, phimax));
986}
987
988//----------------------------------------------------------------------
989/// helper for selecting on both rapidity and azimuthal angle
990class SW_RapPhiRange : public SW_And{
991public:
992 SW_RapPhiRange(double rapmin, double rapmax, double phimin, double phimax)
993 : SW_And(SelectorRapRange(rapmin, rapmax), SelectorPhiRange(phimin, phimax)){
994 _known_area = ((phimax-phimin > twopi) ? twopi : phimax-phimin) * (rapmax-rapmin);
995 }
996
997 /// if it has a computable area, return it
998 virtual double known_area() const{
999 return _known_area;
1000 }
1001
1002protected:
1003 double _known_area;
1004};
1005
1006Selector SelectorRapPhiRange(double rapmin, double rapmax, double phimin, double phimax) {
1007 return Selector(new SW_RapPhiRange(rapmin, rapmax, phimin, phimax));
1008}
1009
1010
1011//----------------------------------------------------------------------
1012/// helper for selecting the n hardest jets
1013class SW_NHardest : public SelectorWorker {
1014public:
1015 /// ctor with specification of the number of objects to keep
1016 SW_NHardest(unsigned int n) : _n(n) {};
1017
1018 /// pass makes no sense here normally the parent selector will throw
1019 /// an error but for internal use in the SW, we'll throw one from
1020 /// here by security
1021 virtual bool pass(const PseudoJet &) const {
1022 if (!applies_jet_by_jet())
1023 throw Error("Cannot apply this selector worker to an individual jet");
1024 return false;
1025 }
1026
1027 /// For each jet that does not pass the cuts, this routine sets the
1028 /// pointer to 0.
1029 virtual void terminator(vector<const PseudoJet *> & jets) const {
1030 // nothing to do if the size is too small
1031 if (jets.size() < _n) return;
1032
1033 // do we want to first chech if things are already ordered before
1034 // going through the ordering process? For now, no. Maybe carry
1035 // out timing tests at some point to establish the optimal
1036 // strategy.
1037
1038 vector<double> minus_pt2(jets.size());
1039 vector<unsigned int> indices(jets.size());
1040
1041 for (unsigned int i=0; i<jets.size(); i++){
1042 indices[i] = i;
1043
1044 // we need to make sure that the object has not already been
1045 // nullified. Note that if we have less than _n jets, this
1046 // whole n-hardest selection will not have any effect.
1047 minus_pt2[i] = jets[i] ? -jets[i]->perp2() : 0.0;
1048 }
1049
1050 IndexedSortHelper sort_helper(& minus_pt2);
1051
1052 partial_sort(indices.begin(), indices.begin()+_n, indices.end(), sort_helper);
1053
1054 for (unsigned int i=_n; i<jets.size(); i++)
1055 jets[indices[i]] = NULL;
1056 }
1057
1058 /// returns true if this can be applied jet by jet
1059 virtual bool applies_jet_by_jet() const {return false;}
1060
1061 /// returns a description of the worker
1062 virtual string description() const {
1063 ostringstream ostr;
1064 ostr << _n << " hardest";
1065 return ostr.str();
1066 }
1067
1068protected:
1069 unsigned int _n;
1070};
1071
1072
1073// returns a selector for the n hardest jets
1074Selector SelectorNHardest(unsigned int n) {
1075 return Selector(new SW_NHardest(n));
1076}
1077
1078
1079
1080//----------------------------------------------------------------------
1081// selector and workers for geometric ranges
1082//----------------------------------------------------------------------
1083
1084//----------------------------------------------------------------------
1085/// a generic class for objects that contain a position
1086class SW_WithReference : public SelectorWorker{
1087public:
1088 /// ctor
1089 SW_WithReference() : _is_initialised(false){};
1090
1091 /// returns true if the worker takes a reference jet
1092 virtual bool takes_reference() const { return true;}
1093
1094 /// sets the reference jet
1095 virtual void set_reference(const PseudoJet &centre){
1096 _is_initialised = true;
1097 _reference = centre;
1098 }
1099
1100protected:
1101 PseudoJet _reference;
1102 bool _is_initialised;
1103};
1104
1105//----------------------------------------------------------------------
1106/// helper for selecting on objects within a distance 'radius' of a reference
1107class SW_Circle : public SW_WithReference {
1108public:
1109 SW_Circle(const double radius) : _radius2(radius*radius) {}
1110
1111 /// return a copy of the current object
1112 virtual SelectorWorker* copy(){ return new SW_Circle(*this);}
1113
1114 /// returns true if a given object passes the selection criterium
1115 /// this has to be overloaded by derived workers
1116 virtual bool pass(const PseudoJet & jet) const {
1117 // make sure the centre is initialised
1118 if (! _is_initialised)
1119 throw Error("To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1120
1121 return jet.squared_distance(_reference) <= _radius2;
1122 }
1123
1124 /// returns a description of the worker
1125 virtual string description() const {
1126 ostringstream ostr;
1127 ostr << "distance from the centre <= " << sqrt(_radius2);
1128 return ostr.str();
1129 }
1130
1131 /// returns the rapidity range for which it may return "true"
1132 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{
1133 // make sure the centre is initialised
1134 if (! _is_initialised)
1135 throw Error("To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)");
1136
1137 rapmax = _reference.rap()+sqrt(_radius2);
1138 rapmin = _reference.rap()-sqrt(_radius2);
1139 }
1140
1141 virtual bool is_geometric() const { return true;} ///< implies a finite area
1142 virtual bool has_finite_area() const { return true;} ///< regardless of the reference
1143 virtual bool has_known_area() const { return true;} ///< the area is analytically known
1144 virtual double known_area() const {
1145 return pi * _radius2;
1146 }
1147
1148protected:
1149 double _radius2;
1150};
1151
1152
1153// select on objets within a distance 'radius' of a variable location
1154Selector SelectorCircle(const double radius) {
1155 return Selector(new SW_Circle(radius));
1156}
1157
1158
1159//----------------------------------------------------------------------
1160/// helper for selecting on objects with a distance to a reference
1161/// betwene 'radius_in' and 'radius_out'
1162class SW_Doughnut : public SW_WithReference {
1163public:
1164 SW_Doughnut(const double radius_in, const double radius_out)
1165 : _radius_in2(radius_in*radius_in), _radius_out2(radius_out*radius_out) {}
1166
1167 /// return a copy of the current object
1168 virtual SelectorWorker* copy(){ return new SW_Doughnut(*this);}
1169
1170 /// returns true if a given object passes the selection criterium
1171 /// this has to be overloaded by derived workers
1172 virtual bool pass(const PseudoJet & jet) const {
1173 // make sure the centre is initialised
1174 if (! _is_initialised)
1175 throw Error("To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1176
1177 double distance2 = jet.squared_distance(_reference);
1178
1179 return (distance2 <= _radius_out2) && (distance2 >= _radius_in2);
1180 }
1181
1182 /// returns a description of the worker
1183 virtual string description() const {
1184 ostringstream ostr;
1185 ostr << sqrt(_radius_in2) << " <= distance from the centre <= " << sqrt(_radius_out2);
1186 return ostr.str();
1187 }
1188
1189 /// returns the rapidity range for which it may return "true"
1190 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{
1191 // make sure the centre is initialised
1192 if (! _is_initialised)
1193 throw Error("To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)");
1194
1195 rapmax = _reference.rap()+sqrt(_radius_out2);
1196 rapmin = _reference.rap()-sqrt(_radius_out2);
1197 }
1198
1199 virtual bool is_geometric() const { return true;} ///< implies a finite area
1200 virtual bool has_finite_area() const { return true;} ///< regardless of the reference
1201 virtual bool has_known_area() const { return true;} ///< the area is analytically known
1202 virtual double known_area() const {
1203 return pi * (_radius_out2-_radius_in2);
1204 }
1205
1206protected:
1207 double _radius_in2, _radius_out2;
1208};
1209
1210
1211
1212// select on objets with distance from the centre is between 'radius_in' and 'radius_out'
1213Selector SelectorDoughnut(const double radius_in, const double radius_out) {
1214 return Selector(new SW_Doughnut(radius_in, radius_out));
1215}
1216
1217
1218//----------------------------------------------------------------------
1219/// helper for selecting on objects with rapidity within a distance 'delta' of a reference
1220class SW_Strip : public SW_WithReference {
1221public:
1222 SW_Strip(const double delta) : _delta(delta) {}
1223
1224 /// return a copy of the current object
1225 virtual SelectorWorker* copy(){ return new SW_Strip(*this);}
1226
1227 /// returns true if a given object passes the selection criterium
1228 /// this has to be overloaded by derived workers
1229 virtual bool pass(const PseudoJet & jet) const {
1230 // make sure the centre is initialised
1231 if (! _is_initialised)
1232 throw Error("To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1233
1234 return abs(jet.rap()-_reference.rap()) <= _delta;
1235 }
1236
1237 /// returns a description of the worker
1238 virtual string description() const {
1239 ostringstream ostr;
1240 ostr << "|rap - rap_reference| <= " << _delta;
1241 return ostr.str();
1242 }
1243
1244 /// returns the rapidity range for which it may return "true"
1245 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{
1246 // make sure the centre is initialised
1247 if (! _is_initialised)
1248 throw Error("To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)");
1249
1250 rapmax = _reference.rap()+_delta;
1251 rapmin = _reference.rap()-_delta;
1252 }
1253
1254 virtual bool is_geometric() const { return true;} ///< implies a finite area
1255 virtual bool has_finite_area() const { return true;} ///< regardless of the reference
1256 virtual bool has_known_area() const { return true;} ///< the area is analytically known
1257 virtual double known_area() const {
1258 return twopi * 2 * _delta;
1259 }
1260
1261protected:
1262 double _delta;
1263};
1264
1265
1266// select on objets within a distance 'radius' of a variable location
1267Selector SelectorStrip(const double half_width) {
1268 return Selector(new SW_Strip(half_width));
1269}
1270
1271
1272//----------------------------------------------------------------------
1273/// helper for selecting on objects with rapidity within a distance
1274/// 'delta_rap' of a reference and phi within a distanve delta_phi of
1275/// a reference
1276class SW_Rectangle : public SW_WithReference {
1277public:
1278 SW_Rectangle(const double delta_rap, const double delta_phi)
1279 : _delta_rap(delta_rap), _delta_phi(delta_phi) {}
1280
1281 /// return a copy of the current object
1282 virtual SelectorWorker* copy(){ return new SW_Rectangle(*this);}
1283
1284 /// returns true if a given object passes the selection criterium
1285 /// this has to be overloaded by derived workers
1286 virtual bool pass(const PseudoJet & jet) const {
1287 // make sure the centre is initialised
1288 if (! _is_initialised)
1289 throw Error("To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1290
1291 return (abs(jet.rap()-_reference.rap()) <= _delta_rap) && (abs(jet.delta_phi_to(_reference)) <= _delta_phi);
1292 }
1293
1294 /// returns a description of the worker
1295 virtual string description() const {
1296 ostringstream ostr;
1297 ostr << "|rap - rap_reference| <= " << _delta_rap << " && |phi - phi_reference| <= " << _delta_phi ;
1298 return ostr.str();
1299 }
1300
1301 /// returns the rapidity range for which it may return "true"
1302 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{
1303 // make sure the centre is initialised
1304 if (! _is_initialised)
1305 throw Error("To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)");
1306
1307 rapmax = _reference.rap()+_delta_rap;
1308 rapmin = _reference.rap()-_delta_rap;
1309 }
1310
1311 virtual bool is_geometric() const { return true;} ///< implies a finite area
1312 virtual bool has_finite_area() const { return true;} ///< regardless of the reference
1313 virtual bool has_known_area() const { return true;} ///< the area is analytically known
1314 virtual double known_area() const {
1315 return 4 * _delta_rap * _delta_phi;
1316 }
1317
1318protected:
1319 double _delta_rap, _delta_phi;
1320};
1321
1322
1323// select on objets within a distance 'radius' of a variable location
1324Selector SelectorRectangle(const double half_rap_width, const double half_phi_width) {
1325 return Selector(new SW_Rectangle(half_rap_width, half_phi_width));
1326}
1327
1328
1329//----------------------------------------------------------------------
1330/// helper for selecting the jets that carry at least a given fraction
1331/// of the reference jet
1332class SW_PtFractionMin : public SW_WithReference {
1333public:
1334 /// ctor with specification of the number of objects to keep
1335 SW_PtFractionMin(double fraction) : _fraction2(fraction*fraction){}
1336
1337 /// return a copy of the current object
1338 virtual SelectorWorker* copy(){ return new SW_PtFractionMin(*this);}
1339
1340 /// return true if the jet carries a large enough fraction of the reference.
1341 /// Throw an error if the reference is not initialised.
1342 virtual bool pass(const PseudoJet & jet) const {
1343 // make sure the centre is initialised
1344 if (! _is_initialised)
1345 throw Error("To use a SelectorPtFractionMin (or any selector that requires a reference), you first have to call set_reference(...)");
1346
1347 // otherwise, just call that method on the jet
1348 return (jet.perp2() >= _fraction2*_reference.perp2());
1349 }
1350
1351 /// returns a description of the worker
1352 virtual string description() const {
1353 ostringstream ostr;
1354 ostr << "pt >= " << sqrt(_fraction2) << "* pt_ref";
1355 return ostr.str();
1356 }
1357
1358protected:
1359 double _fraction2;
1360};
1361
1362
1363// select objects that carry at least a fraction "fraction" of the reference jet
1364// (Note that this selectir takes a reference)
1365Selector SelectorPtFractionMin(double fraction){
1366 return Selector(new SW_PtFractionMin(fraction));
1367}
1368
1369
1370//----------------------------------------------------------------------
1371// additional (mostly helper) selectors
1372//----------------------------------------------------------------------
1373
1374//----------------------------------------------------------------------
1375/// helper for selecting the 0-momentum jets
1376class SW_IsZero : public SelectorWorker {
1377public:
1378 /// ctor
1379 SW_IsZero(){}
1380
1381 /// return true if the jet has zero momentum
1382 virtual bool pass(const PseudoJet & jet) const {
1383 return jet==0;
1384 }
1385
1386 /// rereturns a description of the worker
1387 virtual string description() const { return "zero";}
1388};
1389
1390
1391// select objects with zero momentum
1392Selector SelectorIsZero(){
1393 return Selector(new SW_IsZero());
1394}
1395
1396
1397//----------------------------------------------------------------------
1398#ifndef __FJCORE__
1399/// helper for selecting the pure ghost
1400class SW_IsPureGhost : public SelectorWorker {
1401public:
1402 /// ctor
1403 SW_IsPureGhost(){}
1404
1405 /// return true if the jet is a pure-ghost jet
1406 virtual bool pass(const PseudoJet & jet) const {
1407 // if the jet has no area support then it's certainly not a ghost
1408 if (!jet.has_area()) return false;
1409
1410 // otherwise, just call that method on the jet
1411 return jet.is_pure_ghost();
1412 }
1413
1414 /// rereturns a description of the worker
1415 virtual string description() const { return "pure ghost";}
1416};
1417
1418
1419// select objects that are (or are only made of) ghosts
1420Selector SelectorIsPureGhost(){
1421 return Selector(new SW_IsPureGhost());
1422}
1423
1424//----------------------------------------------------------------------
1425// Selector and workers for obtaining a Selector from an old
1426// RangeDefinition
1427//
1428// This is mostly intended for backward compatibility and is likely to
1429// be removed in a future major release of FastJet
1430//----------------------------------------------------------------------
1431
1432//----------------------------------------------------------------------
1433/// helper for selecting on both rapidity and azimuthal angle
1434class SW_RangeDefinition : public SelectorWorker{
1435public:
1436 /// ctor from a RangeDefinition
1437 SW_RangeDefinition(const RangeDefinition &range) : _range(&range){}
1438
1439 /// transfer the selection creterium to the underlying RangeDefinition
1440 virtual bool pass(const PseudoJet & jet) const {
1441 return _range->is_in_range(jet);
1442 }
1443
1444 /// returns a description of the worker
1445 virtual string description() const {
1446 return _range->description();
1447 }
1448
1449 /// returns the rapidity range for which it may return "true"
1450 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{
1451 _range->get_rap_limits(rapmin, rapmax);
1452 }
1453
1454 /// check if it has a finite area
1455 virtual bool is_geometric() const { return true;}
1456
1457 /// check if it has an analytically computable area
1458 virtual bool has_known_area() const { return true;}
1459
1460 /// if it has a computable area, return it
1461 virtual double known_area() const{
1462 return _range->area();
1463 }
1464
1465protected:
1466 const RangeDefinition *_range;
1467};
1468
1469
1470// ctor from a RangeDefinition
1471//----------------------------------------------------------------------
1472//
1473// This is provided for backward compatibility and will be removed in
1474// a future major release of FastJet
1475Selector::Selector(const RangeDefinition &range) {
1476 _worker.reset(new SW_RangeDefinition(range));
1477}
1478#endif // __FJCORE__
1479
1480
1481// operators applying directly on a Selector
1482//----------------------------------------------------------------------
1483
1484// operator &=
1485// For 2 Selectors a and b, a &= b is eauivalent to a = a & b;
1486Selector & Selector::operator &=(const Selector & b){
1487 _worker.reset(new SW_And(*this, b));
1488 return *this;
1489}
1490
1491// operator &=
1492// For 2 Selectors a and b, a &= b is eauivalent to a = a & b;
1493Selector & Selector::operator |=(const Selector & b){
1494 _worker.reset(new SW_Or(*this, b));
1495 return *this;
1496}
1497
1498FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
Note: See TracBrowser for help on using the repository browser.