test
CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
RecoProducerFP420.cc
Go to the documentation of this file.
3 
7 
8 #include "CLHEP/Vector/LorentzVector.h"
9 
10 #include <math.h>
11 
12 //#include <iostream>
13 
15 
16  // Create LHC beam line
17  // double length = param.getParameter<double>("BeamLineLength");
18  // std::cout << " BeamLineLength = " << length << std::endl;
19 
20  length = conf_.getParameter<double>("BeamLineLength" );//m
21  verbosity = conf_.getUntrackedParameter<int>("VerbosityLevel");
24 
25  edm::LogInfo ("RecoProducerFP420") << "RecoProducerFP420 parameters: \n"
26  << " Verbosity: " << verbosity << "\n"
27  << " length: " << length << "\n";
28  if (verbosity > 1) {
29  std::cout << " RecoProducerFP420: constructor " << std::endl;
30  std::cout << " BeamLineLength: " << length << std::endl;
31  }
32 
33  // edm::FileInPath b1("SimTransport/HectorData/twiss_ip5_b1_v6.5.txt");
34  // edm::FileInPath b2("SimTransport/HectorData/twiss_ip5_b2_v6.5.txt");
35 
36  edm::FileInPath b1(beam1filename.c_str());
37  edm::FileInPath b2(beam2filename.c_str());
38 
39 
40  m_beamline1 = new H_BeamLine( 1, length + 0.1 ); // (direction, length)
41  m_beamline2 = new H_BeamLine( -1, length + 0.1 ); //
42 
43  try {
44  m_beamline1->fill( b1.fullPath(), 1, "IP5" );
45  m_beamline2->fill( b2.fullPath(), -1, "IP5" );
46  } catch ( const edm::Exception& e ) {
47  std::string msg = e.what();
48  msg += " caught in RecoProducerFP420... \nERROR: Could not locate SimTransport/HectorData data files.";
49  edm::LogError ("DataNotFound") << msg;
50  }
51 
52  m_beamline1->offsetElements( 120, -0.097 );
53  m_beamline2->offsetElements( 120, +0.097 );
54 
55  m_beamline1->calcMatrix();
56  m_beamline2->calcMatrix();
57 
58  edm::LogInfo ("RecoProducerFP420") << "==============================\n";
59 
60 }
61 
63 
64 std::vector<RecoFP420> RecoProducerFP420::reconstruct(int direction, double x1_420, double y1_420, double x2_420, double y2_420, double z1_420, double z2_420){
65  // ==================================================================================
66  // ==================================================================================
67  std::vector<RecoFP420> rhits;
68  int restracks = 10;// max # tracks
69  rhits.reserve(restracks);
70  rhits.clear();
71  // ==================================================================================
72 // trivial (TM), angle compensation (AM) and position compensation (PM) methods
73 // #define TM 1 #define AM 2 #define PM 3
74  m_tx0=-100., m_ty0=-100., m_x0=-100., m_y0=-100.;
75  if ( direction == 1 ) {
76  m_rp420_f = new H_RecRPObject( z1_420*0.001, z2_420*0.001, *m_beamline1 );// m
77  if ( m_rp420_f ) {
78  if (verbosity >1) {
79  std::cout << " RecoProducerFP420: input coord. in um " << std::endl;
80  std::cout << " x1_420: " << x1_420 << " y1_420: " << y1_420 << std::endl;
81  std::cout << " x2_420: " << x2_420 << " y2_420: " << y2_420 << std::endl;
82  }
83  m_rp420_f->setPositions( x1_420, y1_420, x2_420 ,y2_420 );//input coord. in um
84  m_e = m_rp420_f->getE( AM );// GeV
85  // std::cout << " m_e1: " << m_rp420_f->getE( TM ) << std::endl;
86  // std::cout << " m_e2: " << m_rp420_f->getE( AM ) << std::endl;
87  m_tx0 = m_rp420_f->getTXIP();// urad
88  m_ty0 = m_rp420_f->getTYIP();// urad
89  m_x0 = m_rp420_f->getX0();// um
90  m_y0 = m_rp420_f->getY0();// um
91  m_q2 = m_rp420_f->getQ2();// GeV^2
92  }// if ( m_rp420_f
93  }// if ( dire
94  else if ( direction == 2 ) {
95  m_rp420_b = new H_RecRPObject( z1_420*0.001, z2_420*0.001, *m_beamline2 );// m
96  if ( m_rp420_b ) {
97  m_rp420_b->setPositions( x1_420, y1_420, x2_420 ,y2_420 );// input coord. in um
98  m_e = m_rp420_b->getE( AM );// GeV
99  m_tx0 = m_rp420_b->getTXIP();// urad
100  m_ty0 = m_rp420_b->getTYIP();// urad
101  m_x0 = m_rp420_b->getX0();// um
102  m_y0 = m_rp420_b->getY0();// um
103  m_q2 = m_rp420_b->getQ2();// GeV^2
104  }// if ( m_rp420_b
105  }
106  else{
107  return rhits;
108  }
109 
110  // ==============================
111  if (verbosity > 1) {
112  std::cout << " RecoProducerFP420: rhits.push_back " << std::endl;
113  std::cout << " m_e: " << m_e << std::endl;
114  std::cout << " m_x0: " << m_x0 << std::endl;
115  std::cout << " m_y0: " << m_y0 << std::endl;
116  std::cout << " m_tx0: " << m_tx0 << std::endl;
117  std::cout << " m_ty0: " << m_ty0 << std::endl;
118  std::cout << " m_q2: " << m_q2 << std::endl;
119  std::cout << " direction: " << direction << std::endl;
120  }
121  rhits.push_back( RecoFP420(m_e,m_x0,m_y0,m_tx0,m_ty0,m_q2,direction) );
122  // ==============================
124  return rhits;
125  //============
126 }
127 
128 
129 
virtual char const * what() const
Definition: Exception.cc:141
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
std::vector< RecoFP420 > rhits
std::string beam1filename
H_BeamLine * m_beamline2
edm::ParameterSet conf_
RecoProducerFP420(const edm::ParameterSet &conf)
H_RecRPObject * m_rp420_f
std::string beam2filename
std::vector< RecoFP420 > reconstruct(int, double, double, double, double, double, double)
H_RecRPObject * m_rp420_b
tuple cout
Definition: gather_cfg.py:145
H_BeamLine * m_beamline1
virtual ~RecoProducerFP420()