#include <iostream>
#include <sstream>
#include <valarray>
#include <vector>
#include <cstddef>
#include "CmdLine.hh"
#include "fastjet/internal/numconsts.hh"
#include "KtJet/KtEvent.h"
#include "KtJet/KtLorentzVector.h"
Include dependency graph for ktjet_timing.cc:
Go to the source code of this file.
Functions | |
double | pow2 (const double x) |
int | main (int argc, char **argv) |
a program to test and time the kt algorithm as implemented in ktjet |
int main | ( | int | argc, | |
char ** | argv | |||
) |
a program to test and time the kt algorithm as implemented in ktjet
Retrieve the final state jets from KtEvent sorted by Pt
Print out jets 4-momentum and Pt
Definition at line 27 of file ktjet_timing.cc.
References CmdLine::double_val(), CmdLine::int_val(), pow2(), CmdLine::present(), and fastjet::twopi.
00027 { 00028 00029 CmdLine cmdline(argc,argv); 00030 //bool clever = !cmdline.present("-dumb"); 00031 int repeat = cmdline.int_val("-repeat",1); 00032 int combine = cmdline.int_val("-combine",1); 00033 bool write = cmdline.present("-write"); 00034 double ktR = cmdline.double_val("-r",1.0); 00035 double inclkt = cmdline.double_val("-incl",-1.0); 00036 int excln = cmdline.int_val ("-excln",-1); 00037 double excld = cmdline.double_val("-excld",-1.0); 00038 int nev = cmdline.int_val("-nev",1); 00039 bool massless = cmdline.present("-massless"); 00040 00041 for (int iev = 0; iev < nev; iev++) { 00042 vector<KtJet::KtLorentzVector> jets; 00043 string line; 00044 int ndone = 0; 00045 while (getline(cin, line)) { 00046 //cout << line<<endl; 00047 istringstream linestream(line); 00048 if (line == "#END") { 00049 ndone += 1; 00050 if (ndone == combine) {break;} 00051 } 00052 if (line.substr(0,1) == "#") {continue;} 00053 valarray<double> fourvec(4); 00054 linestream >> fourvec[0] >> fourvec[1] >> fourvec[2] >> fourvec[3]; 00055 if (massless) { 00056 linestream >> fourvec[0] >> fourvec[1] >> fourvec[2]; 00057 fourvec[3] = sqrt(pow2(fourvec[0])+pow2(fourvec[1])+pow2(fourvec[2]));} 00058 else { 00059 linestream >> fourvec[0] >> fourvec[1] >> fourvec[2] >> fourvec[3]; 00060 } 00061 KtJet::KtLorentzVector p(fourvec[0],fourvec[1],fourvec[2],fourvec[3]); 00062 jets.push_back(p); 00063 } 00064 00065 00066 // set KtEvent flags 00067 int type = 4; // PP 00068 int angle = 2; // delta R 00069 int recom = 1; // E 00070 //double rparameter = 1.0; 00071 00072 for (int i = 0; i < repeat ; i++) { 00073 // Construct the KtEvent object 00074 KtJet::KtEvent ev(jets,type,angle,recom,ktR); 00075 00076 if (i!=0) {continue;} 00077 cout << "Number of particles = "<< jets.size() << endl; 00078 cout << "Algorithm: KtJet (long.inv.) with R = " << ktR << endl; 00079 00080 // Print out the number of final state jets 00081 //std::cout << "Number of final state jets: " << ev.getNJets() << std::endl; 00082 if (write) { 00084 std::vector<KtJet::KtLorentzVector> jets = ev.getJetsPt(); 00085 00087 std::vector<KtJet::KtLorentzVector>::const_iterator itr = jets.begin(); 00088 for( ; itr != jets.end() ; ++itr) { 00089 std::cout << "Jets Pt2: " << pow2((*itr).perp()) << std::endl; 00090 } 00091 } 00092 00093 if (inclkt >= 0.0) { 00094 // Retrieve the final state jets from KtEvent sorted by Pt 00095 std::vector<KtJet::KtLorentzVector> jets = ev.getJetsPt(); 00096 00097 // Print out index, rap, phi, pt 00098 for (size_t j = 0; j < jets.size(); j++) { 00099 if (jets[j].perp() < inclkt) {break;} 00100 double phi = jets[j].phi(); 00101 if (phi < 0.0) {phi += fastjet::twopi;} 00102 printf("%5u %15.8f %15.8f %15.8f\n",j,jets[j].rapidity(),phi,jets[j].perp()); 00103 } 00104 } 00105 00106 if (excln > 0) { 00107 ev.findJetsN(excln); 00108 vector<KtJet::KtLorentzVector> jets = ev.getJetsPt(); 00109 cout << "Printing "<<excln<<" exclusive jets\n"; 00110 for (size_t j = 0; j < jets.size(); j++) { 00111 printf("%5u %15.8f %15.8f %15.8f\n",j, 00112 jets[j].rapidity(),jets[j].phi(),jets[j].perp()); 00113 } 00114 } 00115 00116 if (excld > 0.0) { 00117 ev.findJetsD(excld); 00118 vector<KtJet::KtLorentzVector> jets = ev.getJetsPt(); 00119 cout << "Printing exclusive jets for d = "<<excld<<"\n"; 00120 for (size_t j = 0; j < jets.size(); j++) { 00121 printf("%5u %15.8f %15.8f %15.8f\n",j, 00122 jets[j].rapidity(),jets[j].phi(),jets[j].perp()); 00123 } 00124 } 00125 00126 00127 } 00128 00129 } 00130 }
double pow2 | ( | const double | x | ) | [inline] |