Fork me on GitHub

source: svn/trunk/src/SmearUtil.cc@ 523

Last change on this file since 523 was 494, checked in by Xavier Rouby, 15 years ago

new: calorimeter endcaps

File size: 93.2 KB
Line 
1/***********************************************************************
2** **
3** /----------------------------------------------\ **
4** | Delphes, a framework for the fast simulation | **
5** | of a generic collider experiment | **
6** \------------- arXiv:0903.2225v1 ------------/ **
7** **
8** **
9** This package uses: **
10** ------------------ **
11** ROOT: Nucl. Inst. & Meth. in Phys. Res. A389 (1997) 81-86 **
12** FastJet algorithm: Phys. Lett. B641 (2006) [hep-ph/0512210] **
13** Hector: JINST 2:P09005 (2007) [physics.acc-ph:0707.1198v2] **
14** FROG: [hep-ex/0901.2718v1] **
15** HepMC: Comput. Phys. Commun.134 (2001) 41 **
16** **
17** ------------------------------------------------------------------ **
18** **
19** Main authors: **
20** ------------- **
21** **
22** Severine Ovyn Xavier Rouby **
23** severine.ovyn@uclouvain.be xavier.rouby@cern **
24** **
25** Center for Particle Physics and Phenomenology (CP3) **
26** Universite catholique de Louvain (UCL) **
27** Louvain-la-Neuve, Belgium **
28** **
29** Copyright (C) 2008-2009, **
30** All rights reserved. **
31** **
32***********************************************************************/
33
34/// \file SmearUtil.cc
35/// \brief RESOLution class, and some generic definitions
36
37
38#include "SmearUtil.h"
39#include "TRandom.h"
40#include "TStopwatch.h"
41
42#include <iostream>
43#include <fstream>
44#include <sstream>
45#include <iomanip>
46#include <map>
47#include <vector>
48using namespace std;
49
50//------------------------------------------------------------------------------
51
52RESOLution::RESOLution() {
53
54 // Detector characteristics
55 CEN_max_tracker = 2.5; // Maximum tracker coverage
56 CEN_max_calo_cen = 1.7; // central calorimeter coverage
57 CEN_max_calo_ec = 3.0; // calorimeter endcap coverage
58 CEN_max_calo_fwd = 5.0; // forward calorimeter pseudorapidity coverage
59 CEN_max_mu = 2.4; // muon chambers pseudorapidity coverage
60
61 // Energy resolution for electron/photon
62 // \sigma/E = C + N/E + S/\sqrt{E}
63 ELG_Scen = 0.05; // S term for central ECAL
64 ELG_Ncen = 0.25; // N term
65 ELG_Ccen = 0.005; // C term
66 ELG_Sec = 0.05; // S term for central ECAL endcap
67 ELG_Nec = 0.25; // S term
68 ELG_Cec = 0.005; // S term
69 ELG_Sfwd = 2.084; // S term for FCAL
70 ELG_Nfwd = 0.0; // N term
71 ELG_Cfwd = 0.107; // C term
72 ELG_Szdc = 0.70; // S term for ZDC
73 ELG_Nzdc = 0.0; // N term
74 ELG_Czdc = 0.08; // C term
75
76 // Energy resolution for hadrons in ecal/hcal/fwd
77 // \sigma/E = C + N/E + S/\sqrt{E}
78 HAD_Scen = 1.5; // S term for central HCAL
79 HAD_Ncen = 0.; // N term
80 HAD_Ccen = 0.05; // C term
81 HAD_Sec = 1.5; // S term for HCAL endcap
82 HAD_Nec = 0.; // N term
83 HAD_Cec = 0.05; // C term
84 HAD_Sfwd = 2.7; // S term for FCAL
85 HAD_Nfwd = 0.; // N term
86 HAD_Cfwd = 0.13; // C term
87 HAD_Szdc = 1.38; // S term for ZDC
88 HAD_Nzdc = 0.; // N term
89 HAD_Czdc = 0.13; // C term
90
91 // Muon smearing
92 MU_SmearPt = 0.01;
93
94 // time resolution
95 ZDC_T_resolution = 0; // resolution for time measurement [s]
96 RP220_T_resolution = 0;
97 RP420_T_resolution = 0;
98
99 // Tracking efficiencies
100 TRACK_ptmin = 0.9; // minimal pt needed to reach the calorimeter in GeV
101 TRACK_eff = 100; // efficiency associated to the tracking
102
103 // Calorimetric towers
104 TOWER_number = 40;
105 const float tower_eta_edges[41] = {
106 0., 0.087, 0.174, 0.261, 0.348, 0.435, 0.522, 0.609, 0.696, 0.783, 0.870, 0.957, 1.044, 1.131, 1.218, 1.305, 1.392, 1.479, 1.566,
107 1.653, 1.740, 1.830, 1.930, 2.043, 2.172, 2.322, 2.500, 2.650, 2.868, 2.950, 3.125, 3.300, 3.475, 3.650, 3.825, 4.000, 4.175,
108 4.350, 4.525, 4.700, 5.000}; // temporary object
109 TOWER_eta_edges = new float[TOWER_number+1];
110 for(unsigned int i=0; i<TOWER_number +1; i++) TOWER_eta_edges[i] = tower_eta_edges[i];
111
112 const float tower_dphi[40] = {
113 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 10,
114 10,10,10,10,10, 10,10,10,10,10, 10,10,10,10,10, 10,10,10,20, 20 }; // temporary object
115 TOWER_dphi = new float[TOWER_number];
116 for(unsigned int i=0; i<TOWER_number; i++) TOWER_dphi[i] = tower_dphi[i];
117
118
119 // Thresholds for reconstructed objetcs (GeV)
120 PTCUT_elec = 10.0;
121 PTCUT_muon = 10.0;
122 PTCUT_jet = 20.0;
123 PTCUT_gamma = 10.0;
124 PTCUT_taujet = 10.0;
125
126 ZDC_gamma_E = 20; // GeV
127 ZDC_n_E = 50; // GeV
128
129 // Isolation
130 ISOL_PT = 2.0; //minimal pt of tracks for isolation criteria
131 ISOL_Cone = 0.5; //Cone for isolation criteria
132 ISOL_Calo_ET = 1E99; //minimal tower energy for isolation criteria. Default off = 1E99
133 ISOL_Calo_Grid = 3; //Grid size (N x N) for calorimetric isolation -- should be odd
134
135 // General jet variable
136 JET_coneradius = 0.7; // generic jet radius ; not for tau's !!!
137 JET_jetalgo = 1; // 1 for Cone algorithm, 2 for MidPoint algorithm, 3 for SIScone algorithm, 4 for kt algorithm
138 JET_seed = 1.0; // minimum seed to start jet reconstruction
139 JET_Eflow = 1; // 1 for Energy flow in jets reco ; 0 if not
140
141 // Tagging definition
142 BTAG_b = 40;
143 BTAG_mistag_c = 10;
144 BTAG_mistag_l = 1;
145
146 // FLAGS
147 FLAG_bfield = 1; //1 to run the bfield propagation else 0
148 FLAG_vfd = 1; //1 to run the very forward detectors else 0
149 FLAG_RP = 1; //1 to run the zero degree calorimeter else 0
150 FLAG_trigger = 1; //1 to run the trigger selection else 0
151 FLAG_frog = 1; //1 to run the FROG event display
152 FLAG_lhco = 1;
153
154 // In case BField propagation allowed
155 TRACK_radius = 129; //radius of the BField coverage
156 TRACK_length = 300; //length of the BField coverage
157 TRACK_bfield_x = 0; //X composant of the BField
158 TRACK_bfield_y = 0; //Y composant of the BField
159 TRACK_bfield_z = 3.8; //Z composant of the BField
160
161 // In case Very forward detectors allowed
162 VFD_min_calo_vfd = 5.2; // very forward calorimeter (if any) like CASTOR
163 VFD_max_calo_vfd = 6.6;
164 VFD_min_zdc = 8.3;
165 VFD_s_zdc = 140; // distance of the Zero Degree Calorimeter, from the Interaction poin, in [m]
166
167 RP_220_s = 220; // distance of the RP to the IP, in meters
168 RP_220_x = 0.002; // distance of the RP to the beam, in meters
169 RP_420_s = 420; // distance of the RP to the IP, in meters
170 RP_420_x = 0.004; // distance of the RP to the beam, in meters
171 RP_IP_name = "IP5";
172 RP_beam1Card = "data/LHCB1IR5_v6.500.tfs";
173 RP_beam2Card = "data/LHCB1IR5_v6.500.tfs";
174
175 // In case FROG event display allowed
176 NEvents_Frog = 10;
177
178 // Number of events to be processed
179 NEvents = -1;
180
181 //********************************************
182 //jet stuffs not defined in the input datacard
183 //********************************************
184
185 JET_overlap = 0.75;
186 // MidPoint algorithm definition
187 JET_M_coneareafraction = 0.25;
188 JET_M_maxpairsize = 2;
189 JET_M_maxiterations = 100;
190 // Define Cone algorithm.
191 JET_C_adjacencycut = 2;
192 JET_C_maxiterations = 100;
193 JET_C_iratch = 1;
194 //Define SISCone algorithm.
195 JET_S_npass = 0;
196 JET_S_protojet_ptmin= 0.0;
197
198 //For Tau-jet definition
199 TAU_energy_scone = 0.15; // radius R of the cone for tau definition, based on energy threshold
200 TAU_track_scone = 0.4; // radius R of the cone for tau definition, based on track number
201 TAU_track_pt = 2; // minimal pt [GeV] for tracks to be considered in tau definition
202 TAU_energy_frac = 0.95; // fraction of energy required in the central part of the cone, for tau jets
203
204 PT_QUARKS_MIN = 2.0 ; // minimal pt needed by quarks to do b-tag
205
206 //for very forward detectors
207 RP_offsetEl_s = 120; // distance of beam separation point, from IP
208 RP_offsetEl_x = -0.097; // half distance of separation of beams in horizontal plan, in m
209 RP_offsetEl_y = 0; // half distance of separation of beams in vertical plan, in m
210 RP_cross_x = -500; // IP offset in horizontal plane, in micrometers
211 RP_cross_y = 0.0; // IP offset in vertical plane, in micrometers
212 RP_cross_ang_x = 142.5; // half-crossing angle in horizontal plane, in microrad
213 RP_cross_ang_y = 0.0; // half-crossing angle in vertical plane, in microrad
214
215
216 PdgTableFilename = "data/particle.tbl";
217 inputfilelist = "";
218 detectorcard = "";
219 triggercard = "";
220}
221
222
223RESOLution::RESOLution(const RESOLution & DET) {
224 // Detector characteristics
225 CEN_max_tracker = DET.CEN_max_tracker;
226 CEN_max_calo_cen = DET.CEN_max_calo_cen;
227 CEN_max_calo_ec = DET.CEN_max_calo_ec;
228 CEN_max_calo_fwd = DET.CEN_max_calo_fwd;
229 CEN_max_mu = DET.CEN_max_mu;
230
231 // Energy resolution for electron/photon
232 ELG_Scen = DET.ELG_Scen;
233 ELG_Ncen = DET.ELG_Ncen;
234 ELG_Ccen = DET.ELG_Ccen;
235 ELG_Sec = DET.ELG_Sec;
236 ELG_Nec = DET.ELG_Nec;
237 ELG_Cec = DET.ELG_Cec;
238 ELG_Cfwd = DET.ELG_Cfwd;
239 ELG_Sfwd = DET.ELG_Sfwd;
240 ELG_Nfwd = DET.ELG_Nfwd;
241 ELG_Czdc = DET.ELG_Czdc;
242 ELG_Szdc = DET.ELG_Szdc;
243 ELG_Nzdc = DET.ELG_Nzdc;
244
245 // Energy resolution for hadrons in ecal/hcal/fwd/zdc
246 HAD_Scen = DET.HAD_Scen;
247 HAD_Ncen = DET.HAD_Ncen;
248 HAD_Ccen = DET.HAD_Ccen;
249 HAD_Sec = DET.HAD_Sec;
250 HAD_Nec = DET.HAD_Nec;
251 HAD_Cec = DET.HAD_Cec;
252 HAD_Sfwd = DET.HAD_Sfwd;
253 HAD_Nfwd = DET.HAD_Nfwd;
254 HAD_Cfwd = DET.HAD_Cfwd;
255 HAD_Szdc = DET.HAD_Szdc;
256 HAD_Nzdc = DET.HAD_Nzdc;
257 HAD_Czdc = DET.HAD_Czdc;
258
259 // time resolution
260 ZDC_T_resolution = DET.ZDC_T_resolution; // resolution for time measurement [s]
261 RP220_T_resolution = DET.RP220_T_resolution;
262 RP420_T_resolution = DET.RP420_T_resolution;
263
264 // Muon smearing
265 MU_SmearPt = DET.MU_SmearPt;
266
267 // Tracking efficiencies
268 TRACK_ptmin = DET.TRACK_ptmin;
269 TRACK_eff = DET.TRACK_eff;
270
271 // Calorimetric towers
272 TOWER_number = DET.TOWER_number;
273 TOWER_eta_edges = new float[TOWER_number+1];
274 for(unsigned int i=0; i<TOWER_number +1; i++) TOWER_eta_edges[i] = DET.TOWER_eta_edges[i];
275
276 TOWER_dphi = new float[TOWER_number];
277 for(unsigned int i=0; i<TOWER_number; i++) TOWER_dphi[i] = DET.TOWER_dphi[i];
278
279 // Thresholds for reconstructed objetcs
280 PTCUT_elec = DET.PTCUT_elec;
281 PTCUT_muon = DET.PTCUT_muon;
282 PTCUT_jet = DET.PTCUT_jet;
283 PTCUT_gamma = DET.PTCUT_gamma;
284 PTCUT_taujet = DET.PTCUT_taujet;
285
286 ZDC_gamma_E = DET.ZDC_gamma_E;
287 ZDC_n_E = DET.ZDC_n_E;
288
289 // Isolation
290 ISOL_PT = DET.ISOL_PT; // tracking isolation
291 ISOL_Cone = DET.ISOL_Cone;
292 ISOL_Calo_ET = DET.ISOL_Calo_ET; // calorimeter isolation, defaut off
293 ISOL_Calo_Grid = DET.ISOL_Calo_Grid;
294
295
296 // General jet variable
297 JET_coneradius = DET.JET_coneradius;
298 JET_jetalgo = DET.JET_jetalgo;
299 JET_seed = DET.JET_seed;
300 JET_Eflow = DET.JET_Eflow;
301
302 // Tagging definition
303 BTAG_b = DET.BTAG_b;
304 BTAG_mistag_c = DET.BTAG_mistag_c;
305 BTAG_mistag_l = DET.BTAG_mistag_l;
306
307 // FLAGS
308 FLAG_bfield = DET.FLAG_bfield;
309 FLAG_vfd = DET.FLAG_vfd;
310 FLAG_RP = DET.FLAG_RP;
311 FLAG_trigger = DET.FLAG_trigger;
312 FLAG_frog = DET.FLAG_frog;
313 FLAG_lhco = DET.FLAG_lhco;
314
315 // In case BField propagation allowed
316 TRACK_radius = DET.TRACK_radius;
317 TRACK_length = DET.TRACK_length;
318 TRACK_bfield_x = DET.TRACK_bfield_x;
319 TRACK_bfield_y = DET.TRACK_bfield_y;
320 TRACK_bfield_z = DET.TRACK_bfield_z;
321
322 // In case Very forward detectors allowed
323 VFD_min_calo_vfd = DET.VFD_min_calo_vfd;
324 VFD_max_calo_vfd = DET.VFD_max_calo_vfd;
325 VFD_min_zdc = DET.VFD_min_zdc;
326 VFD_s_zdc = DET.VFD_s_zdc;
327
328 RP_220_s = DET.RP_220_s;
329 RP_220_x = DET.RP_220_x;
330 RP_420_s = DET.RP_420_s;
331 RP_420_x = DET.RP_420_x;
332 RP_beam1Card = DET.RP_beam1Card;
333 RP_beam2Card = DET.RP_beam2Card;
334 RP_offsetEl_s = DET.RP_offsetEl_s;
335 RP_offsetEl_x = DET.RP_offsetEl_x;
336 RP_offsetEl_y = DET.RP_offsetEl_y;
337 RP_cross_x = DET.RP_cross_x;
338 RP_cross_y = DET.RP_cross_y;
339 RP_cross_ang_x = DET.RP_cross_ang_x;
340 RP_cross_ang_y = DET.RP_cross_ang_y;
341 RP_IP_name = DET.RP_IP_name;
342
343 // In case FROG event display allowed
344 NEvents_Frog = DET.NEvents_Frog;
345
346 // Number of events to be processed
347 NEvents = DET.NEvents;
348
349 JET_overlap = DET.JET_overlap;
350 // MidPoint algorithm definition
351 JET_M_coneareafraction = DET.JET_M_coneareafraction;
352 JET_M_maxpairsize = DET.JET_M_maxpairsize;
353 JET_M_maxiterations = DET.JET_M_maxiterations;
354 // Define Cone algorithm.
355 JET_C_adjacencycut = DET.JET_C_adjacencycut;
356 JET_C_maxiterations = DET.JET_C_maxiterations;
357 JET_C_iratch = DET.JET_C_iratch;
358 //Define SISCone algorithm.
359 JET_S_npass = DET.JET_S_npass;
360 JET_S_protojet_ptmin = DET.JET_S_protojet_ptmin;
361
362 //For Tau-jet definition
363 TAU_energy_scone = DET.TAU_energy_scone;
364 TAU_track_scone = DET.TAU_track_scone;
365 TAU_track_pt = DET.TAU_track_pt;
366 TAU_energy_frac = DET.TAU_energy_frac;
367
368 PT_QUARKS_MIN = DET.PT_QUARKS_MIN;
369 PdgTableFilename = DET.PdgTableFilename;
370 PDGtable = DET.PDGtable;
371 inputfilelist = DET.inputfilelist;
372 detectorcard = DET.detectorcard;
373 triggercard = DET.triggercard;
374}
375
376RESOLution& RESOLution::operator=(const RESOLution& DET) {
377 if(this==&DET) return *this;
378 // Detector characteristics
379 CEN_max_tracker = DET.CEN_max_tracker;
380 CEN_max_calo_cen = DET.CEN_max_calo_cen;
381 CEN_max_calo_ec = DET.CEN_max_calo_ec;
382 CEN_max_calo_fwd = DET.CEN_max_calo_fwd;
383 CEN_max_mu = DET.CEN_max_mu;
384
385 // Energy resolution for electron/photon
386 ELG_Scen = DET.ELG_Scen;
387 ELG_Ncen = DET.ELG_Ncen;
388 ELG_Ccen = DET.ELG_Ccen;
389 ELG_Sec = DET.ELG_Sec;
390 ELG_Nec = DET.ELG_Nec;
391 ELG_Cec = DET.ELG_Cec;
392 ELG_Cfwd = DET.ELG_Cfwd;
393 ELG_Sfwd = DET.ELG_Sfwd;
394 ELG_Nfwd = DET.ELG_Nfwd;
395 ELG_Czdc = DET.ELG_Czdc;
396 ELG_Szdc = DET.ELG_Szdc;
397 ELG_Nzdc = DET.ELG_Nzdc;
398
399 // Energy resolution for hadrons in ecal/hcal/fwd/zdc
400 HAD_Scen = DET.HAD_Scen ;
401 HAD_Ncen = DET.HAD_Ncen;
402 HAD_Ccen = DET.HAD_Ccen;
403 HAD_Sec = DET.HAD_Sec;
404 HAD_Nec = DET.HAD_Nec;
405 HAD_Cec = DET.HAD_Cec;
406 HAD_Sfwd = DET.HAD_Sfwd;
407 HAD_Nfwd = DET.HAD_Nfwd;
408 HAD_Cfwd = DET.HAD_Cfwd;
409 HAD_Szdc = DET.HAD_Szdc;
410 HAD_Nzdc = DET.HAD_Nzdc;
411 HAD_Czdc = DET.HAD_Czdc;
412
413 // time resolution
414 ZDC_T_resolution = DET.ZDC_T_resolution; // resolution for time measurement [s]
415 RP220_T_resolution = DET.RP220_T_resolution;
416 RP420_T_resolution = DET.RP420_T_resolution;
417
418 // Muon smearing
419 MU_SmearPt = DET.MU_SmearPt;
420
421 // Tracking efficiencies
422 TRACK_ptmin = DET.TRACK_ptmin;
423 TRACK_eff = DET.TRACK_eff;
424
425 // Calorimetric towers
426 TOWER_number = DET.TOWER_number;
427 TOWER_eta_edges = new float[TOWER_number+1];
428 for(unsigned int i=0; i<TOWER_number +1; i++) TOWER_eta_edges[i] = DET.TOWER_eta_edges[i];
429
430 TOWER_dphi = new float[TOWER_number];
431 for(unsigned int i=0; i<TOWER_number; i++) TOWER_dphi[i] = DET.TOWER_dphi[i];
432
433 // Thresholds for reconstructed objetcs
434 PTCUT_elec = DET.PTCUT_elec;
435 PTCUT_muon = DET.PTCUT_muon;
436 PTCUT_jet = DET.PTCUT_jet;
437 PTCUT_gamma = DET.PTCUT_gamma;
438 PTCUT_taujet = DET.PTCUT_taujet;
439
440 ZDC_gamma_E = DET.ZDC_gamma_E;
441 ZDC_n_E = DET.ZDC_n_E;
442
443 // Isolation
444 ISOL_PT = DET.ISOL_PT; // tracking isolation
445 ISOL_Cone = DET.ISOL_Cone;
446 ISOL_Calo_ET = DET.ISOL_Calo_ET; // calorimeter isolation, defaut off
447 ISOL_Calo_Grid = DET.ISOL_Calo_Grid;
448
449 // General jet variable
450 JET_coneradius = DET.JET_coneradius;
451 JET_jetalgo = DET.JET_jetalgo;
452 JET_seed = DET.JET_seed;
453 JET_Eflow = DET.JET_Eflow;
454
455 // Tagging definition
456 BTAG_b = DET.BTAG_b;
457 BTAG_mistag_c = DET.BTAG_mistag_c;
458 BTAG_mistag_l = DET.BTAG_mistag_l;
459
460 // FLAGS
461 FLAG_bfield = DET.FLAG_bfield;
462 FLAG_vfd = DET.FLAG_vfd;
463 FLAG_RP = DET.FLAG_RP;
464 FLAG_trigger = DET.FLAG_trigger;
465 FLAG_frog = DET.FLAG_frog;
466 FLAG_lhco = DET.FLAG_lhco;
467
468 // In case BField propagation allowed
469 TRACK_radius = DET.TRACK_radius;
470 TRACK_length = DET.TRACK_length;
471 TRACK_bfield_x = DET.TRACK_bfield_x;
472 TRACK_bfield_y = DET.TRACK_bfield_y;
473 TRACK_bfield_z = DET.TRACK_bfield_z;
474
475 // In case Very forward detectors allowed
476 VFD_min_calo_vfd = DET.VFD_min_calo_vfd;
477 VFD_max_calo_vfd = DET.VFD_max_calo_vfd;
478 VFD_min_zdc = DET.VFD_min_zdc;
479 VFD_s_zdc = DET.VFD_s_zdc;
480
481 RP_220_s = DET.RP_220_s;
482 RP_220_x = DET.RP_220_x;
483 RP_420_s = DET.RP_420_s;
484 RP_420_x = DET.RP_420_x;
485 RP_offsetEl_s = DET.RP_offsetEl_s;
486 RP_offsetEl_x = DET.RP_offsetEl_x;
487 RP_offsetEl_y = DET.RP_offsetEl_y;
488 RP_beam1Card = DET.RP_beam1Card;
489 RP_beam2Card = DET.RP_beam2Card;
490 RP_cross_x = DET.RP_cross_x;
491 RP_cross_y = DET.RP_cross_y;
492 RP_cross_ang_x = DET.RP_cross_ang_x;
493 RP_cross_ang_y = DET.RP_cross_ang_y;
494 RP_IP_name = DET.RP_IP_name;
495
496
497 // In case FROG event display allowed
498 NEvents_Frog = DET.NEvents_Frog;
499
500 // Number of events to be processed
501 NEvents = DET.NEvents;
502
503 JET_overlap = DET.JET_overlap;
504 // MidPoint algorithm definition
505 JET_M_coneareafraction = DET.JET_M_coneareafraction;
506 JET_M_maxpairsize = DET.JET_M_maxpairsize;
507 JET_M_maxiterations = DET.JET_M_maxiterations;
508 // Define Cone algorithm.
509 JET_C_adjacencycut = DET.JET_C_adjacencycut;
510 JET_C_maxiterations = DET.JET_C_maxiterations;
511 JET_C_iratch = DET.JET_C_iratch;
512 //Define SISCone algorithm.
513 JET_S_npass = DET.JET_S_npass;
514 JET_S_protojet_ptmin = DET.JET_S_protojet_ptmin;
515
516 //For Tau-jet definition
517 TAU_energy_scone = DET.TAU_energy_scone;
518 TAU_track_scone = DET.TAU_track_scone;
519 TAU_track_pt = DET.TAU_track_pt;
520 TAU_energy_frac = DET.TAU_energy_frac;
521
522 PT_QUARKS_MIN = DET.PT_QUARKS_MIN;
523
524 PdgTableFilename = DET.PdgTableFilename;
525 PDGtable = DET.PDGtable;
526
527 inputfilelist = DET.inputfilelist;
528 detectorcard = DET.detectorcard;
529 triggercard = DET.triggercard;
530
531 return *this;
532}
533
534void RESOLution::setNames(const string& list, const string& det, const string& trig) {
535 inputfilelist = list;
536 detectorcard = det;
537 triggercard = trig;
538}
539
540//------------------------------------------------------------------------------
541void RESOLution::ReadDataCard(const string datacard) {
542
543 string temp_string;
544 istringstream curstring;
545
546 ifstream fichier_a_lire(datacard.c_str());
547 if(!fichier_a_lire.good()) {
548 cout <<"** WARNING: Datadard not found, use default values **" << endl;
549 return;
550 }
551 bool CEN_max_calo_ec_flag = false;
552
553 while (getline(fichier_a_lire,temp_string)) {
554 curstring.clear(); // needed when using several times istringstream::str(string)
555 curstring.str(temp_string);
556 string varname;
557 float value; int ivalue; string svalue;
558
559 if(strstr(temp_string.c_str(),"#")) { }
560 else if(strstr(temp_string.c_str(),"CEN_max_tracker")) {curstring >> varname >> value; CEN_max_tracker = value;}
561 else if(strstr(temp_string.c_str(),"CEN_max_calo_cen")) {curstring >> varname >> value; CEN_max_calo_cen = value;}
562 else if(strstr(temp_string.c_str(),"CEN_max_calo_ec")) {CEN_max_calo_ec_flag=true; curstring >> varname >> value; CEN_max_calo_ec = value;}
563 else if(strstr(temp_string.c_str(),"CEN_max_calo_fwd")) {curstring >> varname >> value; CEN_max_calo_fwd = value;}
564 else if(strstr(temp_string.c_str(),"CEN_max_mu")) {curstring >> varname >> value; CEN_max_mu = value;}
565
566 else if(strstr(temp_string.c_str(),"VFD_min_calo_vfd")) {curstring >> varname >> value; VFD_min_calo_vfd = value;}
567 else if(strstr(temp_string.c_str(),"VFD_max_calo_vfd")) {curstring >> varname >> value; VFD_max_calo_vfd = value;}
568 else if(strstr(temp_string.c_str(),"VFD_min_zdc")) {curstring >> varname >> value; VFD_min_zdc = value;}
569 else if(strstr(temp_string.c_str(),"VFD_s_zdc")) {curstring >> varname >> value; VFD_s_zdc = value;}
570
571 else if(strstr(temp_string.c_str(),"RP_220_s")) {curstring >> varname >> value; RP_220_s = value;}
572 else if(strstr(temp_string.c_str(),"RP_220_x")) {curstring >> varname >> value; RP_220_x = value;}
573 else if(strstr(temp_string.c_str(),"RP_420_s")) {curstring >> varname >> value; RP_420_s = value;}
574 else if(strstr(temp_string.c_str(),"RP_420_x")) {curstring >> varname >> value; RP_420_x = value;}
575 else if(strstr(temp_string.c_str(),"RP_beam1Card")) {curstring >> varname >> svalue;RP_beam1Card = svalue;}
576 else if(strstr(temp_string.c_str(),"RP_beam2Card")) {curstring >> varname >> svalue;RP_beam2Card = svalue;}
577 else if(strstr(temp_string.c_str(),"RP_IP_name")) {curstring >> varname >> svalue;RP_IP_name = svalue;}
578
579 else if(strstr(temp_string.c_str(),"RP_offsetEl_s")) {curstring >> varname >> value; RP_offsetEl_s = value;}
580 else if(strstr(temp_string.c_str(),"RP_offsetEl_x")) {curstring >> varname >> value; RP_offsetEl_x = value;}
581 else if(strstr(temp_string.c_str(),"RP_offsetEl_y")) {curstring >> varname >> value; RP_offsetEl_y = value;}
582 else if(strstr(temp_string.c_str(),"RP_cross_x")) {curstring >> varname >> value; RP_cross_x = value;}
583 else if(strstr(temp_string.c_str(),"RP_cross_y")) {curstring >> varname >> value; RP_cross_y = value;}
584 else if(strstr(temp_string.c_str(),"RP_cross_ang_x")) {curstring >> varname >> value; RP_cross_ang_x = value;}
585 else if(strstr(temp_string.c_str(),"RP_cross_ang_y")) {curstring >> varname >> value; RP_cross_ang_y = value;}
586
587 else if(strstr(temp_string.c_str(),"ELG_Scen")) {curstring >> varname >> value; ELG_Scen = value;}
588 else if(strstr(temp_string.c_str(),"ELG_Ncen")) {curstring >> varname >> value; ELG_Ncen = value;}
589 else if(strstr(temp_string.c_str(),"ELG_Ccen")) {curstring >> varname >> value; ELG_Ccen = value;}
590 else if(strstr(temp_string.c_str(),"ELG_Sec")) {curstring >> varname >> value; ELG_Sec = value;}
591 else if(strstr(temp_string.c_str(),"ELG_Nec")) {curstring >> varname >> value; ELG_Nec = value;}
592 else if(strstr(temp_string.c_str(),"ELG_Cec")) {curstring >> varname >> value; ELG_Cec = value;}
593 else if(strstr(temp_string.c_str(),"ELG_Sfwd")) {curstring >> varname >> value; ELG_Sfwd = value;}
594 else if(strstr(temp_string.c_str(),"ELG_Cfwd")) {curstring >> varname >> value; ELG_Cfwd = value;}
595 else if(strstr(temp_string.c_str(),"ELG_Nfwd")) {curstring >> varname >> value; ELG_Nfwd = value;}
596 else if(strstr(temp_string.c_str(),"ELG_Szdc")) {curstring >> varname >> value; ELG_Szdc = value;}
597 else if(strstr(temp_string.c_str(),"ELG_Czdc")) {curstring >> varname >> value; ELG_Czdc = value;}
598 else if(strstr(temp_string.c_str(),"ELG_Nzdc")) {curstring >> varname >> value; ELG_Nzdc = value;}
599
600 else if(strstr(temp_string.c_str(),"HAD_Shcal")) {warning("HAD_Shcal","HAD_Scen"); curstring >> varname >> value; HAD_Scen = value;}
601 else if(strstr(temp_string.c_str(),"HAD_Nhcal")) {warning("HAD_Nhcal","HAD_Ncen"); curstring >> varname >> value; HAD_Ncen = value;}
602 else if(strstr(temp_string.c_str(),"HAD_Chcal")) {warning("HAD_Chcal","HAD_Ccen"); curstring >> varname >> value; HAD_Ccen = value;}
603 else if(strstr(temp_string.c_str(),"HAD_Shf")) {warning("HAD_Shf","HAD_Sfwd"); curstring >> varname >> value; HAD_Sfwd = value;}
604 else if(strstr(temp_string.c_str(),"HAD_Nhf")) {warning("HAD_Nhf","HAD_Nfwd"); curstring >> varname >> value; HAD_Nfwd = value;}
605 else if(strstr(temp_string.c_str(),"HAD_Chf")) {warning("HAD_Chf","HAD_Cfwd"); curstring >> varname >> value; HAD_Cfwd = value;}
606
607
608
609 else if(strstr(temp_string.c_str(),"HAD_Scen")) {curstring >> varname >> value; HAD_Scen = value;}
610 else if(strstr(temp_string.c_str(),"HAD_Ncen")) {curstring >> varname >> value; HAD_Ncen = value;}
611 else if(strstr(temp_string.c_str(),"HAD_Ccen")) {curstring >> varname >> value; HAD_Ccen = value;}
612 else if(strstr(temp_string.c_str(),"HAD_Sec")) {curstring >> varname >> value; HAD_Sec = value;}
613 else if(strstr(temp_string.c_str(),"HAD_Nec")) {curstring >> varname >> value; HAD_Nec = value;}
614 else if(strstr(temp_string.c_str(),"HAD_Cec")) {curstring >> varname >> value; HAD_Cec = value;}
615 else if(strstr(temp_string.c_str(),"HAD_Sfwd")) {curstring >> varname >> value; HAD_Sfwd = value;}
616 else if(strstr(temp_string.c_str(),"HAD_Nfwd")) {curstring >> varname >> value; HAD_Nfwd = value;}
617 else if(strstr(temp_string.c_str(),"HAD_Cfwd")) {curstring >> varname >> value; HAD_Cfwd = value;}
618 else if(strstr(temp_string.c_str(),"HAD_Szdc")) {curstring >> varname >> value; HAD_Szdc = value;}
619 else if(strstr(temp_string.c_str(),"HAD_Nzdc")) {curstring >> varname >> value; HAD_Nzdc = value;}
620 else if(strstr(temp_string.c_str(),"HAD_Czdc")) {curstring >> varname >> value; HAD_Czdc = value;}
621
622 else if(strstr(temp_string.c_str(),"ZDC_T_resolution")) {curstring >> varname >> value; ZDC_T_resolution = value;}
623 else if(strstr(temp_string.c_str(),"RP220_T_resolution")) {curstring >> varname >> value; RP220_T_resolution = value;}
624 else if(strstr(temp_string.c_str(),"RP420_T_resolution")) {curstring >> varname >> value; RP420_T_resolution = value;}
625 else if(strstr(temp_string.c_str(),"MU_SmearPt")) {curstring >> varname >> value; MU_SmearPt = value;}
626
627 else if(strstr(temp_string.c_str(),"TRACK_radius")) {curstring >> varname >> ivalue;TRACK_radius = ivalue;}
628 else if(strstr(temp_string.c_str(),"TRACK_length")) {curstring >> varname >> ivalue;TRACK_length = ivalue;}
629 else if(strstr(temp_string.c_str(),"TRACK_bfield_x")) {curstring >> varname >> value; TRACK_bfield_x = value;}
630 else if(strstr(temp_string.c_str(),"TRACK_bfield_y")) {curstring >> varname >> value; TRACK_bfield_y = value;}
631 else if(strstr(temp_string.c_str(),"TRACK_bfield_z")) {curstring >> varname >> value; TRACK_bfield_z = value;}
632 else if(strstr(temp_string.c_str(),"FLAG_bfield")) {curstring >> varname >> ivalue; FLAG_bfield = ivalue;}
633 else if(strstr(temp_string.c_str(),"TRACK_ptmin")) {curstring >> varname >> value; TRACK_ptmin = value;}
634 else if(strstr(temp_string.c_str(),"TRACK_eff")) {curstring >> varname >> ivalue;TRACK_eff = ivalue;}
635
636 else if(strstr(temp_string.c_str(),"TOWER_number")) {curstring >> varname >> ivalue;TOWER_number = ivalue;}
637 else if(strstr(temp_string.c_str(),"TOWER_eta_edges")){
638 curstring >> varname; for(unsigned int i=0; i<TOWER_number+1; i++) {curstring >> value; TOWER_eta_edges[i] = value;} }
639 else if(strstr(temp_string.c_str(),"TOWER_dphi")){
640 curstring >> varname; for(unsigned int i=0; i<TOWER_number; i++) {curstring >> value; TOWER_dphi[i] = value;} }
641
642 else if(strstr(temp_string.c_str(),"PTCUT_elec")) {curstring >> varname >> value; PTCUT_elec = value;}
643 else if(strstr(temp_string.c_str(),"PTCUT_muon")) {curstring >> varname >> value; PTCUT_muon = value;}
644 else if(strstr(temp_string.c_str(),"PTCUT_jet")) {curstring >> varname >> value; PTCUT_jet = value;}
645 else if(strstr(temp_string.c_str(),"PTCUT_gamma")) {curstring >> varname >> value; PTCUT_gamma = value;}
646 else if(strstr(temp_string.c_str(),"PTCUT_taujet")) {curstring >> varname >> value; PTCUT_taujet = value;}
647 else if(strstr(temp_string.c_str(),"ZDC_gamma_E")) {curstring >> varname >> value; ZDC_gamma_E = value;}
648 else if(strstr(temp_string.c_str(),"ZDC_n_E")) {curstring >> varname >> value; ZDC_n_E = value;}
649
650 else if(strstr(temp_string.c_str(),"ISOL_PT")) {curstring >> varname >> value; ISOL_PT = value;}
651 else if(strstr(temp_string.c_str(),"ISOL_Cone")) {curstring >> varname >> value; ISOL_Cone = value;}
652 else if(strstr(temp_string.c_str(),"ISOL_Calo_ET")) {curstring >> varname >> value; ISOL_Calo_ET = value;}
653 else if(strstr(temp_string.c_str(),"ISOL_Calo_Grid")) {curstring >> varname >> ivalue; ISOL_Calo_Grid = ivalue;}
654
655 else if(strstr(temp_string.c_str(),"JET_coneradius")) {curstring >> varname >> value; JET_coneradius = value;}
656 else if(strstr(temp_string.c_str(),"JET_jetalgo")) {curstring >> varname >> ivalue;JET_jetalgo = ivalue;}
657 else if(strstr(temp_string.c_str(),"JET_seed")) {curstring >> varname >> value; JET_seed = value;}
658 else if(strstr(temp_string.c_str(),"JET_Eflow")) {curstring >> varname >> ivalue; JET_Eflow = ivalue;}
659
660 else if(strstr(temp_string.c_str(),"BTAG_b")) {curstring >> varname >> ivalue;BTAG_b = ivalue;}
661 else if(strstr(temp_string.c_str(),"BTAG_mistag_c")) {curstring >> varname >> ivalue;BTAG_mistag_c = ivalue;}
662 else if(strstr(temp_string.c_str(),"BTAG_mistag_l")) {curstring >> varname >> ivalue;BTAG_mistag_l = ivalue;}
663
664 else if(strstr(temp_string.c_str(),"FLAG_vfd")) {curstring >> varname >> ivalue; FLAG_vfd = ivalue;}
665 else if(strstr(temp_string.c_str(),"FLAG_RP")) {curstring >> varname >> ivalue; FLAG_RP = ivalue;}
666 else if(strstr(temp_string.c_str(),"FLAG_trigger")) {curstring >> varname >> ivalue; FLAG_trigger = ivalue;}
667 else if(strstr(temp_string.c_str(),"FLAG_frog")) {curstring >> varname >> ivalue; FLAG_frog = ivalue;}
668 else if(strstr(temp_string.c_str(),"FLAG_lhco")) {curstring >> varname >> ivalue; FLAG_lhco = ivalue;}
669 else if(strstr(temp_string.c_str(),"NEvents_Frog")) {curstring >> varname >> ivalue; NEvents_Frog = ivalue;}
670 else if(strstr(temp_string.c_str(),"NEvents")) {curstring >> varname >> ivalue; NEvents = ivalue;}
671
672 else if(strstr(temp_string.c_str(),"PdgTableFilename")) {curstring >> varname >> svalue; PdgTableFilename = svalue;}
673 }
674
675 // for compatibility with old data cards
676 if(!CEN_max_calo_ec_flag) {
677 cout << "** Warning \'CEN_max_calo_ec\' not found in datacard. **"<< endl;
678 cout << "** Same values will be applied for calorimeter endcaps **"<< endl;
679 cout << "** as for central calorimeters **"<< endl;
680 string text = "** Please update your card ("+ datacard +")";
681 cout << left << setw(67) << text << right << setw(2) << "**" << endl;
682 cout << "** This change is 100\% backward compatibly with older DetectorCard. **" << endl;
683 cout << "** However, please update your DetectorCard. **" << endl;
684 CEN_max_calo_ec = CEN_max_calo_cen;
685 CEN_max_calo_cen = CEN_max_calo_cen/2;
686 ELG_Sec = ELG_Scen;
687 ELG_Nec = ELG_Ncen;
688 ELG_Cec = ELG_Ccen;
689 HAD_Sec = HAD_Scen;
690 HAD_Nec = HAD_Ncen;
691 HAD_Cec = HAD_Ccen;
692 }
693
694 if(ISOL_Calo_Grid%2 ==0) {
695 ISOL_Calo_Grid++;
696 cout <<"** WARNING: ISOL_Calo_Grid is not odd. Set it to "<< ISOL_Calo_Grid << " **" << endl;
697 }
698
699 //jet stuffs not defined in the input datacard
700 JET_overlap = 0.75;
701 // MidPoint algorithm definition
702 JET_M_coneareafraction = 0.25;
703 JET_M_maxpairsize = 2;
704 JET_M_maxiterations = 100;
705 // Define Cone algorithm.
706 JET_C_adjacencycut = 2;
707 JET_C_maxiterations = 100;
708 JET_C_iratch = 1;
709 //Define SISCone algorithm.
710 JET_S_npass = 0;
711 JET_S_protojet_ptmin= 0.0;
712
713 //For Tau-jet definition
714 TAU_energy_scone = 0.15; // radius R of the cone for tau definition, based on energy threshold
715 TAU_track_scone = 0.4; // radius R of the cone for tau definition, based on track number
716 TAU_track_pt = 2; // minimal pt [GeV] for tracks to be considered in tau definition
717 TAU_energy_frac = 0.95; // fraction of energy required in the central part of the cone, for tau jets
718
719}
720
721void RESOLution::Logfile(const string& LogName) {
722
723 // creates the list of good input files
724 // this list is vector<string> inputfiles.
725 ifstream infile(inputfilelist.c_str());
726 vector<string> inputfiles;
727 string filename;
728 while(1) {
729 infile >> filename; // reads the first line of the list
730 if(!infile.good()) break; // quits when at the end of the list
731 ifstream checking_the_file(filename.c_str()); // try to open the file
732 if(!checking_the_file.good()) continue; // skips bad/unknown files
733 else checking_the_file.close(); // close file if found
734 inputfiles.push_back(filename); // append the name to the vector
735 }
736 infile.close();
737
738 ofstream f_out(LogName.c_str());
739
740 f_out <<"**********************************************************************"<< endl;
741 f_out <<"**********************************************************************"<< endl;
742 f_out <<"** **"<< endl;
743 f_out <<"** Welcome to **"<< endl;
744 f_out <<"** **"<< endl;
745 f_out <<"** **"<< endl;
746 f_out <<"** .ddddddd- lL hH **"<< endl;
747 f_out <<"** -Dd` `dD: Ll hH` **"<< endl;
748 f_out <<"** dDd dDd eeee. lL .pp+pp Hh+hhh` -eeee- `sssss **"<< endl;
749 f_out <<"** -Dd `DD ee. ee Ll .Pp. PP Hh. HH. ee. ee sSs **"<< endl;
750 f_out <<"** dD` dDd eEeee: lL. pP. pP hH hH` eEeee:` -sSSSs. **"<< endl;
751 f_out <<"** .Dd :dd eE. LlL PpppPP Hh Hh eE sSS **"<< endl;
752 f_out <<"** dddddd:. eee+: lL. pp. hh. hh eee+ sssssS **"<< endl;
753 f_out <<"** Pp **"<< endl;
754 f_out <<"** **"<< endl;
755 f_out <<"** Delphes, a framework for the fast simulation **"<< endl;
756 f_out <<"** of a generic collider experiment **"<< endl;
757 f_out <<"** **"<< endl;
758 f_out <<"** --- Version 1.7 of Delphes --- **"<< endl;
759 f_out <<"** Last date of change: 7 May 2009 **"<< endl;
760 f_out <<"** **"<< endl;
761 f_out <<"** **"<< endl;
762 f_out <<"** This package uses: **"<< endl;
763 f_out <<"** ------------------ **"<< endl;
764 f_out <<"** FastJet algorithm: Phys. Lett. B641 (2006) [hep-ph/0512210] **"<< endl;
765 f_out <<"** Hector: JINST 2:P09005 (2007) [physics.acc-ph:0707.1198v2] **"<< endl;
766 f_out <<"** FROG: L. Quertenmont, V. Roberfroid [hep-ex/0901.2718v1] **"<< endl;
767 f_out <<"** **"<< endl;
768 f_out <<"** ---------------------------------------------------------------- **"<< endl;
769 f_out <<"** **"<< endl;
770 f_out <<"** Main authors: **"<< endl;
771 f_out <<"** ------------- **"<< endl;
772 f_out <<"** **"<< endl;
773 f_out <<"** Séverine Ovyn Xavier Rouby **"<< endl;
774 f_out <<"** severine.ovyn@uclouvain.be xavier.rouby@cern **"<< endl;
775 f_out <<"** Center for Particle Physics and Phenomenology (CP3) **"<< endl;
776 f_out <<"** Universite Catholique de Louvain (UCL) **"<< endl;
777 f_out <<"** Louvain-la-Neuve, Belgium **"<< endl;
778 f_out <<"** **"<< endl;
779 f_out <<"** ---------------------------------------------------------------- **"<< endl;
780 f_out <<"** **"<< endl;
781 f_out <<"** Former Delphes versions and documentation can be found on : **"<< endl;
782 f_out <<"** http://www.fynu.ucl.ac.be/delphes.html **"<< endl;
783 f_out <<"** **"<< endl;
784 f_out <<"** **"<< endl;
785 f_out <<"** Disclaimer: this program is a beta version of Delphes and **"<< endl;
786 f_out <<"** therefore comes without guarantees. Beware of errors and please **"<< endl;
787 f_out <<"** give us your feedbacks about potential bugs **"<< endl;
788 f_out <<"** **"<< endl;
789 f_out <<"**********************************************************************"<< endl;
790 f_out <<"** **"<< endl;
791 f_out<<"#>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"<<"\n";
792 f_out <<"* *"<< endl;
793 f_out <<"#******************************** *"<<"\n";
794 f_out <<"# Input files *"<<"\n";
795 f_out <<"#******************************** *"<<"\n";
796 f_out << left << setw(22) <<"* Input list "<<""
797 << left << setw(39) << inputfilelist << "" << right << setw(9) << "*"<<"\n";
798 for (unsigned int i =0; i<inputfiles.size(); i++) {
799 f_out << left << setw(22) <<"* - file "<<""
800 << left << setw(43) << inputfiles[i] << "" << right << setw(5) << "*"<<"\n";
801 }
802 if(detectorcard != "")
803 f_out << left << setw(22) <<"* Detector card "<<""
804 << left << setw(39) << detectorcard << "" << right << setw(9) << "*"<<"\n";
805 if(triggercard != "")
806 f_out << left << setw(22) <<"* Trigger card "<<""
807 << left << setw(39) << triggercard << "" << right << setw(9) << "*"<<"\n";
808 f_out<<"* Beam optics : *"<<"\n";
809 f_out << left << setw(22) <<"* - beam 1 "<<""
810 << left << setw(33) << RP_beam1Card << "" << right << setw(15) << "*"<<"\n";
811 f_out << left << setw(22) <<"* - beam 2 "<<""
812 << left << setw(33) << RP_beam2Card << "" << right << setw(15) << "*"<<"\n";
813 f_out << left << setw(22) <<"* Input PDG table " << ""
814 << left << setw(39) << PdgTableFilename << "" << right << setw(9) << "*"<<"\n";
815
816 f_out<<"* *"<<"\n";
817 f_out<<"* *"<<"\n";
818 f_out<<"#******************************** *"<<"\n";
819 f_out<<"# Central detector caracteristics *"<<"\n";
820 f_out<<"#******************************** *"<<"\n";
821 f_out<<"* *"<<"\n";
822 f_out << left << setw(30) <<"* Maximum tracking system: "<<""
823 << left << setw(10) <<CEN_max_tracker <<""<< right << setw(15)<<"*"<<"\n";
824 f_out << left << setw(30) <<"* Maximum central calorimeter: "<<""
825 << left << setw(10) <<CEN_max_calo_cen <<""<< right << setw(15)<<"*"<<"\n";
826 f_out << left << setw(30) <<"* Maximum endcap calorimeter: "<<""
827 << left << setw(10) <<CEN_max_calo_ec <<""<< right << setw(15)<<"*"<<"\n";
828 f_out << left << setw(30) <<"* Maximum central calorimeter: "<<""
829 << left << setw(10) <<CEN_max_calo_fwd <<""<< right << setw(15)<<"*"<<"\n";
830 f_out << left << setw(30) <<"* Muon chambers coverage: "<<""
831 << left << setw(10) <<CEN_max_mu <<""<< right << setw(15)<<"*"<<"\n";
832 f_out<<"* *"<<"\n";
833 if(FLAG_RP==1){
834 f_out<<"#************************************ *"<<"\n";
835 f_out<<"# Very forward Roman Pots switched on *"<<"\n";
836 f_out<<"#************************************ *"<<"\n";
837 f_out<<"* *"<<"\n";
838 f_out << left << setw(55) <<"* Distance of the 220 RP to the IP in meters:"<<""
839 << left << setw(5) <<RP_220_s <<""<< right << setw(10)<<"*"<<"\n";
840 f_out << left << setw(55) <<"* Distance of the 220 RP to the beam in meters:"<<""
841 << left << setw(5) <<RP_220_x <<""<< right << setw(10)<<"*"<<"\n";
842 f_out << left << setw(55) <<"* Distance of the 420 RP to the IP in meters:"<<""
843 << left << setw(5) <<RP_420_s <<""<< right << setw(10)<<"*"<<"\n";
844 f_out << left << setw(55) <<"* Distance of the 420 RP to the beam in meters:"<<""
845 << left << setw(5) <<RP_420_x <<""<< right << setw(10)<<"*"<<"\n";
846 f_out << left << setw(55) <<"* Interaction point at the LHC named: "<<""
847 << left << setw(5) <<RP_IP_name <<""<< right << setw(10)<<"*"<<"\n";
848 f_out << left << setw(35) <<"* Datacard for beam 1: "<<""
849 << left << setw(25) <<RP_beam1Card <<""<< right << setw(10)<<"*"<<"\n";
850 f_out << left << setw(35) <<"* Datacard for beam 2: "<<""
851 << left << setw(25) <<RP_beam2Card <<""<< right << setw(10)<<"*"<<"\n";
852 f_out << left << setw(54) <<"* Beam separation, in meters(hor):"<<""
853 << left << setw(6) << RP_offsetEl_x <<""<< right << setw(10)<<"*"<<"\n";
854 f_out << left << setw(54) <<"* Beam separation, in meters(ver):"<<""
855 << left << setw(6) << RP_offsetEl_y <<""<< right << setw(10)<<"*"<<"\n";
856 f_out << left << setw(54) <<"* Distance from IP for Beam separation (m):"<<""
857 << left << setw(6) <<RP_offsetEl_s <<""<< right << setw(10)<<"*"<<"\n";
858 f_out << left << setw(54) <<"* X offset of beam crossing in micrometers:"<<""
859 << left << setw(6) <<RP_cross_x <<""<< right << setw(10)<<"*"<<"\n";
860 f_out << left << setw(54) <<"* Y offset of beam crossing in micrometers:"<<""
861 << left << setw(6) <<RP_cross_y <<""<< right << setw(10)<<"*"<<"\n";
862 f_out << left << setw(54) <<"* X Angle of beam crossing:"<<""
863 << left << setw(6) <<RP_cross_ang_x <<""<< right << setw(10)<<"*"<<"\n";
864 f_out << left << setw(54) <<"* Y Angle of beam crossing:"<<""
865 << left << setw(6) <<RP_cross_ang_y <<""<< right << setw(10)<<"*"<<"\n";
866 f_out<<"* *"<<"\n";
867 }
868 else {
869 f_out<<"#************************************* *"<<"\n";
870 f_out<<"# Very forward Roman Pots switched off *"<<"\n";
871 f_out<<"#************************************* *"<<"\n";
872 f_out<<"* *"<<"\n";
873 }
874 if(FLAG_vfd==1){
875 f_out<<"#************************************** *"<<"\n";
876 f_out<<"# Very forward calorimeters switched on *"<<"\n";
877 f_out<<"#************************************** *"<<"\n";
878 f_out<<"* *"<<"\n";
879 f_out << left << setw(55) <<"* Minimum very forward calorimeter: "<<""
880 << left << setw(5) <<VFD_min_calo_vfd <<""<< right << setw(10)<<"*"<<"\n";
881 f_out << left << setw(55) <<"* Maximum very forward calorimeter: "<<""
882 << left << setw(5) <<VFD_max_calo_vfd <<""<< right << setw(10)<<"*"<<"\n";
883 f_out << left << setw(55) <<"* Minimum coverage zero_degree calorimeter "<<""
884 << left << setw(5) <<VFD_min_zdc <<""<< right << setw(10)<<"*"<<"\n";
885 f_out << left << setw(55) <<"* Distance of the ZDC to the IP, in meters: "<<""
886 << left << setw(5) <<VFD_s_zdc <<""<< right << setw(10)<<"*"<<"\n";
887 f_out<<"* *"<<"\n";
888 }
889 else {
890 f_out<<"#*************************************** *"<<"\n";
891 f_out<<"# Very forward calorimeters switched off *"<<"\n";
892 f_out<<"#*************************************** *"<<"\n";
893 f_out<<"* *"<<"\n";
894 }
895
896 f_out<<"#************************************ *"<<"\n";
897 f_out<<"# Electromagnetic smearing parameters *"<<"\n";
898 f_out<<"#************************************ *"<<"\n";
899 f_out<<"* *"<<"\n";
900 //# \sigma/E = C + N/E + S/\sqrt{E}
901 f_out << left << setw(30) <<"* S term for central ECAL: "<<""
902 << left << setw(30) <<ELG_Scen <<""<< right << setw(10)<<"*"<<"\n";
903 f_out << left << setw(30) <<"* N term for central ECAL: "<<""
904 << left << setw(30) <<ELG_Ncen <<""<< right << setw(10)<<"*"<<"\n";
905 f_out << left << setw(30) <<"* C term for central ECAL: "<<""
906 << left << setw(30) <<ELG_Ccen <<""<< right << setw(10)<<"*"<<"\n";
907 f_out << left << setw(30) <<"* S term for ECAL end-cap: "<<""
908 << left << setw(30) <<ELG_Sec <<""<< right << setw(10)<<"*"<<"\n";
909 f_out << left << setw(30) <<"* N term for ECAL end-cap: "<<""
910 << left << setw(30) <<ELG_Nec <<""<< right << setw(10)<<"*"<<"\n";
911 f_out << left << setw(30) <<"* C term for ECAL end-cap: "<<""
912 << left << setw(30) <<ELG_Cec <<""<< right << setw(10)<<"*"<<"\n";
913
914
915 f_out << left << setw(30) <<"* S term for FCAL: "<<""
916 << left << setw(30) <<ELG_Sfwd <<""<< right << setw(10)<<"*"<<"\n";
917 f_out << left << setw(30) <<"* N term for FCAL: "<<""
918 << left << setw(30) <<ELG_Nfwd <<""<< right << setw(10)<<"*"<<"\n";
919 f_out << left << setw(30) <<"* C term for FCAL: "<<""
920 << left << setw(30) <<ELG_Cfwd <<""<< right << setw(10)<<"*"<<"\n";
921 f_out << left << setw(30) <<"* S term for ZDC: "<<""
922 << left << setw(30) <<ELG_Szdc <<""<< right << setw(10)<<"*"<<"\n";
923 f_out << left << setw(30) <<"* N term for ZDC: "<<""
924 << left << setw(30) <<ELG_Nzdc <<""<< right << setw(10)<<"*"<<"\n";
925 f_out << left << setw(30) <<"* C term for ZDC: "<<""
926 << left << setw(30) <<ELG_Czdc <<""<< right << setw(10)<<"*"<<"\n";
927
928 f_out<<"* *"<<"\n";
929 f_out<<"#***************************** *"<<"\n";
930 f_out<<"# Hadronic smearing parameters *"<<"\n";
931 f_out<<"#***************************** *"<<"\n";
932 f_out<<"* *"<<"\n";
933 f_out << left << setw(30) <<"* S term for central HCAL: "<<""
934 << left << setw(30) <<HAD_Scen <<""<< right << setw(10)<<"*"<<"\n";
935 f_out << left << setw(30) <<"* N term for central HCAL: "<<""
936 << left << setw(30) <<HAD_Ncen <<""<< right << setw(10)<<"*"<<"\n";
937 f_out << left << setw(30) <<"* C term for central HCAL: "<<""
938 << left << setw(30) <<HAD_Ccen <<""<< right << setw(10)<<"*"<<"\n";
939 f_out << left << setw(30) <<"* S term for HCAL endcap: "<<""
940 << left << setw(30) <<HAD_Sec <<""<< right << setw(10)<<"*"<<"\n";
941 f_out << left << setw(30) <<"* N term for HCAL endcap: "<<""
942 << left << setw(30) <<HAD_Nec <<""<< right << setw(10)<<"*"<<"\n";
943 f_out << left << setw(30) <<"* C term for HCAL endcap: "<<""
944 << left << setw(30) <<HAD_Cec <<""<< right << setw(10)<<"*"<<"\n";
945 f_out << left << setw(30) <<"* S term for FCAL: "<<""
946 << left << setw(30) <<HAD_Sfwd <<""<< right << setw(10)<<"*"<<"\n";
947 f_out << left << setw(30) <<"* N term for FCAL: "<<""
948 << left << setw(30) <<HAD_Nfwd <<""<< right << setw(10)<<"*"<<"\n";
949 f_out << left << setw(30) <<"* C term for FCAL: "<<""
950 << left << setw(30) <<HAD_Cfwd <<""<< right << setw(10)<<"*"<<"\n";
951 f_out << left << setw(30) <<"* S term for ZDC: "<<""
952 << left << setw(30) <<HAD_Szdc <<""<< right << setw(10)<<"*"<<"\n";
953 f_out << left << setw(30) <<"* N term for ZDC: "<<""
954 << left << setw(30) <<HAD_Nzdc <<""<< right << setw(10)<<"*"<<"\n";
955 f_out << left << setw(30) <<"* C term for ZDC: "<<""
956 << left << setw(30) <<HAD_Czdc <<""<< right << setw(10)<<"*"<<"\n";
957
958 f_out<<"* *"<<"\n";
959 f_out<<"#************************* *"<<"\n";
960 f_out<<"# Time smearing parameters *"<<"\n";
961 f_out<<"#************************* *"<<"\n";
962 f_out<<"* *"<<"\n";
963 f_out << left << setw(55) <<"* Time resolution for ZDC : "<<""
964 << left << setw(5) <<ZDC_T_resolution <<""<< right << setw(10)<<"*"<<"\n";
965 f_out << left << setw(55) <<"* Time resolution for RP220 : "<<""
966 << left << setw(5) <<RP220_T_resolution <<""<< right << setw(10)<<"*"<<"\n";
967 f_out << left << setw(55) <<"* Time resolution for RP420 : "<<""
968 << left << setw(5) <<RP420_T_resolution <<""<< right << setw(10)<<"*"<<"\n";
969 f_out<<"* *"<<"\n";
970
971 f_out<<"* *"<<"\n";
972 f_out<<"#************************* *"<<"\n";
973 f_out<<"# Muon smearing parameters *"<<"\n";
974 f_out<<"#************************* *"<<"\n";
975 f_out<<"* *"<<"\n";
976 f_out << left << setw(55) <<"* PT resolution for muons : "<<""
977 << left << setw(5) <<MU_SmearPt <<""<< right << setw(10)<<"*"<<"\n";
978 f_out<<"* *"<<"\n";
979 if(FLAG_bfield==1){
980 f_out<<"#*************************** *"<<"\n";
981 f_out<<"# Magnetic field switched on *"<<"\n";
982 f_out<<"#*************************** *"<<"\n";
983 f_out<<"* *"<<"\n";
984 f_out << left << setw(55) <<"* Radius of the BField coverage: "<<""
985 << left << setw(5) <<TRACK_radius <<""<< right << setw(10)<<"*"<<"\n";
986 f_out << left << setw(55) <<"* Length of the BField coverage: "<<""
987 << left << setw(5) <<TRACK_length <<""<< right << setw(10)<<"*"<<"\n";
988 f_out << left << setw(55) <<"* BField X component: "<<""
989 << left << setw(5) <<TRACK_bfield_x <<""<< right << setw(10)<<"*"<<"\n";
990 f_out << left << setw(55) <<"* BField Y component: "<<""
991 << left << setw(5) <<TRACK_bfield_y <<""<< right << setw(10)<<"*"<<"\n";
992 f_out << left << setw(55) <<"* BField Z component: "<<""
993 << left << setw(5) <<TRACK_bfield_z <<""<< right << setw(10)<<"*"<<"\n";
994 f_out << left << setw(55) <<"* Minimal pT needed to reach the calorimeter [GeV]: "<<""
995 << left << setw(10) <<TRACK_ptmin <<""<< right << setw(5)<<"*"<<"\n";
996 f_out << left << setw(55) <<"* Efficiency associated to the tracking: "<<""
997 << left << setw(10) <<TRACK_eff <<""<< right << setw(5)<<"*"<<"\n";
998 f_out<<"* *"<<"\n";
999 }
1000 else {
1001 f_out<<"#**************************** *"<<"\n";
1002 f_out<<"# Magnetic field switched off *"<<"\n";
1003 f_out<<"#**************************** *"<<"\n";
1004 f_out << left << setw(55) <<"* Minimal pT needed to reach the calorimeter [GeV]: "<<""
1005 << left << setw(10) <<TRACK_ptmin <<""<< right << setw(5)<<"*"<<"\n";
1006 f_out << left << setw(55) <<"* Efficiency associated to the tracking: "<<""
1007 << left << setw(10) <<TRACK_eff <<""<< right << setw(5)<<"*"<<"\n";
1008 f_out<<"* *"<<"\n";
1009 }
1010 f_out<<"#******************** *"<<"\n";
1011 f_out<<"# Calorimetric Towers *"<<"\n";
1012 f_out<<"#******************** *"<<"\n";
1013 f_out << left << setw(55) <<"* Number of calorimetric towers in eta, for eta>0: "<<""
1014 << left << setw(5) << TOWER_number <<""<< right << setw(10)<<"*"<<"\n";
1015 f_out << left << setw(55) <<"* Tower edges in eta, for eta>0: "<<"" << right << setw(15)<<"*"<<"\n";
1016 f_out << "* ";
1017 for (unsigned int i=0; i<TOWER_number+1; i++) {
1018 f_out << left << setw(7) << TOWER_eta_edges[i];
1019 if(!( (i+1) %9 )) f_out << right << setw(3) << "*" << "\n" << "* ";
1020 }
1021 for (unsigned int i=(TOWER_number+1)%9; i<9; i++) f_out << left << setw(7) << "";
1022 f_out << right << setw(3)<<"*"<<"\n";
1023 f_out << left << setw(55) <<"* Tower sizes in phi, for eta>0 [degree]:"<<"" << right << setw(15)<<"*"<<"\n";
1024 f_out << "* ";
1025 for (unsigned int i=0; i<TOWER_number; i++) {
1026 f_out << left << setw(7) << TOWER_dphi[i];
1027 if(!( (i+1) %9 )) f_out << right << setw(3) << "*" << "\n" << "* ";
1028 }
1029 for (unsigned int i=(TOWER_number)%9; i<9; i++) f_out << left << setw(7) << "";
1030 f_out << right << setw(3)<<"*"<<"\n";
1031 f_out<<"* *"<<"\n";
1032 f_out<<"#******************* *"<<"\n";
1033 f_out<<"# Minimum pT's [GeV] *"<<"\n";
1034 f_out<<"#******************* *"<<"\n";
1035 f_out<<"* *"<<"\n";
1036 f_out << left << setw(40) <<"* Minimum pT for electrons: "<<""
1037 << left << setw(20) <<PTCUT_elec <<""<< right << setw(10)<<"*"<<"\n";
1038 f_out << left << setw(40) <<"* Minimum pT for muons: "<<""
1039 << left << setw(20) <<PTCUT_muon <<""<< right << setw(10)<<"*"<<"\n";
1040 f_out << left << setw(40) <<"* Minimum pT for jets: "<<""
1041 << left << setw(20) <<PTCUT_jet <<""<< right << setw(10)<<"*"<<"\n";
1042 f_out << left << setw(40) <<"* Minimum pT for Tau-jets: "<<""
1043 << left << setw(20) <<PTCUT_taujet <<""<< right << setw(10)<<"*"<<"\n";
1044 f_out << left << setw(40) <<"* Minimum pT for photons: "<<""
1045 << left << setw(20) <<PTCUT_gamma <<""<< right << setw(10)<<"*"<<"\n";
1046 f_out << left << setw(40) <<"* Minimum E for photons in ZDC: "<<""
1047 << left << setw(20) <<ZDC_gamma_E <<""<< right << setw(10)<<"*"<<"\n";
1048 f_out << left << setw(40) <<"* Minimum E for neutrons in ZDC: "<<""
1049 << left << setw(20) <<ZDC_n_E <<""<< right << setw(10)<<"*"<<"\n";
1050
1051 f_out<<"* *"<<"\n";
1052 f_out<<"#******************* *"<<"\n";
1053 f_out<<"# Isolation criteria *"<<"\n";
1054 f_out<<"#******************* *"<<"\n";
1055 f_out<<"* *"<<"\n";
1056 f_out << left << setw(40) <<"* Minimum pT for tracks [GeV]: "<<""
1057 << left << setw(20) <<ISOL_PT <<""<< right << setw(10)<<"*"<<"\n";
1058 f_out << left << setw(40) <<"* Cone for isolation criteria: "<<""
1059 << left << setw(20) <<ISOL_Cone <<""<< right << setw(10)<<"*"<<"\n";
1060
1061 if(ISOL_Calo_ET > 1E98) f_out<<"# No Calorimetric isolation applied *"<<"\n";
1062 else {
1063 f_out << left << setw(40) <<"* Minimum ET for towers [GeV]: "<<""
1064 << left << setw(20) <<ISOL_Calo_ET <<""<< right << setw(10)<<"*"<<"\n";
1065 f_out << left << setw(40) <<"* Grid size (NxN) for calorimetric isolation: "<<""
1066 << left << setw(20) <<ISOL_Calo_Grid <<""<< right << setw(4)<<"*"<<"\n";
1067 }
1068
1069
1070 f_out<<"* *"<<"\n";
1071 f_out<<"#*************** *"<<"\n";
1072 f_out<<"# Jet definition *"<<"\n";
1073 f_out<<"#*************** *"<<"\n";
1074 if(JET_Eflow)
1075 {
1076 f_out<<"#*************** *"<<"\n";
1077 f_out<<"#* Running considering perfect energy flow on the tracker coverage *"<<"\n";
1078 }
1079 else
1080 {
1081 f_out<<"#* Running considering no energy flow on the tracker coverage *"<<"\n";
1082 f_out<<"#* --> jet algo applied on the calorimetric towers *"<<"\n";
1083 }
1084 f_out<<"* *"<<"\n";
1085 f_out<<"* Six algorithms are currently available: *"<<"\n";
1086 f_out<<"* - 1) CDF cone algorithm, *"<<"\n";
1087 f_out<<"* - 2) CDF MidPoint algorithm, *"<<"\n";
1088 f_out<<"* - 3) SIScone algorithm, *"<<"\n";
1089 f_out<<"* - 4) kt algorithm, *"<<"\n";
1090 f_out<<"* - 5) Cambrigde/Aachen algorithm, *"<<"\n";
1091 f_out<<"* - 6) Anti-kt algorithm. *"<<"\n";
1092 f_out<<"* *"<<"\n";
1093 f_out<<"* You have chosen *"<<"\n";
1094 switch(JET_jetalgo) {
1095 default:
1096 case 1: {
1097 f_out<<"* CDF JetClu jet algorithm with parameters: *"<<"\n";
1098 f_out << left << setw(40) <<"* - Seed threshold: "<<""
1099 << left << setw(10) <<JET_seed <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1100 f_out << left << setw(40) <<"* - Cone radius: "<<""
1101 << left << setw(10) <<JET_coneradius <<""<< right << setw(20)<<"*"<<"\n";
1102 f_out << left << setw(40) <<"* - Adjacency cut: "<<""
1103 << left << setw(10) <<JET_C_adjacencycut <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1104 f_out << left << setw(40) <<"* - Max iterations: "<<""
1105 << left << setw(10) <<JET_C_maxiterations <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1106 f_out << left << setw(40) <<"* - Iratch: "<<""
1107 << left << setw(10) <<JET_C_iratch <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1108 f_out << left << setw(40) <<"* - Overlap threshold: "<<""
1109 << left << setw(10) <<JET_overlap <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1110 }
1111 break;
1112 case 2: {
1113 f_out<<"* CDF midpoint jet algorithm with parameters: *"<<"\n";
1114 f_out << left << setw(40) <<"* - Seed threshold: "<<""
1115 << left << setw(20) <<JET_seed <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1116 f_out << left << setw(40) <<"* - Cone radius: "<<""
1117 << left << setw(20) <<JET_coneradius <<""<< right << setw(10)<<"*"<<"\n";
1118 f_out << left << setw(40) <<"* - Cone area fraction:"<<""
1119 << left << setw(20) <<JET_M_coneareafraction <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1120 f_out << left << setw(40) <<"* - Maximum pair size: "<<""
1121 << left << setw(20) <<JET_M_maxpairsize <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1122 f_out << left << setw(40) <<"* - Max iterations: "<<""
1123 << left << setw(20) <<JET_M_maxiterations <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1124 f_out << left << setw(40) <<"* - Overlap threshold: "<<""
1125 << left << setw(20) <<JET_overlap <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1126 }
1127 break;
1128 case 3: {
1129 f_out <<"* SISCone jet algorithm with parameters: *"<<"\n";
1130 f_out << left << setw(40) <<"* - Cone radius: "<<""
1131 << left << setw(20) <<JET_coneradius <<""<< right << setw(10)<<"*"<<"\n";
1132 f_out << left << setw(40) <<"* - Overlap threshold: "<<""
1133 << left << setw(20) <<JET_overlap <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1134 f_out << left << setw(40) <<"* - Number pass max: "<<""
1135 << left << setw(20) <<JET_S_npass <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1136 f_out << left << setw(40) <<"* - Minimum pT for protojet: "<<""
1137 << left << setw(20) <<JET_S_protojet_ptmin <<""<< right << setw(10)<<"! not in datacard *"<<"\n";
1138 }
1139 break;
1140 case 4: {
1141 f_out <<"* KT jet algorithm with parameters: *"<<"\n";
1142 f_out << left << setw(40) <<"* - Cone radius: "<<""
1143 << left << setw(20) <<JET_coneradius <<""<< right << setw(10)<<"*"<<"\n";
1144 }
1145 break;
1146 case 5: {
1147 f_out <<"* Cambridge/Aachen jet algorithm with parameters: *"<<"\n";
1148 f_out << left << setw(40) <<"* - Cone radius: "<<""
1149 << left << setw(20) <<JET_coneradius <<""<< right << setw(10)<<"*"<<"\n";
1150 }
1151 break;
1152 case 6: {
1153 f_out <<"* Anti-kt jet algorithm with parameters: *"<<"\n";
1154 f_out << left << setw(40) <<"* - Cone radius: "<<""
1155 << left << setw(20) <<JET_coneradius <<""<< right << setw(10)<<"*"<<"\n";
1156 }
1157 break;
1158 }
1159 f_out<<"* *"<<"\n";
1160 f_out<<"#****************************** *"<<"\n";
1161 f_out<<"# Tau-jet definition parameters *"<<"\n";
1162 f_out<<"#****************************** *"<<"\n";
1163 f_out<<"* *"<<"\n";
1164 f_out << left << setw(45) <<"* Cone radius for calorimeter tagging: "<<""
1165 << left << setw(5) <<TAU_energy_scone <<""<< right << setw(20)<<"*"<<"\n";
1166 f_out << left << setw(45) <<"* Fraction of energy in the small cone: "<<""
1167 << left << setw(5) <<TAU_energy_frac*100 <<""<< right << setw(20)<<"! not in datacard *"<<"\n";
1168 f_out << left << setw(45) <<"* Cone radius for tracking tagging: "<<""
1169 << left << setw(5) <<TAU_track_scone <<""<< right << setw(20)<<"*"<<"\n";
1170 f_out << left << setw(45) <<"* Minimum track pT [GeV]: "<<""
1171 << left << setw(5) <<TAU_track_pt <<""<< right << setw(20)<<"*"<<"\n";
1172 f_out<<"* *"<<"\n";
1173 f_out<<"#*************************** *"<<"\n";
1174 f_out<<"# B-tagging efficiencies [%] *"<<"\n";
1175 f_out<<"#*************************** *"<<"\n";
1176 f_out<<"* *"<<"\n";
1177 f_out << left << setw(50) <<"* Efficiency to tag a \"b\" as a b-jet: "<<""
1178 << left << setw(10) <<BTAG_b <<""<< right << setw(10)<<"*"<<"\n";
1179 f_out << left << setw(50) <<"* Efficiency to mistag a c-jet as a b-jet: "<<""
1180 << left << setw(10) <<BTAG_mistag_c <<""<< right << setw(10)<<"*"<<"\n";
1181 f_out << left << setw(50) <<"* Efficiency to mistag a light jet as a b-jet: "<<""
1182 << left << setw(10) <<BTAG_mistag_l <<""<< right << setw(10)<<"*"<<"\n";
1183 f_out<<"* *"<<"\n";
1184 f_out<<"* *"<<"\n";
1185 f_out<<"#....................................................................*"<<"\n";
1186 f_out<<"#>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"<<"\n";
1187
1188 f_out.close();
1189}
1190
1191// **********Provides the smeared TLorentzVector for the electrons********
1192// Smears the electron energy, and changes the 4-momentum accordingly
1193// different smearing if the electron is central (eta < 2.5) or forward
1194void RESOLution::SmearElectron(TLorentzVector &electron) {
1195 // the 'electron' variable will be changed by the function
1196 float energy = electron.E(); // before smearing
1197 float energyS = 0.0; // after smearing // \sigma/E = C + N/E + S/\sqrt{E}
1198
1199 if(fabs(electron.Eta()) < CEN_max_tracker) { // if the electron is inside the tracker
1200 energyS = gRandom->Gaus(energy, sqrt(
1201 pow(ELG_Ncen,2) +
1202 pow(ELG_Ccen*energy,2) +
1203 pow(ELG_Scen*sqrt(energy),2) ));
1204 }
1205 else if(fabs(electron.Eta()) >= CEN_max_tracker && fabs(electron.Eta()) < CEN_max_calo_ec){
1206 energyS = gRandom->Gaus(energy, sqrt(
1207 pow(ELG_Nec,2) +
1208 pow(ELG_Cec*energy,2) +
1209 pow(ELG_Sec*sqrt(energy),2) ) );
1210 }
1211 else if(fabs(electron.Eta()) >= CEN_max_calo_ec && fabs(electron.Eta()) < CEN_max_calo_fwd){
1212 energyS = gRandom->Gaus(energy, sqrt(
1213 pow(ELG_Nfwd,2) +
1214 pow(ELG_Cfwd*energy,2) +
1215 pow(ELG_Sfwd*sqrt(energy),2) ) );
1216 }
1217 electron.SetPtEtaPhiE(energyS/cosh(electron.Eta()), electron.Eta(), electron.Phi(), energyS);
1218 if(electron.E() < 0)electron.SetPxPyPzE(0,0,0,0); // no negative values after smearing !
1219}
1220
1221
1222// **********Provides the smeared TLorentzVector for the muons********
1223// Smears the muon pT and changes the 4-momentum accordingly
1224void RESOLution::SmearMu(TLorentzVector &muon) {
1225 // the 'muon' variable will be changed by the function
1226 float pt = muon.Pt(); // before smearing
1227 float ptS=pt;
1228
1229 if(fabs(muon.Eta()) < CEN_max_mu )
1230 {
1231 ptS = gRandom->Gaus(pt, MU_SmearPt*pt ); // after smearing // \sigma/E = C + N/E + S/\sqrt{E}
1232 }
1233 muon.SetPtEtaPhiE(ptS, muon.Eta(), muon.Phi(), ptS*cosh(muon.Eta()));
1234
1235 if(muon.E() < 0)muon.SetPxPyPzE(0,0,0,0); // no negative values after smearing !
1236}
1237
1238
1239// **********Provides the smeared TLorentzVector for the hadrons********
1240// Smears the hadron 4-momentum
1241void RESOLution::SmearHadron(TLorentzVector &hadron, const float frac)
1242 // the 'hadron' variable will be changed by the function
1243 // the 'frac' variable describes the long-living particles. Should be 0.7 for K0S and Lambda, 1. otherwise
1244{
1245 float energy = hadron.E(); // before smearing
1246 float energyS = 0.0; // after smearing // \sigma/E = C + N/E + S/\sqrt{E}
1247 float energy_ecal = (1.0 - frac)*energy; // electromagnetic calorimeter
1248 float energy_hcal = frac*energy; // hadronic calorimeter
1249 // frac takes into account the decay of long-living particles, that decay in the calorimeters
1250 // some of the particles decay mostly in the ecal, some mostly in the hcal
1251
1252 float energyS1,energyS2;
1253 if(fabs(hadron.Eta()) < CEN_max_calo_cen) {
1254 energyS1 = gRandom->Gaus(energy_hcal, sqrt(
1255 pow(HAD_Ncen,2) +
1256 pow(HAD_Ccen*energy_hcal,2) +
1257 pow(HAD_Scen*sqrt(energy_hcal),2) )) ;
1258
1259 energyS2 = gRandom->Gaus(energy_ecal, sqrt(
1260 pow(ELG_Ncen,2) +
1261 pow(ELG_Ccen*energy_ecal,2) +
1262 pow(ELG_Scen*sqrt(energy_ecal),2) ) );
1263
1264 energyS = ((energyS1>0)?energyS1:0) + ((energyS2>0)?energyS2:0);
1265 }
1266 else if(fabs(hadron.Eta()) >= CEN_max_calo_cen && fabs(hadron.Eta()) < CEN_max_calo_ec) {
1267 energyS1 = gRandom->Gaus(energy_hcal, sqrt(
1268 pow(HAD_Nec,2) +
1269 pow(HAD_Cec*energy_hcal,2) +
1270 pow(HAD_Sec*sqrt(energy_hcal),2) )) ;
1271
1272 energyS2 = gRandom->Gaus(energy_ecal, sqrt(
1273 pow(ELG_Nec,2) +
1274 pow(ELG_Cec*energy_ecal,2) +
1275 pow(ELG_Sec*sqrt(energy_ecal),2) ) );
1276
1277 energyS = ((energyS1>0)?energyS1:0) + ((energyS2>0)?energyS2:0);
1278 }
1279 else if(fabs(hadron.Eta()) >= CEN_max_calo_ec && fabs(hadron.Eta()) < CEN_max_calo_fwd){
1280 energyS = gRandom->Gaus(energy, sqrt(
1281 pow(HAD_Nfwd,2) +
1282 pow(HAD_Cfwd*energy,2) +
1283 pow(HAD_Sfwd*sqrt(energy),2) ));
1284}
1285
1286
1287
1288 hadron.SetPtEtaPhiE(energyS/cosh(hadron.Eta()),hadron.Eta(), hadron.Phi(), energyS);
1289
1290 if(hadron.E() < 0)hadron.SetPxPyPzE(0,0,0,0);
1291}
1292
1293//******************************************************************************************
1294
1295//void RESOLution::SortedVector(vector<ParticleUtil> &vect)
1296void RESOLution::SortedVector(vector<D_Particle> &vect)
1297{
1298 int i,j = 0;
1299 TLorentzVector tmp;
1300 bool en_desordre = true;
1301 int entries=vect.size();
1302 for(i = 0 ; (i < entries) && en_desordre; i++)
1303 {
1304 en_desordre = false;
1305 for(j = 1 ; j < entries - i ; j++)
1306 {
1307 if ( vect[j].Pt() > vect[j-1].Pt() )
1308 {
1309 //ParticleUtil tmp = vect[j-1];
1310 D_Particle tmp = vect[j-1];
1311 vect[j-1] = vect[j];
1312 vect[j] = tmp;
1313 en_desordre = true;
1314 }
1315 }
1316 }
1317}
1318
1319// **********Provides the energy in the cone of radius TAU_CONE_ENERGY for the tau identification********
1320// to be taken into account, a calo tower should
1321// 1) have a transverse energy \f$ E_T = \sqrt{E_X^2 + E_Y^2} \f$ above a given threshold
1322// 2) be inside a cone with a radius R and the axis defined by (eta,phi)
1323double RESOLution::EnergySmallCone(const vector<PhysicsTower> &towers, const float eta, const float phi) {
1324 double Energie=0;
1325 for(unsigned int i=0; i < towers.size(); i++) {
1326 if(towers[i].fourVector.pt() < JET_seed) continue;
1327 if((DeltaR(phi,eta,towers[i].fourVector.phi(),towers[i].fourVector.eta()) < TAU_energy_scone)) {
1328 Energie += towers[i].fourVector.E;
1329 }
1330 }
1331 return Energie;
1332}
1333
1334
1335// **********Provides the number of tracks in the cone of radius TAU_CONE_TRACKS for the tau identification********
1336// to be taken into account, a track should
1337// 1) avec a transverse momentum \$f p_T \$ above a given threshold
1338// 2) be inside a cone with a radius R and the axis defined by (eta,phi)
1339// IMPORTANT REMARK !!!!!
1340// NEW : "charge" will contain the sum of all charged tracks in the cone TAU_track_scone
1341unsigned int RESOLution::NumTracks(float& charge, const vector<TRootTracks> &tracks, const float pt_track, const float eta, const float phi) {
1342 unsigned int numbtrack=0; // number of track in the tau-jet cone, which is smaller than R;
1343 charge=0;
1344 for(unsigned int i=0; i < tracks.size(); i++) {
1345 if(tracks[i].PT < pt_track ) continue;
1346 //float dr = DeltaR(phi,eta,tracks[i].PhiOuter,tracks[i].EtaOuter);
1347 float dr = DeltaR(phi,eta,tracks[i].Phi,tracks[i].Eta);
1348 if (dr > TAU_track_scone) continue;
1349 numbtrack++;
1350 charge += tracks[i].Charge; // total charge in the cone for Tau-jet
1351 }
1352 return numbtrack;
1353}
1354
1355//*** Returns the PID of the particle with the highest energy, in a cone with a radius CONERADIUS and an axis (eta,phi) *********
1356//used by Btaggedjet
1357///// Attention : bug removed => CONERADIUS/2 -> CONERADIUS !!
1358int RESOLution::Bjets(const TSimpleArray<TRootC::GenParticle> &subarray, const float& eta, const float& phi) {
1359 float emax=0;
1360 int Ppid=0;
1361 if(subarray.GetEntries()>0) {
1362 for(int i=0; i < subarray.GetEntries();i++) { // should have pt>PT_JETMIN and a small cone radius (r<CONE_JET)
1363 float genDeltaR = DeltaR(subarray[i]->Phi,subarray[i]->Eta,phi,eta);
1364 if(genDeltaR < JET_coneradius && subarray[i]->E > emax) {
1365 emax=subarray[i]->E;
1366 Ppid=abs(subarray[i]->PID);
1367 }
1368 }
1369 }
1370 return Ppid;
1371}
1372
1373
1374//******************** Simulates the b-tagging efficiency for real bjet, or the misendentification for other jets****************
1375bool RESOLution::Btaggedjet(const TLorentzVector &JET, const TSimpleArray<TRootC::GenParticle> &subarray) {
1376 if( rand()%100 < (BTAG_b+1) && Bjets(subarray,JET.Eta(),JET.Phi())==pB ) return true; // b-tag of b-jets is 40%
1377 else if( rand()%100 < (BTAG_mistag_c+1) && Bjets(subarray,JET.Eta(),JET.Phi())==pC ) return true; // b-tag of c-jets is 10%
1378 else if( rand()%100 < (BTAG_mistag_l+1) && Bjets(subarray,JET.Eta(),JET.Phi())!=0) return true; // b-tag of light jets is 1%
1379 return false;
1380}
1381
1382//***********************Isolation criteria***********************
1383//****************************************************************
1384bool RESOLution::Isolation(const D_Particle& part, const vector<TRootTracks> &tracks, const float& pt_second_track, const float& isolCone, float& ptiso )
1385{
1386 bool isolated = false;
1387 ptiso = 0; // sum of all track pt in isolation cone
1388 float deltar=1E99; // Initial value; should be high; no further repercussion
1389
1390 // loop on all tracks, with p_t above threshold, close enough from the charged lepton
1391 for(unsigned int i=0; i < tracks.size(); i++) {
1392 if(tracks[i].PT < pt_second_track) continue; // ptcut on tracks
1393 float genDeltaR = DeltaR(part.Phi(),part.Eta(),tracks[i].Phi,tracks[i].Eta);
1394 if(
1395 (genDeltaR > deltar) ||
1396 (genDeltaR==0) // rejets the track of the particle itself
1397 ) continue ;
1398 deltar=genDeltaR; // finds the closest track
1399
1400 // as long as (genDeltaR==0) is put above, the particle itself is not taken into account
1401 if( genDeltaR < ISOL_Cone) ptiso += tracks[i].PT; // dR cut on tracks
1402 }
1403 if(deltar > isolCone) isolated = true;
1404 return isolated;
1405}
1406
1407// ******* Calorimetric isolation
1408float RESOLution::CaloIsolation(const D_Particle& part, const D_CaloTowerList & towers, const float iPhi, const float iEta) {
1409 // etrat, which is a percentage between 00 and 99. It is the ratio of the transverse energy
1410 // in a 3×3 grid surrounding the muon to the pT of the muon. For well-isolated muons, both ptiso and etrat will be small.
1411 if(ISOL_Calo_ET>1E10) return UNDEFINED; // avoid doing anything unreasonable...
1412 float et_sum=0;
1413 // available parameters: ISOL_Calo_ET , ISOL_Calo_Grid
1414 // Get the EtaCalo/PhiCalo of the muon ;
1415 // transform it into iEta/iPhi to get the towers, and their neighbourh (i-1, i-2, etc)
1416
1417 unsigned int N = ISOL_Calo_Grid;
1418 int index= iUNDEFINED; // index of the central tower of the grid in TOWER_eta_edges[.];
1419 // !! TOWER_eta_edges is only with eta>0
1420 // finds the index of the central tower of the NxN grid
1421 for (unsigned int i=1; i< TOWER_number+1; i++) {
1422 if(fabs(iEta) >= TOWER_eta_edges[i-1] && fabs(iEta) < TOWER_eta_edges[i]) {
1423 index = i-1;
1424 break;
1425 }
1426 }
1427 if(index != iUNDEFINED) {
1428 // finds the size in phi of the cells for this eta
1429 float dphi = TOWER_dphi[index]*pi/180.; // in rad
1430
1431 //cout << "Grid " << " ----------\n";
1432 for (unsigned int i_eta=0; i_eta<N; i_eta++) {
1433 unsigned int real_index = (iEta>0) ? index+i_eta-(N-1)/2 : index+1+i_eta-(N-1)/2 ;
1434 float eta_ith_tower = TOWER_eta_edges[real_index];
1435 if(iEta<0) eta_ith_tower *= -1;
1436
1437 for (unsigned int i_phi=0; i_phi<N; i_phi++) {
1438 float phi_ith_tower = iPhi + (float)(i_phi - (N-1.)/2.)*dphi;
1439 D_CaloTower calMuon(towers.getElement(eta_ith_tower,phi_ith_tower));
1440 if(calMuon.getEta() != UNDEFINED && calMuon.getE() > ISOL_Calo_ET) {
1441 et_sum += calMuon.getE();
1442 //cout << "eta/phi = " << eta_ith_tower << "\t" << phi_ith_tower << "\t" << calMuon.getE() << " GeV" << endl;
1443 }
1444 //else cout << "eta/phi = " << eta_ith_tower << "\t" << phi_ith_tower << "\tnot active\n";
1445 }
1446
1447 } // NxN grid
1448 //cout << "-----------" << etrat << endl;
1449 }
1450 else if (CEN_max_mu < CEN_max_calo_fwd)
1451 cout << "** ERROR in RESOLution::CaloIsolation: 'muon'-tower not found! **" << endl;
1452 // should never happen ! this would be a bug
1453 //cout << "etrat = " << et_sum << "\t Pt=" << part.Pt() << endl;
1454
1455 // should return a number between 0 and 99 (due to LHCO definitions)
1456 // which is Pt(muon) / sum(ET)
1457 float etrat = 0.;
1458 if(et_sum==0) etrat = 99.;
1459 else if(et_sum>0) etrat = 100*part.Pt()/et_sum;
1460 if(etrat<0) cout << "Error: negative etrat in CaloIsolation (" << etrat <<")\n";
1461 //else if(etrat>99) cout << "Error: etrat should be in [0;99] in CaloIsolation (" << etrat <<")\n";
1462 if(etrat>99) etrat = 99;
1463 return etrat;
1464}
1465
1466
1467 //********** returns a segmented value for eta and phi, for calo towers *****
1468void RESOLution::BinEtaPhi(const float phi, const float eta, float& iPhi, float& iEta){
1469 iEta = UNDEFINED;
1470 int index= iUNDEFINED;
1471 for (unsigned int i=1; i< TOWER_number+1; i++) {
1472 if(fabs(eta)>TOWER_eta_edges[i-1] && fabs(eta)<=TOWER_eta_edges[i]) {
1473 iEta = (eta>0) ? ((TOWER_eta_edges[i-1]+TOWER_eta_edges[i])/2.0) : -((TOWER_eta_edges[i-1]+TOWER_eta_edges[i])/2.0);
1474 index = i-1;
1475 break;
1476 }
1477 }
1478 if(index==UNDEFINED) return;
1479 iPhi = UNDEFINED;
1480 float dphi = TOWER_dphi[index]*pi/180.;
1481 for (unsigned int i=1; i < 360/TOWER_dphi[index]; i++ ) {
1482 float low = -pi+(i-1)*dphi;
1483 float high= low+dphi;
1484 if(phi > low && phi <= high ){
1485 iPhi = (low+high)/2.0;
1486 break;
1487 }
1488 }
1489 if (phi > pi-dphi) iPhi = pi-dphi;
1490}
1491
1492//**************************** Returns the delta Phi ****************************
1493float DeltaPhi(const float phi1, const float phi2) {
1494 float deltaphi=phi1-phi2; // in here, -pi < phi < pi
1495 if(fabs(deltaphi) > pi) {
1496 deltaphi=2.*pi -fabs(deltaphi);// put deltaphi between 0 and pi
1497 }
1498 else deltaphi=fabs(deltaphi);
1499
1500 return deltaphi;
1501}
1502
1503//**************************** Returns the delta R****************************
1504float DeltaR(const float phi1, const float eta1, const float phi2, const float eta2) {
1505 return sqrt(pow(DeltaPhi(phi1,phi2),2) + pow(eta1-eta2,2));
1506}
1507
1508int sign(const int myint) {
1509 if (myint >0) return 1;
1510 else if (myint <0) return -1;
1511 else return 0;
1512}
1513
1514int sign(const float myfloat) {
1515 if (myfloat >0) return 1;
1516 else if (myfloat <0) return -1;
1517 else return 0;
1518}
1519
1520int ChargeVal(const int pid)
1521{
1522 cout << "ChargeVal :: deprecated function, do not use it anymore" << endl;
1523 int charge;
1524 if(
1525 (pid == pGAMMA) ||
1526 (pid == pPI0) ||
1527 (pid == pK0L) ||
1528 (pid == pN) ||
1529 (pid == pSIGMA0) ||
1530 (pid == pDELTA0) ||
1531 (pid == pK0S) // not charged particles : invisible by tracker
1532 )
1533 charge = 0;
1534 else charge = sign(pid);
1535 return charge;
1536
1537}
1538
1539//------------------------------------------------------------------------------
1540void RESOLution::ReadParticleDataGroupTable() {
1541
1542 string temp_string;
1543 istringstream curstring;
1544
1545 ifstream fichier_a_lire(PdgTableFilename.c_str());
1546 if(!fichier_a_lire.good()) {
1547 cout <<"** ERROR: PDG Table ("<< PdgTableFilename
1548 << ") not found! exit. **" << endl;
1549 exit(1);
1550 return;
1551 }
1552 // first three lines of the file are useless
1553 getline(fichier_a_lire,temp_string);
1554 getline(fichier_a_lire,temp_string);
1555 getline(fichier_a_lire,temp_string);
1556
1557
1558 while (getline(fichier_a_lire,temp_string)) {
1559 curstring.clear(); // needed when using several times istringstream::str(string)
1560 curstring.str(temp_string);
1561 long int ID; std::string name; int charge; float mass; float width; float lifetime;
1562 // ID name chg mass total width lifetime
1563 // 1 d -1 0.33000 0.00000 0.00000E+00
1564 // in the table, the charge is in units of e+/3
1565 // the total width is in GeV
1566 // the lifetime is ctau in mm
1567 curstring >> ID >> name >> charge >> mass >> width >> lifetime;
1568 PdgParticle particle(ID,name,mass,charge/3.,width,lifetime/1000.);
1569 PDGtable.insert(ID,particle);
1570 //PdgTable.insert(pair<int,PdgParticle>(ID,particle));
1571 //cout << PDGtable[ID].name() << "\t" << PDGtable[ID].mass() << "\t" << PDGtable[ID].charge() << endl;
1572 }
1573
1574} // ReadParticleDataGroupTable
1575
1576
1577// to be improved in order to avoid code repetition
1578// sorry, no time to do it right now (XR, 19/05/2009)
1579void time_report(const TStopwatch& global,const TStopwatch& loop,const TStopwatch& trigger,const TStopwatch& frog,const TStopwatch& lhco, const int flag_frog, const int flag_trigger, const int flag_lhco, const string& LogName, const Long64_t allEntries) {
1580
1581TStopwatch globalwatch(global), loopwatch(loop), triggerwatch(trigger), frogwatch(frog), lhcowatch(lhco);
1582
1583 cout <<"** **"<< endl;
1584 cout <<"** ################## Time report ################# **"<< endl;
1585 cout << left << setw(32) <<"** Time report for "<<""
1586 << left << setw(15) << allEntries <<""
1587 << right << setw(22) <<"events **"<<endl;
1588 cout <<"** **"<< endl;
1589 cout << left << setw(10) <<"**"<<""
1590 << left << setw(15) <<"Time (s):"<<""
1591 << right << setw(15) <<"CPU"<<""
1592 << right << setw(15) <<"Real"<<""
1593 << right << setw(14) <<"**"<<endl;
1594 cout << left << setw(10) <<"**"<<""
1595 << left << setw(15) <<" + Global:"<<""
1596 << right << setw(15) <<globalwatch.CpuTime()<<""
1597 << right << setw(15) <<globalwatch.RealTime()<<""
1598 << right << setw(14) <<"**"<<endl;
1599 cout << left << setw(10) <<"**"<<""
1600 << left << setw(15) <<" + Events:"<<""
1601 << right << setw(15) <<loopwatch.CpuTime()<<""
1602 << right << setw(15) <<loopwatch.RealTime()<<""
1603 << right << setw(14) <<"**"<<endl;
1604 if(flag_trigger == 1)
1605 {
1606 cout << left << setw(10) <<"**"<<""
1607 << left << setw(15) <<" + Trigger:"<<""
1608 << right << setw(15) <<triggerwatch.CpuTime()<<""
1609 << right << setw(15) <<triggerwatch.RealTime()<<""
1610 << right << setw(14) <<"**"<<endl;
1611 }
1612 if(flag_frog == 1)
1613 {
1614 cout << left << setw(10) <<"**"<<""
1615 << left << setw(15) <<" + Frog:"<<""
1616 << right << setw(15) <<frogwatch.CpuTime()<<""
1617 << right << setw(15) <<frogwatch.RealTime()<<""
1618 << right << setw(14) <<"**"<<endl;
1619 }
1620 if(flag_lhco == 1)
1621 {
1622 cout << left << setw(10) <<"**"<<""
1623 << left << setw(15) <<" + LHCO:"<<""
1624 << right << setw(15) <<lhcowatch.CpuTime()<<""
1625 << right << setw(15) <<lhcowatch.RealTime()<<""
1626 << right << setw(14) <<"**"<<endl;
1627 }
1628
1629
1630 ofstream f_out(LogName.c_str(),ios_base::app);
1631
1632 f_out <<"** *"<< endl;
1633 f_out <<"** ################## Time report ################# *"<< endl;
1634 f_out << left << setw(32) <<"** Time report for "<<""
1635 << left << setw(15) << allEntries <<""
1636 << right << setw(23) <<"events *"<<endl;
1637 f_out <<"** *"<< endl;
1638 f_out << left << setw(10) <<"**"<<""
1639 << left << setw(15) <<"Time (s):"<<""
1640 << right << setw(15) <<"CPU"<<""
1641 << right << setw(15) <<"Real"<<""
1642 << right << setw(15) <<" *"<<endl;
1643 f_out << left << setw(10) <<"**"<<""
1644 << left << setw(15) <<" + Global:"<<""
1645 << right << setw(15) <<globalwatch.CpuTime()<<""
1646 << right << setw(15) <<globalwatch.RealTime()<<""
1647 << right << setw(15) <<" *"<<endl;
1648 f_out << left << setw(10) <<"**"<<""
1649 << left << setw(15) <<" + Events:"<<""
1650 << right << setw(15) <<loopwatch.CpuTime()<<""
1651 << right << setw(15) <<loopwatch.RealTime()<<""
1652 << right << setw(15) <<" *"<<endl;
1653 if(flag_trigger == 1)
1654 {
1655 f_out << left << setw(10) <<"**"<<""
1656 << left << setw(15) <<" + Trigger:"<<""
1657 << right << setw(15) <<triggerwatch.CpuTime()<<""
1658 << right << setw(15) <<triggerwatch.RealTime()<<""
1659 << right << setw(15) <<" *"<<endl;
1660 }
1661 if(flag_frog == 1)
1662 {
1663 f_out << left << setw(10) <<"**"<<""
1664 << left << setw(15) <<" + Frog:"<<""
1665 << right << setw(15) <<frogwatch.CpuTime()<<""
1666 << right << setw(15) <<frogwatch.RealTime()<<""
1667 << right << setw(15) <<" *"<<endl;
1668 }
1669 if(flag_lhco == 1)
1670 {
1671 f_out << left << setw(10) <<"**"<<""
1672 << left << setw(15) <<" + LHCO:"<<""
1673 << right << setw(15) <<lhcowatch.CpuTime()<<""
1674 << right << setw(15) <<lhcowatch.RealTime()<<""
1675 << right << setw(15) <<" *"<<endl;
1676 }
1677 f_out << left << setw(16) << "** " << ""
1678 << left << setw(52) << get_time_date() << "**" << endl;
1679
1680
1681 f_out<<"* *"<<"\n";
1682 f_out<<"* *"<<"\n";
1683 f_out<<"#....................................................................*"<<"\n";
1684 f_out<<"#>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"<<"\n";
1685
1686 f_out.close();
1687
1688}
1689
1690void print_header() {
1691 cout << endl << endl;
1692
1693 cout <<"*********************************************************************"<< endl;
1694 cout <<"*********************************************************************"<< endl;
1695 cout <<"** **"<< endl;
1696 cout <<"** Welcome to **"<< endl;
1697 cout <<"** **"<< endl;
1698 cout <<"** **"<< endl;
1699 cout <<"** .ddddddd- lL hH **"<< endl;
1700 cout <<"** -Dd` `dD: Ll hH` **"<< endl;
1701 cout <<"** dDd dDd eeee. lL .pp+pp Hh+hhh` -eeee- `sssss **"<< endl;
1702 cout <<"** -Dd `DD ee. ee Ll .Pp. PP Hh. HH. ee. ee sSs **"<< endl;
1703 cout <<"** dD` dDd eEeee: lL. pP. pP hH hH` eEeee:` -sSSSs. **"<< endl;
1704 cout <<"** .Dd :dd eE. LlL PpppPP Hh Hh eE sSS **"<< endl;
1705 cout <<"** dddddd:. eee+: lL. pp. hh. hh eee+ sssssS **"<< endl;
1706 cout <<"** Pp **"<< endl;
1707 cout <<"** **"<< endl;
1708 cout <<"** Delphes, a framework for the fast simulation **"<< endl;
1709 cout <<"** of a generic collider experiment **"<< endl;
1710 cout <<"** arXiv:0903.2225v1 [hep-ph] **"<< endl;
1711 cout <<"** **"<< endl;
1712 cout <<"** --- Version 1.7 of Delphes --- **"<< endl;
1713 cout <<"** Last date of change: 7 May 2009 **"<< endl;
1714 cout <<"** **"<< endl;
1715 cout <<"** **"<< endl;
1716 cout <<"** This package uses: **"<< endl;
1717 cout <<"** ------------------ **"<< endl;
1718 cout <<"** FastJet algorithm: Phys. Lett. B641 (2006) [hep-ph/0512210] **"<< endl;
1719 cout <<"** Hector: JINST 2:P09005 (2007) [physics.acc-ph:0707.1198v2] **"<< endl;
1720 cout <<"** FROG: [hep-ex/0901.2718v1] **"<< endl;
1721 cout <<"** **"<< endl;
1722 cout <<"**-----------------------------------------------------------------**"<< endl;
1723 cout <<"** **"<< endl;
1724 cout <<"** Main authors: **"<< endl;
1725 cout <<"** ------------- **"<< endl;
1726 cout <<"** **"<< endl;
1727 cout <<"** Séverine Ovyn Xavier Rouby **"<< endl;
1728 cout <<"** severine.ovyn@uclouvain.be xavier.rouby@cern **"<< endl;
1729 cout <<"** Center for Particle Physics and Phenomenology (CP3) **"<< endl;
1730 cout <<"** Universite Catholique de Louvain (UCL) **"<< endl;
1731 cout <<"** Louvain-la-Neuve, Belgium **"<< endl;
1732 cout <<"** **"<< endl;
1733 cout <<"**-----------------------------------------------------------------**"<< endl;
1734 cout <<"** **"<< endl;
1735 cout <<"** Former Delphes versions and documentation can be found on : **"<< endl;
1736 cout <<"** http://www.fynu.ucl.ac.be/delphes.html **"<< endl;
1737 cout <<"** **"<< endl;
1738 cout <<"** **"<< endl;
1739 cout <<"** Disclaimer: this program is a beta version of Delphes and **"<< endl;
1740 cout <<"** therefore comes without guarantees. Beware of errors and please **"<< endl;
1741 cout <<"** give us your feedbacks about potential bugs **"<< endl;
1742 cout <<"** **"<< endl;
1743 cout <<"*********************************************************************"<< endl;
1744 cout <<"*********************************************************************"<< endl;
1745
1746}
1747
1748string get_time_date() {
1749 time_t rawtime;
1750 struct tm * timeinfo;
1751
1752 time ( &rawtime );
1753 timeinfo = localtime ( &rawtime );
1754
1755 char temp[100];
1756 sprintf(temp,"%i/%i/%i %i:%i:%i",timeinfo->tm_mday,timeinfo->tm_mon+1,timeinfo->tm_year+1900,timeinfo->tm_hour,timeinfo->tm_min,timeinfo->tm_sec);
1757 string tempstring(temp);
1758 return tempstring;
1759}
1760
1761void warning(const string oldname, const string newname) {
1762 string text = "** Warning in datacard: " + oldname + " deprecated. Use " + newname + " instead.";
1763 cout << left << setw(67) << text
1764 << right << setw(2) <<"**" << endl;
1765}
Note: See TracBrowser for help on using the repository browser.