ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CylinderGeomMicromegas.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file CylinderGeomMicromegas.cc
1 
7 
8 #include <g4main/PHG4Hit.h>
9 
10 #include <TVector3.h>
11 
12 #include <cassert>
13 
14 namespace
15 {
16  // convenient square function
17  template<class T>
18  inline constexpr T square( const T& x ) { return x*x; }
19 
21  template<class T>
22  inline constexpr T get_r( T x, T y ) { return std::sqrt( square(x) + square(y) ); }
23 
24  // bind angle to [-M_PI,+M_PI[. This is useful to avoid edge effects when making the difference between two angles
25  template<class T>
26  inline T bind_angle( const T& angle )
27  {
28  if( angle >= M_PI ) return angle - 2*M_PI;
29  else if( angle < -M_PI ) return angle + 2*M_PI;
30  else return angle;
31  }
32 
33 }
34 
35 //________________________________________________________________________________
36 TVector3 CylinderGeomMicromegas::get_local_from_world_coords( uint tileid, const TVector3& world_coordinates ) const
37 {
38  assert( tileid < m_tiles.size() );
39 
40  // store world coordinates in array
41  std::array<double,3> master;
42  world_coordinates.GetXYZ( &master[0] );
43 
44  // convert to local coordinate
45  std::array<double,3> local;
46  transformation_matrix(tileid).MasterToLocal( &master[0], &local[0] );
47 
48  return TVector3( &local[0] );
49 
50 }
51 
52 //________________________________________________________________________________
53 TVector3 CylinderGeomMicromegas::get_world_from_local_coords( uint tileid, const TVector3& local_coordinates ) const
54 {
55  assert( tileid < m_tiles.size() );
56 
57  // store world coordinates in array
58  std::array<double,3> local;
59  local_coordinates.GetXYZ( &local[0] );
60 
61  // convert to local coordinate
62  std::array<double,3> master;
63  transformation_matrix(tileid).LocalToMaster( &local[0], &master[0] );
64 
65  return TVector3( &master[0] );
66 
67 }
68 
69 //________________________________________________________________________________
70 TVector3 CylinderGeomMicromegas::get_world_from_local_vect( uint tileid, const TVector3& local_vect ) const
71 {
72  assert( tileid < m_tiles.size() );
73 
74  // store world vect in array
75  std::array<double,3> local;
76  local_vect.GetXYZ( &local[0] );
77 
78  // convert to local coordinate
79  std::array<double,3> master;
80  transformation_matrix(tileid).LocalToMasterVect( &local[0], &master[0] );
81 
82  return TVector3( &master[0] );
83 
84 }
85 
86 //________________________________________________________________________________
87 TVector3 CylinderGeomMicromegas::get_local_from_world_vect( uint tileid, const TVector3& world_vect ) const
88 {
89  assert( tileid < m_tiles.size() );
90 
91  // store world vect in array
92  std::array<double,3> master;
93  world_vect.GetXYZ( &master[0] );
94 
95  // convert to local coordinate
96  std::array<double,3> local;
97  transformation_matrix(tileid).MasterToLocalVect( &master[0], &local[0] );
98 
99  return TVector3( &local[0] );
100 
101 }
102 
103 //________________________________________________________________________________
104 int CylinderGeomMicromegas::find_tile_cylindrical( const TVector3& world_coordinates ) const
105 {
106  // check radius
107  if( !check_radius(world_coordinates) ) return -1;
108 
109  // convert to polar coordinates
110  const auto phi = std::atan2( world_coordinates.y(), world_coordinates.x() );
111  const auto z = world_coordinates.z();
112 
113  for( size_t itile = 0; itile < m_tiles.size(); ++itile )
114  {
115  const auto& tile = m_tiles.at(itile);
116 
117  if( std::abs( z - tile.m_centerZ ) > tile.m_sizeZ/2 ) continue;
118  if( std::abs( bind_angle( phi - tile.m_centerPhi ) ) > tile.m_sizePhi/2 ) continue;
119 
120  return itile;
121  }
122 
123  return -1;
124 }
125 
126 //________________________________________________________________________________
127 int CylinderGeomMicromegas::find_tile_planar( const TVector3& world_coordinates ) const
128 {
129 
130  for( size_t itile = 0; itile < m_tiles.size(); ++itile )
131  {
132 
133  // get local coordinates
134  const auto local_coordinates = get_local_from_world_coords( itile, world_coordinates );
135 
136  // store tile struct locally
137  const auto& tile = m_tiles.at(itile);
138 
139  // check against thickness
140  if( std::abs( local_coordinates.y() ) > m_thickness/2 ) continue;
141 
142  // check azimuth
143  if( std::abs( local_coordinates.x() ) > tile.m_sizePhi*reference_radius/2 ) continue;
144 
145  // check z extend
146  if( std::abs( local_coordinates.z() ) > tile.m_sizeZ/2 ) continue;
147 
148  return itile;
149  }
150 
151  return -1;
152 }
153 
154 //________________________________________________________________________________
155 void CylinderGeomMicromegas::convert_to_planar( uint tileid, PHG4Hit* hit ) const
156 {
157 
158  // get extrapolation direction based on in and out G4Hit positions
159  const TVector3 world_direction( hit->get_x(1)-hit->get_x(0), hit->get_y(1)-hit->get_y(0), hit->get_z(1)-hit->get_z(0) );
160 
161  // same treatment is applied to both coordinates of the G4Hit
162  for( int i = 0; i < 2; ++i )
163  {
164 
165  const TVector3 world_coordinates( hit->get_x(i), hit->get_y(i), hit->get_z(i) );
166 
167  // get coordinate radius, and compare to reference radius. This will be use to offset local position along y in local space
168  const double delta_radius = get_r( world_coordinates.x(), world_coordinates.y() ) - m_radius;
169 
170  // convert position and direction to local coordinates
171  const auto local_coordinates = get_local_from_world_coords( tileid, world_coordinates );
172  const auto local_direction = get_local_from_world_vect( tileid, world_direction );
173 
174  // propagate along local_direction to the rigth y
175  const double delta_y = delta_radius - local_coordinates.y();
176  TVector3 new_local_coordinates(
177  local_coordinates.x() + delta_y*local_direction.x()/local_direction.y(),
178  local_coordinates.y() + delta_y,
179  local_coordinates.z() + delta_y*local_direction.z()/local_direction.y() );
180 
181  // convert back to world coordinate and assign to g4hit
182  const auto new_world_coordinates = get_world_from_local_coords( tileid, new_local_coordinates );
183  hit->set_x(i, new_world_coordinates.x() );
184  hit->set_y(i, new_world_coordinates.y() );
185  hit->set_z(i, new_world_coordinates.z() );
186  }
187 }
188 
189 //________________________________________________________________________________
190 int CylinderGeomMicromegas::find_strip_from_world_coords( uint tileid, const TVector3& world_coordinates ) const
191 {
192  // convert to local coordinates
193  const auto local_coordinates = get_local_from_world_coords( tileid, world_coordinates );
194  return find_strip_from_local_coords( tileid, local_coordinates );
195 }
196 
197 //________________________________________________________________________________
198 int CylinderGeomMicromegas::find_strip_from_local_coords( uint tileid, const TVector3& local_coordinates ) const
199 {
200 
201  // check against thickness
202  if( std::abs( local_coordinates.y() ) > m_thickness/2 ) return -1;
203 
204  // get tile
205  const auto& tile = m_tiles.at(tileid);
206 
207  // check azimuth
208  if( std::abs( local_coordinates.x() ) > tile.m_sizePhi*reference_radius/2 ) return -1;
209 
210  // check z extend
211  if( std::abs( local_coordinates.z() ) > tile.m_sizeZ/2 ) return -1;
212 
213  // we found a tile to which the hit belong
214  // calculate strip index, depending on cylinder direction
215  switch( m_segmentation_type )
216  {
218  return (int) std::floor( (local_coordinates.x() + tile.m_sizePhi*reference_radius/2)/m_pitch );
219 
221  return (int) std::floor( (local_coordinates.z() + tile.m_sizeZ/2)/m_pitch );
222  }
223 
224  // unreachable
225  return -1;
226 }
227 
228 //________________________________________________________________________________
229 double CylinderGeomMicromegas::get_strip_length( uint tileid ) const
230 {
231  assert( tileid < m_tiles.size() );
232  switch( m_segmentation_type )
233  {
235  return m_tiles[tileid].m_sizeZ;
236 
238  return m_tiles[tileid].m_sizePhi*reference_radius;
239  }
240 
241  // unreachable
242  return 0;
243 }
244 
245 //________________________________________________________________________________
247 {
248  assert( tileid < m_tiles.size() );
249  switch( m_segmentation_type )
250  {
252  return std::floor( m_tiles[tileid].m_sizePhi*reference_radius/m_pitch );
253 
255  return std::floor( m_tiles[tileid].m_sizeZ/m_pitch );
256  }
257 
258  // unreachable
259  return 0;
260 }
261 
262 //________________________________________________________________________________
263 TVector3 CylinderGeomMicromegas::get_local_coordinates( uint tileid, uint stripnum ) const
264 {
265  assert( tileid < m_tiles.size() );
266 
267  // get tile
268  const auto& tile = m_tiles[tileid];
269 
270  switch( m_segmentation_type )
271  {
273  return TVector3( (0.5+stripnum)*m_pitch - tile.m_sizePhi*reference_radius/2, 0, 0 );
274 
276  return TVector3( 0, 0, (0.5+stripnum)*m_pitch - tile.m_sizeZ/2 );
277  }
278 
279  // unreachable
280  return TVector3();
281 }
282 
283 //________________________________________________________________________________
284 TVector3 CylinderGeomMicromegas::get_world_coordinates( uint tileid, uint stripnum ) const
285 {
286  assert( tileid < m_tiles.size() );
287  return get_world_from_local_coords( tileid, get_local_coordinates( tileid, stripnum ) );
288 }
289 
290 //________________________________________________________________________________
291 void CylinderGeomMicromegas::identify( std::ostream& out ) const
292 {
293  out << "CylinderGeomMicromegas" << std::endl;
294  out << "layer: " << m_layer << std::endl;
295  out << "segmentation_type: " << (m_segmentation_type == MicromegasDefs::SegmentationType::SEGMENTATION_PHI ? "SEGMENTATION_PHI":"SEGMENTATION_Z") << std::endl;
296  out << "drift_direction: " << (m_drift_direction == MicromegasDefs::DriftDirection::INWARD ? "INWARD":"OUTWARD") << std::endl;
297  out << "radius: " << m_radius << "cm" << std::endl;
298  out << "thickness: " << m_thickness << "cm" << std::endl;
299  out << "zmin: " << m_zmin << "cm" << std::endl;
300  out << "zmax: " << m_zmax << "cm" << std::endl;
301  out << "pitch: " << m_pitch << "cm" << std::endl;
302  out << std::endl;
303 }
304 
305 //________________________________________________________________________________
306 bool CylinderGeomMicromegas::check_radius( const TVector3& world_coordinates ) const
307 {
308  const auto radius = get_r( world_coordinates.x(), world_coordinates.y() );
309  return std::abs( radius - m_radius ) <= m_thickness/2;
310 }
311 
312 //________________________________________________________________________________
313 TGeoHMatrix CylinderGeomMicromegas::transformation_matrix( uint tileid ) const
314 {
315  assert( tileid < m_tiles.size() );
316  const auto& tile = m_tiles[tileid];
317 
318  // local to master transformation matrix
319  TGeoHMatrix matrix;
320 
321  // rotate around z axis
322  matrix.RotateZ( tile.m_centerPhi*180./M_PI - 90 );
323 
324  // translate to tile center
325  matrix.SetDx( m_radius*std::cos( tile.m_centerPhi ) );
326  matrix.SetDy( m_radius*std::sin( tile.m_centerPhi ) );
327  matrix.SetDz( tile.m_centerZ );
328 
329  return matrix;
330 }