FastJet 3.0.2
|
00001 //STARTHEADER 00002 // $Id: Selector.cc 2687 2011-11-14 11:17:51Z soyez $ 00003 // 00004 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez 00005 // 00006 //---------------------------------------------------------------------- 00007 // This file is part of FastJet. 00008 // 00009 // FastJet is free software; you can redistribute it and/or modify 00010 // it under the terms of the GNU General Public License as published by 00011 // the Free Software Foundation; either version 2 of the License, or 00012 // (at your option) any later version. 00013 // 00014 // The algorithms that underlie FastJet have required considerable 00015 // development and are described in hep-ph/0512210. If you use 00016 // FastJet as part of work towards a scientific publication, please 00017 // include a citation to the FastJet paper. 00018 // 00019 // FastJet is distributed in the hope that it will be useful, 00020 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00021 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00022 // GNU General Public License for more details. 00023 // 00024 // You should have received a copy of the GNU General Public License 00025 // along with FastJet. If not, see <http://www.gnu.org/licenses/>. 00026 //---------------------------------------------------------------------- 00027 //ENDHEADER 00028 00029 00030 #include <sstream> 00031 #include <algorithm> 00032 #include "fastjet/Selector.hh" 00033 #include "fastjet/GhostedAreaSpec.hh" // for area support 00034 00035 using namespace std; 00036 00037 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh 00038 00039 //---------------------------------------------------------------------- 00040 // implementations of some of the more complex bits of Selector 00041 //---------------------------------------------------------------------- 00042 00043 // implementation of the operator() acting on a vector of jets 00044 std::vector<PseudoJet> Selector::operator()(const std::vector<PseudoJet> & jets) const { 00045 std::vector<PseudoJet> result; 00046 const SelectorWorker * worker_local = validated_worker(); 00047 if (worker_local->applies_jet_by_jet()) { 00048 //if (false) { 00049 // for workers that apply jet by jet, this is more efficient 00050 for (std::vector<PseudoJet>::const_iterator jet = jets.begin(); 00051 jet != jets.end(); jet++) { 00052 if (worker_local->pass(*jet)) result.push_back(*jet); 00053 } 00054 } else { 00055 // for workers that can only be applied to entire vectors, 00056 // go through the following 00057 std::vector<const PseudoJet *> jetptrs(jets.size()); 00058 for (unsigned i = 0; i < jets.size(); i++) { 00059 jetptrs[i] = & jets[i]; 00060 } 00061 worker_local->terminator(jetptrs); 00062 for (unsigned i = 0; i < jetptrs.size(); i++) { 00063 if (jetptrs[i]) result.push_back(jets[i]); 00064 } 00065 } 00066 return result; 00067 } 00068 00069 00070 //---------------------------------------------------------------------- 00071 // count the number of jets that pass the cuts 00072 unsigned int Selector::count(const std::vector<PseudoJet> & jets) const { 00073 unsigned n = 0; 00074 const SelectorWorker * worker_local = validated_worker(); 00075 00076 // separate strategies according to whether the worker applies jet by jet 00077 if (worker_local->applies_jet_by_jet()) { 00078 for (unsigned i = 0; i < jets.size(); i++) { 00079 if (worker_local->pass(jets[i])) n++; 00080 } 00081 } else { 00082 std::vector<const PseudoJet *> jetptrs(jets.size()); 00083 for (unsigned i = 0; i < jets.size(); i++) { 00084 jetptrs[i] = & jets[i]; 00085 } 00086 worker_local->terminator(jetptrs); 00087 for (unsigned i = 0; i < jetptrs.size(); i++) { 00088 if (jetptrs[i]) n++; 00089 } 00090 } 00091 00092 return n; 00093 } 00094 00095 00096 //---------------------------------------------------------------------- 00097 // sift the input jets into two vectors -- those that pass the selector 00098 // and those that do not 00099 void Selector::sift(const std::vector<PseudoJet> & jets, 00100 std::vector<PseudoJet> & jets_that_pass, 00101 std::vector<PseudoJet> & jets_that_fail 00102 ) const { 00103 const SelectorWorker * worker_local = validated_worker(); 00104 00105 jets_that_pass.clear(); 00106 jets_that_fail.clear(); 00107 00108 // separate strategies according to whether the worker applies jet by jet 00109 if (worker_local->applies_jet_by_jet()) { 00110 for (unsigned i = 0; i < jets.size(); i++) { 00111 if (worker_local->pass(jets[i])) { 00112 jets_that_pass.push_back(jets[i]); 00113 } else { 00114 jets_that_fail.push_back(jets[i]); 00115 } 00116 } 00117 } else { 00118 std::vector<const PseudoJet *> jetptrs(jets.size()); 00119 for (unsigned i = 0; i < jets.size(); i++) { 00120 jetptrs[i] = & jets[i]; 00121 } 00122 worker_local->terminator(jetptrs); 00123 for (unsigned i = 0; i < jetptrs.size(); i++) { 00124 if (jetptrs[i]) { 00125 jets_that_pass.push_back(jets[i]); 00126 } else { 00127 jets_that_fail.push_back(jets[i]); 00128 } 00129 } 00130 } 00131 } 00132 00133 // area using default ghost area 00134 double Selector::area() const{ 00135 return area(gas::def_ghost_area); 00136 } 00137 00138 // implementation of the Selector's area function 00139 double Selector::area(double ghost_area) const{ 00140 if (! is_geometric()) throw InvalidArea(); 00141 00142 // has area will already check we've got a valid worker 00143 if (_worker->has_known_area()) return _worker->known_area(); 00144 00145 // generate a set of "ghosts" 00146 double rapmin, rapmax; 00147 get_rapidity_extent(rapmin, rapmax); 00148 GhostedAreaSpec ghost_spec(rapmin, rapmax, 1, ghost_area); 00149 std::vector<PseudoJet> ghosts; 00150 ghost_spec.add_ghosts(ghosts); 00151 00152 // check what passes the selection 00153 return ghost_spec.ghost_area() * ((*this)(ghosts)).size(); 00154 } 00155 00156 00157 //---------------------------------------------------------------------- 00158 // implementations of some of the more complex bits of SelectorWorker 00159 //---------------------------------------------------------------------- 00160 // check if it has a finite area 00161 bool SelectorWorker::has_finite_area() const { 00162 if (! is_geometric()) return false; 00163 double rapmin, rapmax; 00164 get_rapidity_extent(rapmin, rapmax); 00165 return (rapmax != std::numeric_limits<double>::infinity()) 00166 && (-rapmin != std::numeric_limits<double>::infinity()); 00167 } 00168 00169 00170 00171 //---------------------------------------------------------------------- 00172 // very basic set of selectors (at the moment just the identity!) 00173 //---------------------------------------------------------------------- 00174 00175 //---------------------------------------------------------------------- 00176 /// helper for selecting the n hardest jets 00177 class SW_Identity : public SelectorWorker { 00178 public: 00179 /// ctor with specification of the number of objects to keep 00180 SW_Identity(){} 00181 00182 /// just let everything pass 00183 virtual bool pass(const PseudoJet &) const { 00184 return true; 00185 } 00186 00187 /// For each jet that does not pass the cuts, this routine sets the 00188 /// pointer to 0. 00189 virtual void terminator(vector<const PseudoJet *> &) const { 00190 // everything passes, hence nothing to nullify 00191 return; 00192 } 00193 00194 /// returns a description of the worker 00195 virtual string description() const { return "Identity";} 00196 00197 /// strictly speaking, this is geometric 00198 virtual bool is_geometric() const { return true;} 00199 }; 00200 00201 00202 // returns an "identity" selector that lets everything pass 00203 Selector SelectorIdentity() { 00204 return Selector(new SW_Identity); 00205 } 00206 00207 00208 //---------------------------------------------------------------------- 00209 // selector and workers for operators 00210 //---------------------------------------------------------------------- 00211 00212 //---------------------------------------------------------------------- 00213 /// helper for combining selectors with a logical not 00214 class SW_Not : public SelectorWorker { 00215 public: 00216 /// ctor 00217 SW_Not(const Selector & s) : _s(s) {} 00218 00219 /// return a copy of the current object 00220 virtual SelectorWorker* copy(){ return new SW_Not(*this);} 00221 00222 /// returns true if a given object passes the selection criterium 00223 /// this has to be overloaded by derived workers 00224 virtual bool pass(const PseudoJet & jet) const { 00225 // make sure that the "pass" can be applied on both selectors 00226 if (!applies_jet_by_jet()) 00227 throw Error("Cannot apply this selector worker to an individual jet"); 00228 00229 return ! _s.pass(jet); 00230 } 00231 00232 /// returns true if this can be applied jet by jet 00233 virtual bool applies_jet_by_jet() const {return _s.applies_jet_by_jet();} 00234 00235 /// select the jets in the list that pass both selectors 00236 virtual void terminator(vector<const PseudoJet *> & jets) const { 00237 // if we can apply the selector jet-by-jet, call the base selector 00238 // that does exactly that 00239 if (applies_jet_by_jet()){ 00240 SelectorWorker::terminator(jets); 00241 return; 00242 } 00243 00244 // check the effect of the selector we want to negate 00245 vector<const PseudoJet *> s_jets = jets; 00246 _s.worker()->terminator(s_jets); 00247 00248 // now apply the negation: all the jets that pass the base 00249 // selector (i.e. are not NULL) have to be set to NULL 00250 for (unsigned int i=0; i<s_jets.size(); i++){ 00251 if (s_jets[i]) jets[i] = NULL; 00252 } 00253 } 00254 00255 /// returns a description of the worker 00256 virtual string description() const { 00257 ostringstream ostr; 00258 ostr << "!(" << _s.description() << ")"; 00259 return ostr.str(); 00260 } 00261 00262 /// is geometric if the underlying selector is 00263 virtual bool is_geometric() const { return _s.is_geometric();} 00264 00265 /// returns true if the worker can be set_referenced 00266 virtual bool takes_reference() const { return _s.takes_reference();} 00267 00268 /// set the reference jet for this selector 00269 virtual void set_reference(const PseudoJet &ref) { _s.set_reference(ref);} 00270 00271 protected: 00272 Selector _s; 00273 }; 00274 00275 00276 // logical not applied on a selector 00277 Selector operator!(const Selector & s) { 00278 return Selector(new SW_Not(s)); 00279 } 00280 00281 00282 //---------------------------------------------------------------------- 00283 /// Base class for binary operators 00284 class SW_BinaryOperator: public SelectorWorker { 00285 public: 00286 /// ctor 00287 SW_BinaryOperator(const Selector & s1, const Selector & s2) : _s1(s1), _s2(s2) { 00288 // stores info for more efficient access to the selector's properties 00289 00290 // we can apply jet by jet only if this is the case for both sub-selectors 00291 _applies_jet_by_jet = _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet(); 00292 00293 // the selector takes a reference if either of the sub-selectors does 00294 _takes_reference = _s1.takes_reference() || _s2.takes_reference(); 00295 00296 // we have a well-defined area provided the two objects have one 00297 _is_geometric = _s1.is_geometric() && _s2.is_geometric(); 00298 } 00299 00300 /// returns true if this can be applied jet by jet 00301 virtual bool applies_jet_by_jet() const {return _applies_jet_by_jet;} 00302 00303 /// returns true if this takes a reference jet 00304 virtual bool takes_reference() const{ 00305 return _takes_reference; 00306 } 00307 00308 /// sets the reference jet 00309 virtual void set_reference(const PseudoJet ¢re){ 00310 _s1.set_reference(centre); 00311 _s2.set_reference(centre); 00312 } 00313 00314 /// check if it has a finite area 00315 virtual bool is_geometric() const { return _is_geometric;} 00316 00317 protected: 00318 Selector _s1, _s2; 00319 bool _applies_jet_by_jet; 00320 bool _takes_reference; 00321 bool _is_geometric; 00322 }; 00323 00324 00325 00326 //---------------------------------------------------------------------- 00327 /// helper for combining selectors with a logical and 00328 class SW_And: public SW_BinaryOperator { 00329 public: 00330 /// ctor 00331 SW_And(const Selector & s1, const Selector & s2) : SW_BinaryOperator(s1,s2){} 00332 00333 /// return a copy of this 00334 virtual SelectorWorker* copy(){ return new SW_And(*this);} 00335 00336 /// returns true if a given object passes the selection criterium 00337 /// this has to be overloaded by derived workers 00338 virtual bool pass(const PseudoJet & jet) const { 00339 // make sure that the "pass" can be applied on both selectors 00340 if (!applies_jet_by_jet()) 00341 throw Error("Cannot apply this selector worker to an individual jet"); 00342 00343 return _s1.pass(jet) && _s2.pass(jet); 00344 } 00345 00346 /// select the jets in the list that pass both selectors 00347 virtual void terminator(vector<const PseudoJet *> & jets) const { 00348 // if we can apply the selector jet-by-jet, call the base selector 00349 // that does exactly that 00350 if (applies_jet_by_jet()){ 00351 SelectorWorker::terminator(jets); 00352 return; 00353 } 00354 00355 // check the effect of the first selector 00356 vector<const PseudoJet *> s1_jets = jets; 00357 _s1.worker()->terminator(s1_jets); 00358 00359 // apply the second 00360 _s2.worker()->terminator(jets); 00361 00362 // terminate the jets that wiuld be terminated by _s1 00363 for (unsigned int i=0; i<jets.size(); i++){ 00364 if (! s1_jets[i]) jets[i] = NULL; 00365 } 00366 } 00367 00368 /// returns the rapidity range for which it may return "true" 00369 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const { 00370 double s1min, s1max, s2min, s2max; 00371 _s1.get_rapidity_extent(s1min, s1max); 00372 _s2.get_rapidity_extent(s2min, s2max); 00373 rapmax = min(s1max, s2max); 00374 rapmin = max(s1min, s2min); 00375 } 00376 00377 /// returns a description of the worker 00378 virtual string description() const { 00379 ostringstream ostr; 00380 ostr << "(" << _s1.description() << " && " << _s2.description() << ")"; 00381 return ostr.str(); 00382 } 00383 }; 00384 00385 00386 // logical and between two selectors 00387 Selector operator&&(const Selector & s1, const Selector & s2) { 00388 return Selector(new SW_And(s1,s2)); 00389 } 00390 00391 00392 00393 //---------------------------------------------------------------------- 00394 /// helper for combining selectors with a logical or 00395 class SW_Or: public SW_BinaryOperator { 00396 public: 00397 /// ctor 00398 SW_Or(const Selector & s1, const Selector & s2) : SW_BinaryOperator(s1,s2) {} 00399 00400 /// return a copy of this 00401 virtual SelectorWorker* copy(){ return new SW_Or(*this);} 00402 00403 /// returns true if a given object passes the selection criterium 00404 /// this has to be overloaded by derived workers 00405 virtual bool pass(const PseudoJet & jet) const { 00406 // make sure that the "pass" can be applied on both selectors 00407 if (!applies_jet_by_jet()) 00408 throw Error("Cannot apply this selector worker to an individual jet"); 00409 00410 return _s1.pass(jet) || _s2.pass(jet); 00411 } 00412 00413 /// returns true if this can be applied jet by jet 00414 virtual bool applies_jet_by_jet() const { 00415 // watch out, even though it's the "OR" selector, to be applied jet 00416 // by jet, both the base selectors need to be jet-by-jet-applicable, 00417 // so the use of a && in the line below 00418 return _s1.applies_jet_by_jet() && _s2.applies_jet_by_jet(); 00419 } 00420 00421 /// select the jets in the list that pass both selectors 00422 virtual void terminator(vector<const PseudoJet *> & jets) const { 00423 // if we can apply the selector jet-by-jet, call the base selector 00424 // that does exactly that 00425 if (applies_jet_by_jet()){ 00426 SelectorWorker::terminator(jets); 00427 return; 00428 } 00429 00430 // check the effect of the first selector 00431 vector<const PseudoJet *> s1_jets = jets; 00432 _s1.worker()->terminator(s1_jets); 00433 00434 // apply the second 00435 _s2.worker()->terminator(jets); 00436 00437 // resurrect any jet that has been terminated by the second one 00438 // and not by the first one 00439 for (unsigned int i=0; i<jets.size(); i++){ 00440 if (s1_jets[i]) jets[i] = s1_jets[i]; 00441 } 00442 } 00443 00444 /// returns a description of the worker 00445 virtual string description() const { 00446 ostringstream ostr; 00447 ostr << "(" << _s1.description() << " || " << _s2.description() << ")"; 00448 return ostr.str(); 00449 } 00450 00451 /// returns the rapidity range for which it may return "true" 00452 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const { 00453 double s1min, s1max, s2min, s2max; 00454 _s1.get_rapidity_extent(s1min, s1max); 00455 _s2.get_rapidity_extent(s2min, s2max); 00456 rapmax = max(s1max, s2max); 00457 rapmin = min(s1min, s2min); 00458 } 00459 }; 00460 00461 00462 // logical or between two selectors 00463 Selector operator ||(const Selector & s1, const Selector & s2) { 00464 return Selector(new SW_Or(s1,s2)); 00465 } 00466 00467 //---------------------------------------------------------------------- 00468 /// helper for multiplying two selectors (in an operator-like way) 00469 class SW_Mult: public SW_And { 00470 public: 00471 /// ctor 00472 SW_Mult(const Selector & s1, const Selector & s2) : SW_And(s1,s2) {} 00473 00474 /// return a copy of this 00475 virtual SelectorWorker* copy(){ return new SW_Mult(*this);} 00476 00477 /// select the jets in the list that pass both selectors 00478 virtual void terminator(vector<const PseudoJet *> & jets) const { 00479 // if we can apply the selector jet-by-jet, call the base selector 00480 // that does exactly that 00481 if (applies_jet_by_jet()){ 00482 SelectorWorker::terminator(jets); 00483 return; 00484 } 00485 00486 // first apply _s2 00487 _s2.worker()->terminator(jets); 00488 00489 // then apply _s1 00490 _s1.worker()->terminator(jets); 00491 } 00492 00493 /// returns a description of the worker 00494 virtual string description() const { 00495 ostringstream ostr; 00496 ostr << "(" << _s1.description() << " * " << _s2.description() << ")"; 00497 return ostr.str(); 00498 } 00499 }; 00500 00501 00502 // logical and between two selectors 00503 Selector operator*(const Selector & s1, const Selector & s2) { 00504 return Selector(new SW_Mult(s1,s2)); 00505 } 00506 00507 00508 //---------------------------------------------------------------------- 00509 // selector and workers for kinematic cuts 00510 //---------------------------------------------------------------------- 00511 00512 //---------------------------------------------------------------------- 00513 // a series of basic classes that allow easy implementations of 00514 // min, max and ranges on a quantity-to-be-defined 00515 00516 // generic holder for a quantity 00517 class QuantityBase{ 00518 public: 00519 QuantityBase(double q) : _q(q){} 00520 virtual ~QuantityBase(){} 00521 virtual double operator()(const PseudoJet & jet ) const =0; 00522 virtual string description() const =0; 00523 virtual bool is_geometric() const { return false;} 00524 virtual double comparison_value() const {return _q;} 00525 virtual double description_value() const {return comparison_value();} 00526 protected: 00527 double _q; 00528 }; 00529 00530 // generic holder for a squared quantity 00531 class QuantitySquareBase : public QuantityBase{ 00532 public: 00533 QuantitySquareBase(double sqrtq) : QuantityBase(sqrtq*sqrtq), _sqrtq(sqrtq){} 00534 virtual double description_value() const {return _sqrtq;} 00535 protected: 00536 double _sqrtq; 00537 }; 00538 00539 // generic_quantity >= minimum 00540 template<typename QuantityType> 00541 class SW_QuantityMin : public SelectorWorker{ 00542 public: 00543 /// detfault ctor (initialises the pt cut) 00544 SW_QuantityMin(double qmin) : _qmin(qmin) {} 00545 00546 /// returns true is the given object passes the selection pt cut 00547 virtual bool pass(const PseudoJet & jet) const {return _qmin(jet) >= _qmin.comparison_value();} 00548 00549 /// returns a description of the worker 00550 virtual string description() const { 00551 ostringstream ostr; 00552 ostr << _qmin.description() << " >= " << _qmin.description_value(); 00553 return ostr.str(); 00554 } 00555 00556 virtual bool is_geometric() const { return _qmin.is_geometric();} 00557 00558 protected: 00559 QuantityType _qmin; ///< the cut 00560 }; 00561 00562 00563 // generic_quantity <= maximum 00564 template<typename QuantityType> 00565 class SW_QuantityMax : public SelectorWorker { 00566 public: 00567 /// detfault ctor (initialises the pt cut) 00568 SW_QuantityMax(double qmax) : _qmax(qmax) {} 00569 00570 /// returns true is the given object passes the selection pt cut 00571 virtual bool pass(const PseudoJet & jet) const {return _qmax(jet) <= _qmax.comparison_value();} 00572 00573 /// returns a description of the worker 00574 virtual string description() const { 00575 ostringstream ostr; 00576 ostr << _qmax.description() << " <= " << _qmax.description_value(); 00577 return ostr.str(); 00578 } 00579 00580 virtual bool is_geometric() const { return _qmax.is_geometric();} 00581 00582 protected: 00583 QuantityType _qmax; ///< the cut 00584 }; 00585 00586 00587 // generic quantity in [minimum:maximum] 00588 template<typename QuantityType> 00589 class SW_QuantityRange : public SelectorWorker { 00590 public: 00591 /// detfault ctor (initialises the pt cut) 00592 SW_QuantityRange(double qmin, double qmax) : _qmin(qmin), _qmax(qmax) {} 00593 00594 /// returns true is the given object passes the selection pt cut 00595 virtual bool pass(const PseudoJet & jet) const { 00596 double q = _qmin(jet); // we could identically use _qmax 00597 return (q >= _qmin.comparison_value()) && (q <= _qmax.comparison_value()); 00598 } 00599 00600 /// returns a description of the worker 00601 virtual string description() const { 00602 ostringstream ostr; 00603 ostr << _qmin.description_value() << " <= " << _qmin.description() << " <= " << _qmax.description_value(); 00604 return ostr.str(); 00605 } 00606 00607 virtual bool is_geometric() const { return _qmin.is_geometric();} 00608 00609 protected: 00610 QuantityType _qmin; // the lower cut 00611 QuantityType _qmax; // the upper cut 00612 }; 00613 00614 00615 //---------------------------------------------------------------------- 00616 /// helper class for selecting on pt 00617 class QuantityPt2 : public QuantitySquareBase{ 00618 public: 00619 QuantityPt2(double pt) : QuantitySquareBase(pt){} 00620 virtual double operator()(const PseudoJet & jet ) const { return jet.perp2();} 00621 virtual string description() const {return "pt";} 00622 }; 00623 00624 // returns a selector for a minimum pt 00625 Selector SelectorPtMin(double ptmin) { 00626 return Selector(new SW_QuantityMin<QuantityPt2>(ptmin)); 00627 } 00628 00629 // returns a selector for a maximum pt 00630 Selector SelectorPtMax(double ptmax) { 00631 return Selector(new SW_QuantityMax<QuantityPt2>(ptmax)); 00632 } 00633 00634 // returns a selector for a pt range 00635 Selector SelectorPtRange(double ptmin, double ptmax) { 00636 return Selector(new SW_QuantityRange<QuantityPt2>(ptmin, ptmax)); 00637 } 00638 00639 00640 //---------------------------------------------------------------------- 00641 /// helper class for selecting on transverse energy 00642 class QuantityEt2 : public QuantitySquareBase{ 00643 public: 00644 QuantityEt2(double Et) : QuantitySquareBase(Et){} 00645 virtual double operator()(const PseudoJet & jet ) const { return jet.Et2();} 00646 virtual string description() const {return "Et";} 00647 }; 00648 00649 // returns a selector for a minimum Et 00650 Selector SelectorEtMin(double Etmin) { 00651 return Selector(new SW_QuantityMin<QuantityEt2>(Etmin)); 00652 } 00653 00654 // returns a selector for a maximum Et 00655 Selector SelectorEtMax(double Etmax) { 00656 return Selector(new SW_QuantityMax<QuantityEt2>(Etmax)); 00657 } 00658 00659 // returns a selector for a Et range 00660 Selector SelectorEtRange(double Etmin, double Etmax) { 00661 return Selector(new SW_QuantityRange<QuantityEt2>(Etmin, Etmax)); 00662 } 00663 00664 00665 //---------------------------------------------------------------------- 00666 /// helper class for selecting on energy 00667 class QuantityE : public QuantityBase{ 00668 public: 00669 QuantityE(double E) : QuantityBase(E){} 00670 virtual double operator()(const PseudoJet & jet ) const { return jet.E();} 00671 virtual string description() const {return "E";} 00672 }; 00673 00674 // returns a selector for a minimum E 00675 Selector SelectorEMin(double Emin) { 00676 return Selector(new SW_QuantityMin<QuantityE>(Emin)); 00677 } 00678 00679 // returns a selector for a maximum E 00680 Selector SelectorEMax(double Emax) { 00681 return Selector(new SW_QuantityMax<QuantityE>(Emax)); 00682 } 00683 00684 // returns a selector for a E range 00685 Selector SelectorERange(double Emin, double Emax) { 00686 return Selector(new SW_QuantityRange<QuantityE>(Emin, Emax)); 00687 } 00688 00689 00690 //---------------------------------------------------------------------- 00691 /// helper class for selecting on mass 00692 class QuantityM2 : public QuantitySquareBase{ 00693 public: 00694 QuantityM2(double m) : QuantitySquareBase(m){} 00695 virtual double operator()(const PseudoJet & jet ) const { return jet.m2();} 00696 virtual string description() const {return "mass";} 00697 }; 00698 00699 // returns a selector for a minimum mass 00700 Selector SelectorMassMin(double mmin) { 00701 return Selector(new SW_QuantityMin<QuantityM2>(mmin)); 00702 } 00703 00704 // returns a selector for a maximum mass 00705 Selector SelectorMassMax(double mmax) { 00706 return Selector(new SW_QuantityMax<QuantityM2>(mmax)); 00707 } 00708 00709 // returns a selector for a mass range 00710 Selector SelectorMassRange(double mmin, double mmax) { 00711 return Selector(new SW_QuantityRange<QuantityM2>(mmin, mmax)); 00712 } 00713 00714 00715 00716 //---------------------------------------------------------------------- 00717 /// helper for selecting on rapidities: quantity 00718 class QuantityRap : public QuantityBase{ 00719 public: 00720 QuantityRap(double rap) : QuantityBase(rap){} 00721 virtual double operator()(const PseudoJet & jet ) const { return jet.rap();} 00722 virtual string description() const {return "rap";} 00723 virtual bool is_geometric() const { return true;} 00724 }; 00725 00726 00727 /// helper for selecting on rapidities: min 00728 class SW_RapMin : public SW_QuantityMin<QuantityRap>{ 00729 public: 00730 SW_RapMin(double rapmin) : SW_QuantityMin<QuantityRap>(rapmin){} 00731 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{ 00732 rapmax = std::numeric_limits<double>::max(); 00733 rapmin = _qmin.comparison_value(); 00734 } 00735 }; 00736 00737 /// helper for selecting on rapidities: max 00738 class SW_RapMax : public SW_QuantityMax<QuantityRap>{ 00739 public: 00740 SW_RapMax(double rapmax) : SW_QuantityMax<QuantityRap>(rapmax){} 00741 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{ 00742 rapmax = _qmax.comparison_value(); 00743 rapmin = -std::numeric_limits<double>::max(); 00744 } 00745 }; 00746 00747 /// helper for selecting on rapidities: range 00748 class SW_RapRange : public SW_QuantityRange<QuantityRap>{ 00749 public: 00750 SW_RapRange(double rapmin, double rapmax) : SW_QuantityRange<QuantityRap>(rapmin, rapmax){ 00751 assert(rapmin<=rapmax); 00752 } 00753 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{ 00754 rapmax = _qmax.comparison_value(); 00755 rapmin = _qmin.comparison_value(); 00756 } 00757 virtual bool has_known_area() const { return true;} ///< the area is analytically known 00758 virtual double known_area() const { 00759 return twopi * (_qmax.comparison_value()-_qmin.comparison_value()); 00760 } 00761 }; 00762 00763 // returns a selector for a minimum rapidity 00764 Selector SelectorRapMin(double rapmin) { 00765 return Selector(new SW_RapMin(rapmin)); 00766 } 00767 00768 // returns a selector for a maximum rapidity 00769 Selector SelectorRapMax(double rapmax) { 00770 return Selector(new SW_RapMax(rapmax)); 00771 } 00772 00773 // returns a selector for a rapidity range 00774 Selector SelectorRapRange(double rapmin, double rapmax) { 00775 return Selector(new SW_RapRange(rapmin, rapmax)); 00776 } 00777 00778 00779 //---------------------------------------------------------------------- 00780 /// helper for selecting on |rapidities| 00781 class QuantityAbsRap : public QuantityBase{ 00782 public: 00783 QuantityAbsRap(double absrap) : QuantityBase(absrap){} 00784 virtual double operator()(const PseudoJet & jet ) const { return abs(jet.rap());} 00785 virtual string description() const {return "|rap|";} 00786 virtual bool is_geometric() const { return true;} 00787 }; 00788 00789 00790 /// helper for selecting on |rapidities|: max 00791 class SW_AbsRapMax : public SW_QuantityMax<QuantityAbsRap>{ 00792 public: 00793 SW_AbsRapMax(double absrapmax) : SW_QuantityMax<QuantityAbsRap>(absrapmax){} 00794 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{ 00795 rapmax = _qmax.comparison_value(); 00796 rapmin = -_qmax.comparison_value(); 00797 } 00798 virtual bool has_known_area() const { return true;} ///< the area is analytically known 00799 virtual double known_area() const { 00800 return twopi * 2 * _qmax.comparison_value(); 00801 } 00802 }; 00803 00804 /// helper for selecting on |rapidities|: max 00805 class SW_AbsRapRange : public SW_QuantityRange<QuantityAbsRap>{ 00806 public: 00807 SW_AbsRapRange(double absrapmin, double absrapmax) : SW_QuantityRange<QuantityAbsRap>(absrapmin, absrapmax){} 00808 virtual void get_rapidity_extent(double &rapmin, double & rapmax) const{ 00809 rapmax = _qmax.comparison_value(); 00810 rapmin = -_qmax.comparison_value(); 00811 } 00812 virtual bool has_known_area() const { return true;} ///< the area is analytically known 00813 virtual double known_area() const { 00814 return twopi * 2 * (_qmax.comparison_value()-max(_qmin.comparison_value(),0.0)); // this should handle properly absrapmin<0 00815 } 00816 }; 00817 00818 // returns a selector for a minimum |rapidity| 00819 Selector SelectorAbsRapMin(double absrapmin) { 00820 return Selector(new SW_QuantityMin<QuantityAbsRap>(absrapmin)); 00821 } 00822 00823 // returns a selector for a maximum |rapidity| 00824 Selector SelectorAbsRapMax(double absrapmax) { 00825 return Selector(new SW_AbsRapMax(absrapmax)); 00826 } 00827 00828 // returns a selector for a |rapidity| range 00829 Selector SelectorAbsRapRange(double rapmin, double rapmax) { 00830 return Selector(new SW_AbsRapRange(rapmin, rapmax)); 00831 } 00832 00833 00834 //---------------------------------------------------------------------- 00835 /// helper for selecting on pseudo-rapidities 00836 class QuantityEta : public QuantityBase{ 00837 public: 00838 QuantityEta(double eta) : QuantityBase(eta){} 00839 virtual double operator()(const PseudoJet & jet ) const { return jet.eta();} 00840 virtual string description() const {return "eta";} 00841 // virtual bool is_geometric() const { return true;} // not strictly only y and phi-dependent 00842 }; 00843 00844 // returns a selector for a pseudo-minimum rapidity 00845 Selector SelectorEtaMin(double etamin) { 00846 return Selector(new SW_QuantityMin<QuantityEta>(etamin)); 00847 } 00848 00849 // returns a selector for a pseudo-maximum rapidity 00850 Selector SelectorEtaMax(double etamax) { 00851 return Selector(new SW_QuantityMax<QuantityEta>(etamax)); 00852 } 00853 00854 // returns a selector for a pseudo-rapidity range 00855 Selector SelectorEtaRange(double etamin, double etamax) { 00856 return Selector(new SW_QuantityRange<QuantityEta>(etamin, etamax)); 00857 } 00858 00859 00860 //---------------------------------------------------------------------- 00861 /// helper for selecting on |pseudo-rapidities| 00862 class QuantityAbsEta : public QuantityBase{ 00863 public: 00864 QuantityAbsEta(double abseta) : QuantityBase(abseta){} 00865 virtual double operator()(const PseudoJet & jet ) const { return abs(jet.eta());} 00866 virtual string description() const {return "|eta|";} 00867 virtual bool is_geometric() const { return true;} 00868 }; 00869 00870 // returns a selector for a minimum |pseudo-rapidity| 00871 Selector SelectorAbsEtaMin(double absetamin) { 00872 return Selector(new SW_QuantityMin<QuantityAbsEta>(absetamin)); 00873 } 00874 00875 // returns a selector for a maximum |pseudo-rapidity| 00876 Selector SelectorAbsEtaMax(double absetamax) { 00877 return Selector(new SW_QuantityMax<QuantityAbsEta>(absetamax)); 00878 } 00879 00880 // returns a selector for a |pseudo-rapidity| range 00881 Selector SelectorAbsEtaRange(double absetamin, double absetamax) { 00882 return Selector(new SW_QuantityRange<QuantityAbsEta>(absetamin, absetamax)); 00883 } 00884 00885 00886 //---------------------------------------------------------------------- 00887 /// helper for selecting on azimuthal angle 00888 /// 00889 /// Note that the bounds have to be specified as min<max 00890 /// phimin has to be > -2pi 00891 /// phimax has to be < 4pi 00892 class SW_PhiRange : public SelectorWorker { 00893 public: 00894 /// detfault ctor (initialises the pt cut) 00895 SW_PhiRange(double phimin, double phimax) : _phimin(phimin), _phimax(phimax){ 00896 assert(_phimin<_phimax); 00897 assert(_phimin>-twopi); 00898 assert(_phimax<2*twopi); 00899 00900 _phispan = _phimax - _phimin; 00901 } 00902 00903 /// returns true is the given object passes the selection pt cut 00904 virtual bool pass(const PseudoJet & jet) const { 00905 double dphi=jet.phi()-_phimin; 00906 if (dphi >= twopi) dphi -= twopi; 00907 if (dphi < 0) dphi += twopi; 00908 return (dphi <= _phispan); 00909 } 00910 00911 /// returns a description of the worker 00912 virtual string description() const { 00913 ostringstream ostr; 00914 ostr << _phimin << " <= phi <= " << _phimax; 00915 return ostr.str(); 00916 } 00917 00918 virtual bool is_geometric() const { return true;} 00919 00920 protected: 00921 double _phimin; // the lower cut 00922 double _phimax; // the upper cut 00923 double _phispan; // the span of the range 00924 }; 00925 00926 00927 // returns a selector for a phi range 00928 Selector SelectorPhiRange(double phimin, double phimax) { 00929 return Selector(new SW_PhiRange(phimin, phimax)); 00930 } 00931 00932 //---------------------------------------------------------------------- 00933 /// helper for selecting on both rapidity and azimuthal angle 00934 class SW_RapPhiRange : public SW_And{ 00935 public: 00936 SW_RapPhiRange(double rapmin, double rapmax, double phimin, double phimax) 00937 : SW_And(SelectorRapRange(rapmin, rapmax), SelectorPhiRange(phimin, phimax)){ 00938 _known_area = ((phimax-phimin > twopi) ? twopi : phimax-phimin) * (rapmax-rapmin); 00939 } 00940 00941 /// if it has a computable area, return it 00942 virtual double known_area() const{ 00943 return _known_area; 00944 } 00945 00946 protected: 00947 double _known_area; 00948 }; 00949 00950 Selector SelectorRapPhiRange(double rapmin, double rapmax, double phimin, double phimax) { 00951 return Selector(new SW_RapPhiRange(rapmin, rapmax, phimin, phimax)); 00952 } 00953 00954 00955 //---------------------------------------------------------------------- 00956 /// helper for selecting the n hardest jets 00957 class SW_NHardest : public SelectorWorker { 00958 public: 00959 /// ctor with specification of the number of objects to keep 00960 SW_NHardest(unsigned int n) : _n(n) {}; 00961 00962 /// pass makes no sense here normally the parent selector will throw 00963 /// an error but for internal use in the SW, we'll throw one from 00964 /// here by security 00965 virtual bool pass(const PseudoJet &) const { 00966 if (!applies_jet_by_jet()) 00967 throw Error("Cannot apply this selector worker to an individual jet"); 00968 return false; 00969 } 00970 00971 /// For each jet that does not pass the cuts, this routine sets the 00972 /// pointer to 0. 00973 virtual void terminator(vector<const PseudoJet *> & jets) const { 00974 // nothing to do if the size is too small 00975 if (jets.size() < _n) return; 00976 00977 // do we want to first chech if things are already ordered before 00978 // going through the ordering process? For now, no. Maybe carry 00979 // out timing tests at some point to establish the optimal 00980 // strategy. 00981 00982 vector<double> minus_pt2(jets.size()); 00983 vector<unsigned int> indices(jets.size()); 00984 00985 for (unsigned int i=0; i<jets.size(); i++){ 00986 indices[i] = i; 00987 00988 // we need to make sure that the object has not already been 00989 // nullified. Note that if we have less than _n jets, this 00990 // whole n-hardest selection will not have any effect. 00991 minus_pt2[i] = jets[i] ? -jets[i]->perp2() : 0.0; 00992 } 00993 00994 IndexedSortHelper sort_helper(& minus_pt2); 00995 00996 partial_sort(indices.begin(), indices.begin()+_n, indices.end(), sort_helper); 00997 00998 for (unsigned int i=_n; i<jets.size(); i++) 00999 jets[indices[i]] = NULL; 01000 } 01001 01002 /// returns true if this can be applied jet by jet 01003 virtual bool applies_jet_by_jet() const {return false;} 01004 01005 /// returns a description of the worker 01006 virtual string description() const { 01007 ostringstream ostr; 01008 ostr << _n << " hardest"; 01009 return ostr.str(); 01010 } 01011 01012 protected: 01013 unsigned int _n; 01014 }; 01015 01016 01017 // returns a selector for the n hardest jets 01018 Selector SelectorNHardest(unsigned int n) { 01019 return Selector(new SW_NHardest(n)); 01020 } 01021 01022 01023 01024 //---------------------------------------------------------------------- 01025 // selector and workers for geometric ranges 01026 //---------------------------------------------------------------------- 01027 01028 //---------------------------------------------------------------------- 01029 /// a generic class for objects that contain a position 01030 class SW_WithReference : public SelectorWorker{ 01031 public: 01032 /// ctor 01033 SW_WithReference() : _is_initialised(false){}; 01034 01035 /// returns true if the worker takes a reference jet 01036 virtual bool takes_reference() const { return true;} 01037 01038 /// sets the reference jet 01039 virtual void set_reference(const PseudoJet ¢re){ 01040 _is_initialised = true; 01041 _reference = centre; 01042 } 01043 01044 protected: 01045 PseudoJet _reference; 01046 bool _is_initialised; 01047 }; 01048 01049 //---------------------------------------------------------------------- 01050 /// helper for selecting on objects within a distance 'radius' of a reference 01051 class SW_Circle : public SW_WithReference { 01052 public: 01053 SW_Circle(const double &radius) : _radius2(radius*radius) {} 01054 01055 /// return a copy of the current object 01056 virtual SelectorWorker* copy(){ return new SW_Circle(*this);} 01057 01058 /// returns true if a given object passes the selection criterium 01059 /// this has to be overloaded by derived workers 01060 virtual bool pass(const PseudoJet & jet) const { 01061 // make sure the centre is initialised 01062 if (! _is_initialised) 01063 throw Error("To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)"); 01064 01065 return jet.squared_distance(_reference) <= _radius2; 01066 } 01067 01068 /// returns a description of the worker 01069 virtual string description() const { 01070 ostringstream ostr; 01071 ostr << "distance from the centre <= " << sqrt(_radius2); 01072 return ostr.str(); 01073 } 01074 01075 /// returns the rapidity range for which it may return "true" 01076 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{ 01077 // make sure the centre is initialised 01078 if (! _is_initialised) 01079 throw Error("To use a SelectorCircle (or any selector that requires a reference), you first have to call set_reference(...)"); 01080 01081 rapmax = _reference.rap()+sqrt(_radius2); 01082 rapmin = _reference.rap()-sqrt(_radius2); 01083 } 01084 01085 virtual bool is_geometric() const { return true;} ///< implies a finite area 01086 virtual bool has_finite_area() const { return true;} ///< regardless of the reference 01087 virtual bool has_known_area() const { return true;} ///< the area is analytically known 01088 virtual double known_area() const { 01089 return pi * _radius2; 01090 } 01091 01092 protected: 01093 double _radius2; 01094 }; 01095 01096 01097 // select on objets within a distance 'radius' of a variable location 01098 Selector SelectorCircle(const double & radius) { 01099 return Selector(new SW_Circle(radius)); 01100 } 01101 01102 01103 //---------------------------------------------------------------------- 01104 /// helper for selecting on objects with a distance to a reference 01105 /// betwene 'radius_in' and 'radius_out' 01106 class SW_Doughnut : public SW_WithReference { 01107 public: 01108 SW_Doughnut(const double &radius_in, const double &radius_out) 01109 : _radius_in2(radius_in*radius_in), _radius_out2(radius_out*radius_out) {} 01110 01111 /// return a copy of the current object 01112 virtual SelectorWorker* copy(){ return new SW_Doughnut(*this);} 01113 01114 /// returns true if a given object passes the selection criterium 01115 /// this has to be overloaded by derived workers 01116 virtual bool pass(const PseudoJet & jet) const { 01117 // make sure the centre is initialised 01118 if (! _is_initialised) 01119 throw Error("To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)"); 01120 01121 double distance2 = jet.squared_distance(_reference); 01122 01123 return (distance2 <= _radius_out2) && (distance2 >= _radius_in2); 01124 } 01125 01126 /// returns a description of the worker 01127 virtual string description() const { 01128 ostringstream ostr; 01129 ostr << sqrt(_radius_in2) << " <= distance from the centre <= " << sqrt(_radius_out2); 01130 return ostr.str(); 01131 } 01132 01133 /// returns the rapidity range for which it may return "true" 01134 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{ 01135 // make sure the centre is initialised 01136 if (! _is_initialised) 01137 throw Error("To use a SelectorDoughnut (or any selector that requires a reference), you first have to call set_reference(...)"); 01138 01139 rapmax = _reference.rap()+sqrt(_radius_out2); 01140 rapmin = _reference.rap()-sqrt(_radius_out2); 01141 } 01142 01143 virtual bool is_geometric() const { return true;} ///< implies a finite area 01144 virtual bool has_finite_area() const { return true;} ///< regardless of the reference 01145 virtual bool has_known_area() const { return true;} ///< the area is analytically known 01146 virtual double known_area() const { 01147 return pi * (_radius_out2-_radius_in2); 01148 } 01149 01150 protected: 01151 double _radius_in2, _radius_out2; 01152 }; 01153 01154 01155 01156 // select on objets with distance from the centre is between 'radius_in' and 'radius_out' 01157 Selector SelectorDoughnut(const double & radius_in, const double & radius_out) { 01158 return Selector(new SW_Doughnut(radius_in, radius_out)); 01159 } 01160 01161 01162 //---------------------------------------------------------------------- 01163 /// helper for selecting on objects with rapidity within a distance 'delta' of a reference 01164 class SW_Strip : public SW_WithReference { 01165 public: 01166 SW_Strip(const double &delta) : _delta(delta) {} 01167 01168 /// return a copy of the current object 01169 virtual SelectorWorker* copy(){ return new SW_Strip(*this);} 01170 01171 /// returns true if a given object passes the selection criterium 01172 /// this has to be overloaded by derived workers 01173 virtual bool pass(const PseudoJet & jet) const { 01174 // make sure the centre is initialised 01175 if (! _is_initialised) 01176 throw Error("To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)"); 01177 01178 return abs(jet.rap()-_reference.rap()) <= _delta; 01179 } 01180 01181 /// returns a description of the worker 01182 virtual string description() const { 01183 ostringstream ostr; 01184 ostr << "|rap - rap_reference| <= " << _delta; 01185 return ostr.str(); 01186 } 01187 01188 /// returns the rapidity range for which it may return "true" 01189 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{ 01190 // make sure the centre is initialised 01191 if (! _is_initialised) 01192 throw Error("To use a SelectorStrip (or any selector that requires a reference), you first have to call set_reference(...)"); 01193 01194 rapmax = _reference.rap()+_delta; 01195 rapmin = _reference.rap()-_delta; 01196 } 01197 01198 virtual bool is_geometric() const { return true;} ///< implies a finite area 01199 virtual bool has_finite_area() const { return true;} ///< regardless of the reference 01200 virtual bool has_known_area() const { return true;} ///< the area is analytically known 01201 virtual double known_area() const { 01202 return twopi * 2 * _delta; 01203 } 01204 01205 protected: 01206 double _delta; 01207 }; 01208 01209 01210 // select on objets within a distance 'radius' of a variable location 01211 Selector SelectorStrip(const double & half_width) { 01212 return Selector(new SW_Strip(half_width)); 01213 } 01214 01215 01216 //---------------------------------------------------------------------- 01217 /// helper for selecting on objects with rapidity within a distance 01218 /// 'delta_rap' of a reference and phi within a distanve delta_phi of 01219 /// a reference 01220 class SW_Rectangle : public SW_WithReference { 01221 public: 01222 SW_Rectangle(const double &delta_rap, const double &delta_phi) 01223 : _delta_rap(delta_rap), _delta_phi(delta_phi) {} 01224 01225 /// return a copy of the current object 01226 virtual SelectorWorker* copy(){ return new SW_Rectangle(*this);} 01227 01228 /// returns true if a given object passes the selection criterium 01229 /// this has to be overloaded by derived workers 01230 virtual bool pass(const PseudoJet & jet) const { 01231 // make sure the centre is initialised 01232 if (! _is_initialised) 01233 throw Error("To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)"); 01234 01235 return (abs(jet.rap()-_reference.rap()) <= _delta_rap) && (abs(jet.delta_phi_to(_reference)) <= _delta_phi); 01236 } 01237 01238 /// returns a description of the worker 01239 virtual string description() const { 01240 ostringstream ostr; 01241 ostr << "|rap - rap_reference| <= " << _delta_rap << " && |phi - phi_reference| <= " << _delta_phi ; 01242 return ostr.str(); 01243 } 01244 01245 /// returns the rapidity range for which it may return "true" 01246 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{ 01247 // make sure the centre is initialised 01248 if (! _is_initialised) 01249 throw Error("To use a SelectorRectangle (or any selector that requires a reference), you first have to call set_reference(...)"); 01250 01251 rapmax = _reference.rap()+_delta_rap; 01252 rapmin = _reference.rap()-_delta_rap; 01253 } 01254 01255 virtual bool is_geometric() const { return true;} ///< implies a finite area 01256 virtual bool has_finite_area() const { return true;} ///< regardless of the reference 01257 virtual bool has_known_area() const { return true;} ///< the area is analytically known 01258 virtual double known_area() const { 01259 return 4 * _delta_rap * _delta_phi; 01260 } 01261 01262 protected: 01263 double _delta_rap, _delta_phi; 01264 }; 01265 01266 01267 // select on objets within a distance 'radius' of a variable location 01268 Selector SelectorRectangle(const double & half_rap_width, const double & half_phi_width) { 01269 return Selector(new SW_Rectangle(half_rap_width, half_phi_width)); 01270 } 01271 01272 01273 //---------------------------------------------------------------------- 01274 /// helper for selecting the jets that carry at least a given fraction 01275 /// of the reference jet 01276 class SW_PtFractionMin : public SW_WithReference { 01277 public: 01278 /// ctor with specification of the number of objects to keep 01279 SW_PtFractionMin(double fraction) : _fraction2(fraction*fraction){} 01280 01281 /// return a copy of the current object 01282 virtual SelectorWorker* copy(){ return new SW_PtFractionMin(*this);} 01283 01284 /// return true if the jet carries a large enough fraction of the reference. 01285 /// Throw an error if the reference is not initialised. 01286 virtual bool pass(const PseudoJet & jet) const { 01287 // make sure the centre is initialised 01288 if (! _is_initialised) 01289 throw Error("To use a SelectorPtFractionMin (or any selector that requires a reference), you first have to call set_reference(...)"); 01290 01291 // otherwise, just call that method on the jet 01292 return (jet.perp2() >= _fraction2*_reference.perp2()); 01293 } 01294 01295 /// returns a description of the worker 01296 virtual string description() const { 01297 ostringstream ostr; 01298 ostr << "pt >= " << sqrt(_fraction2) << "* pt_ref"; 01299 return ostr.str(); 01300 } 01301 01302 protected: 01303 double _fraction2; 01304 }; 01305 01306 01307 // select objects that carry at least a fraction "fraction" of the reference jet 01308 // (Note that this selectir takes a reference) 01309 Selector SelectorPtFractionMin(double fraction){ 01310 return Selector(new SW_PtFractionMin(fraction)); 01311 } 01312 01313 01314 //---------------------------------------------------------------------- 01315 // additional (mostly helper) selectors 01316 //---------------------------------------------------------------------- 01317 01318 //---------------------------------------------------------------------- 01319 /// helper for selecting the 0-momentum jets 01320 class SW_IsZero : public SelectorWorker { 01321 public: 01322 /// ctor 01323 SW_IsZero(){} 01324 01325 /// return true if the jet has zero momentum 01326 virtual bool pass(const PseudoJet & jet) const { 01327 return jet==0; 01328 } 01329 01330 /// rereturns a description of the worker 01331 virtual string description() const { return "zero";} 01332 }; 01333 01334 01335 // select objects with zero momentum 01336 Selector SelectorIsZero(){ 01337 return Selector(new SW_IsZero()); 01338 } 01339 01340 01341 //---------------------------------------------------------------------- 01342 /// helper for selecting the pure ghost 01343 class SW_IsPureGhost : public SelectorWorker { 01344 public: 01345 /// ctor 01346 SW_IsPureGhost(){} 01347 01348 /// return true if the jet is a pure-ghost jet 01349 virtual bool pass(const PseudoJet & jet) const { 01350 // if the jet has no area support then it's certainly not a ghost 01351 if (!jet.has_area()) return false; 01352 01353 // otherwise, just call that method on the jet 01354 return jet.is_pure_ghost(); 01355 } 01356 01357 /// rereturns a description of the worker 01358 virtual string description() const { return "pure ghost";} 01359 }; 01360 01361 01362 // select objects that are (or are only made of) ghosts 01363 Selector SelectorIsPureGhost(){ 01364 return Selector(new SW_IsPureGhost()); 01365 } 01366 01367 01368 //---------------------------------------------------------------------- 01369 // Selector and workers for obtaining a Selector from an old 01370 // RangeDefinition 01371 // 01372 // This is mostly intended for backward compatibility and is likely to 01373 // be removed in a future major release of FastJet 01374 //---------------------------------------------------------------------- 01375 01376 //---------------------------------------------------------------------- 01377 /// helper for selecting on both rapidity and azimuthal angle 01378 class SW_RangeDefinition : public SelectorWorker{ 01379 public: 01380 /// ctor from a RangeDefinition 01381 SW_RangeDefinition(const RangeDefinition &range) : _range(&range){} 01382 01383 /// transfer the selection creterium to the underlying RangeDefinition 01384 virtual bool pass(const PseudoJet & jet) const { 01385 return _range->is_in_range(jet); 01386 } 01387 01388 /// returns a description of the worker 01389 virtual string description() const { 01390 return _range->description(); 01391 } 01392 01393 /// returns the rapidity range for which it may return "true" 01394 virtual void get_rapidity_extent(double & rapmin, double & rapmax) const{ 01395 _range->get_rap_limits(rapmin, rapmax); 01396 } 01397 01398 /// check if it has a finite area 01399 virtual bool is_geometric() const { return true;} 01400 01401 /// check if it has an analytically computable area 01402 virtual bool has_known_area() const { return true;} 01403 01404 /// if it has a computable area, return it 01405 virtual double known_area() const{ 01406 return _range->area(); 01407 } 01408 01409 protected: 01410 const RangeDefinition *_range; 01411 }; 01412 01413 01414 // ctor from a RangeDefinition 01415 //---------------------------------------------------------------------- 01416 // 01417 // This is provided for backward compatibility and will be removed in 01418 // a future major release of FastJet 01419 Selector::Selector(const RangeDefinition &range) { 01420 _worker.reset(new SW_RangeDefinition(range)); 01421 } 01422 01423 01424 // operators applying directly on a Selector 01425 //---------------------------------------------------------------------- 01426 01427 // operator &= 01428 // For 2 Selectors a and b, a &= b is eauivalent to a = a & b; 01429 Selector & Selector::operator &=(const Selector & b){ 01430 _worker.reset(new SW_And(*this, b)); 01431 return *this; 01432 } 01433 01434 // operator &= 01435 // For 2 Selectors a and b, a &= b is eauivalent to a = a & b; 01436 Selector & Selector::operator |=(const Selector & b){ 01437 _worker.reset(new SW_Or(*this, b)); 01438 return *this; 01439 } 01440 01441 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh