FastJet 3.0.2
|
00001 //---------------------------------------------------------------------- 00002 /// \file 00003 /// \page Example08 08 - using the Selector tool 00004 /// 00005 /// fastjet sample program to illustrate the use of fastjet::Selector 00006 /// 00007 /// run it with : ./08-selector < data/single-event.dat 00008 /// 00009 /// Source code: 08-selector.cc 00010 //---------------------------------------------------------------------- 00011 00012 //STARTHEADER 00013 // $Id: 08-selector.cc 2684 2011-11-14 07:41:44Z soyez $ 00014 // 00015 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez 00016 // 00017 //---------------------------------------------------------------------- 00018 // This file is part of FastJet. 00019 // 00020 // FastJet is free software; you can redistribute it and/or modify 00021 // it under the terms of the GNU General Public License as published by 00022 // the Free Software Foundation; either version 2 of the License, or 00023 // (at your option) any later version. 00024 // 00025 // The algorithms that underlie FastJet have required considerable 00026 // development and are described in hep-ph/0512210. If you use 00027 // FastJet as part of work towards a scientific publication, please 00028 // include a citation to the FastJet paper. 00029 // 00030 // FastJet is distributed in the hope that it will be useful, 00031 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00032 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00033 // GNU General Public License for more details. 00034 // 00035 // You should have received a copy of the GNU General Public License 00036 // along with FastJet. If not, see <http://www.gnu.org/licenses/>. 00037 //---------------------------------------------------------------------- 00038 //ENDHEADER 00039 00040 #include "fastjet/ClusterSequence.hh" 00041 #include "fastjet/Selector.hh" 00042 #include <iostream> // needed for io 00043 #include <cstdio> // needed for io 00044 00045 using namespace std; 00046 using namespace fastjet; 00047 00048 int main(){ 00049 00050 // read in input particles 00051 //---------------------------------------------------------- 00052 vector<PseudoJet> input_particles; 00053 00054 double px, py , pz, E; 00055 while (cin >> px >> py >> pz >> E) { 00056 // create a PseudoJet with these components and put it onto 00057 // the back of the input_particles vector 00058 input_particles.push_back(PseudoJet(px,py,pz,E)); 00059 } 00060 00061 // Selector application #1: keep particles within a given acceptance 00062 // e.g. all particles with 1<|y|<2.5, and all particles with pt>1 for |y|<1 00063 // (we include the (redundant) fastjet:: prefix just as a reminder that this 00064 // is the namespace where Selector is to be found). 00065 //---------------------------------------------------------- 00066 fastjet::Selector particle_selector = fastjet::SelectorAbsRapRange(1.0,2.5) 00067 || (fastjet::SelectorAbsRapMax(1.0) && fastjet::SelectorPtMin(1.0)); 00068 cout << input_particles.size() << " particles before selector" << endl; 00069 input_particles = particle_selector(input_particles); 00070 cout << input_particles.size() << " particles after selector" << endl; 00071 00072 // create a jet definition: 00073 // a jet algorithm with a given radius parameter 00074 //---------------------------------------------------------- 00075 double R = 0.6; 00076 JetDefinition jet_def(kt_algorithm, R); 00077 00078 00079 // run the jet clustering with the above jet definition 00080 //---------------------------------------------------------- 00081 ClusterSequence clust_seq(input_particles, jet_def); 00082 00083 00084 // get the 5 hardest jets within |y|<2 00085 // 00086 // Note that this nicely illustrates that you should watch out that 00087 // Selectors do not necessarily commute. 00088 // 00089 // The && operator behaves like a logical and, i.e. it keeps objects 00090 // that satisfy both criteria (independently). It does commute. 00091 // 00092 // The * operator applies Selectors successively (starting from the 00093 // rightmost one as in a usual operator product). Here, order may matter. 00094 //---------------------------------------------------------- 00095 Selector jet_selector = SelectorNHardest(5) * SelectorAbsRapMax(2.0); 00096 vector<PseudoJet> inclusive_jets = sorted_by_pt(jet_selector(clust_seq.inclusive_jets())); 00097 00098 00099 // tell the user what was done 00100 // - the description of the algorithm used 00101 // - the description of teh selectors used 00102 // - the jets 00103 // show the output as 00104 // {index, rap, phi, pt} 00105 //---------------------------------------------------------- 00106 double rapmin, rapmax; 00107 00108 cout << "Ran " << jet_def.description() << endl; 00109 00110 cout << "Selected particles: " << particle_selector.description() << endl; 00111 particle_selector.get_rapidity_extent(rapmin, rapmax); 00112 cout << " with a total rapidity range of [" << rapmin << ", " << rapmax << "]" << endl; 00113 00114 cout << "Selected jets: " << jet_selector.description() << endl; 00115 jet_selector.get_rapidity_extent(rapmin, rapmax); 00116 cout << " with a total rapidity range of [" << rapmin << ", " << rapmax << "]" << endl; 00117 00118 // label the columns 00119 printf("%5s %15s %15s %15s\n","jet #", "rapidity", "phi", "pt"); 00120 00121 // print out the details for each jet 00122 for (unsigned int i = 0; i < inclusive_jets.size(); i++) { 00123 printf("%5u %15.8f %15.8f %15.8f\n", 00124 i, inclusive_jets[i].rap(), inclusive_jets[i].phi(), 00125 inclusive_jets[i].perp()); 00126 } 00127 00128 return 0; 00129 }