00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00026
00027 #include <stdio.h>
00028 #include <iostream>
00029 #include "siscone/momentum.h"
00030 #include "siscone/siscone.h"
00031 #include "options.h"
00032
00033 using namespace std;
00034 using namespace siscone;
00035
00036 int main(int argc, char *argv[]){
00037 vector<Cmomentum> particles;
00038 Csiscone siscone;
00039 int i,N;
00040 double px,py,pz,E;
00041 Coptions opts;
00042 char fline[512];
00043
00044 if (opts.parse_options(argc, argv))
00045 exit(1);
00046
00047
00048 if (opts.help_flag){
00049 opts.print_help();
00050 exit(0);
00051 }
00052
00053
00054 if (opts.version_flag){
00055 opts.print_version();
00056 exit(0);
00057 }
00058
00059
00060 FILE *flux;
00061 FILE *fpart;
00062
00063
00064 if (opts.verbose_flag) cout << "reading particles" << endl;
00065 flux = fopen(opts.ev_name, "r");
00066 if (flux==NULL){
00067 cerr << "cannot read event '" << opts.ev_name << "'" << endl;
00068 cerr << "specify the event to read using the -e option" << endl;
00069 return 1;
00070 }
00071
00072 N=0;
00073 fpart = fopen("particles.dat", "w+");
00074 while ((opts.N_stop!=0) && (fgets(fline, 512, flux)!=NULL)){
00075 if (fline[0]!='#'){
00076 if (sscanf(fline, "%le%le%le%le", &px, &py, &pz, &E)==4){
00077 particles.push_back(Cmomentum(px, py, pz, E));
00078 fprintf(fpart, "%le\t%le\n", particles[N].eta, particles[N].phi);
00079 N++;
00080 opts.N_stop--;
00081 } else {
00082 cout << "error in reading event file Giving up." << endl;
00083 fclose(flux);
00084 fclose(fpart);
00085 exit(2);
00086 }
00087 }
00088 }
00089 fclose(flux);
00090 fclose(fpart);
00091 if (opts.verbose_flag)
00092 cout << " working with " << N << " particles" << endl;
00093
00094
00095 if (opts.verbose_flag) cout << "computing jet contents" << endl;
00096 i=siscone.compute_jets(particles, opts.R, opts.f, opts.npass, opts.ptmin, opts.SM_var);
00097 if (opts.verbose_flag){
00098 unsigned int pass;
00099 for (pass=0;pass<siscone.protocones_list.size();pass++)
00100 cout << " pass " << pass << " found " << siscone.protocones_list[pass].size()
00101 << " stable cones" << endl;
00102 cout << " Final result: " << i << " jets found" << endl;
00103 }
00104
00105
00106 if (opts.verbose_flag)
00107 cout << "saving result" << endl;
00108 flux = fopen("jets.dat", "w+");
00109 siscone.save_contents(flux);
00110 fclose(flux);
00111
00112 if (opts.verbose_flag)
00113 cout << "bye..." << endl;
00114
00115 return 0;
00116 }