ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ParameterisationTrd.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4ParameterisationTrd.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // G4ParameterisationTrd[X/Y/Z] implementation
27 //
28 // 26.05.03 - P.Arce, Initial version
29 // 08.04.04 - I.Hrivnacova, Implemented reflection
30 // 21.04.10 - M.Asai, Added gaps
31 // --------------------------------------------------------------------
32 
33 #include "G4ParameterisationTrd.hh"
34 
35 #include <iomanip>
36 #include "G4ThreeVector.hh"
37 #include "G4RotationMatrix.hh"
38 #include "G4VPhysicalVolume.hh"
39 #include "G4LogicalVolume.hh"
40 #include "G4ReflectedSolid.hh"
41 #include "G4Trd.hh"
42 #include "G4Trap.hh"
43 
44 //--------------------------------------------------------------------------
45 G4VParameterisationTrd::
46 G4VParameterisationTrd( EAxis axis, G4int nDiv, G4double width,
47  G4double offset, G4VSolid* msolid,
48  DivisionType divType )
49  : G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
50 {
51  G4Trd* msol = (G4Trd*)(msolid);
52  if (msolid->GetEntityType() == "G4ReflectedSolid")
53  {
54  // Get constituent solid
55  G4VSolid* mConstituentSolid
56  = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
57  msol = (G4Trd*)(mConstituentSolid);
58 
59  // Create a new solid with inversed parameters
60  G4Trd* newSolid
61  = new G4Trd(msol->GetName(),
62  msol->GetXHalfLength2(), msol->GetXHalfLength1(),
63  msol->GetYHalfLength2(), msol->GetYHalfLength1(),
64  msol->GetZHalfLength());
65  msol = newSolid;
66  fmotherSolid = newSolid;
67  fReflectedSolid = true;
68  fDeleteSolid = true;
69  }
70 }
71 
72 //------------------------------------------------------------------------
73 G4VParameterisationTrd::~G4VParameterisationTrd()
74 {
75 }
76 
77 //------------------------------------------------------------------------
78 G4ParameterisationTrdX::
79 G4ParameterisationTrdX( EAxis axis, G4int nDiv,
80  G4double width, G4double offset,
81  G4VSolid* msolid, DivisionType divType )
82  : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
83 {
84  CheckParametersValidity();
85  SetType( "DivisionTrdX" );
86 
87  G4Trd* msol = (G4Trd*)(fmotherSolid);
88  if( divType == DivWIDTH )
89  {
90  fnDiv = CalculateNDiv( msol->GetXHalfLength1()+msol->GetXHalfLength2(),
91  width, offset );
92  }
93  else if( divType == DivNDIV )
94  {
95  fwidth = CalculateWidth( msol->GetXHalfLength1()+msol->GetXHalfLength2(),
96  nDiv, offset );
97  }
98 
99 #ifdef G4DIVDEBUG
100  if( verbose >= 1 )
101  {
102  G4cout << " G4ParameterisationTrdX - ## divisions " << fnDiv << " = "
103  << nDiv << G4endl
104  << " Offset " << foffset << " = " << offset << G4endl
105  << " Width " << fwidth << " = " << width << G4endl;
106  }
107 #endif
108 
109  G4double mpDx1 = msol->GetXHalfLength1();
110  G4double mpDx2 = msol->GetXHalfLength2();
111  if( std::fabs(mpDx1 - mpDx2) > kCarTolerance )
112  {
113  bDivInTrap = true;
114  }
115 }
116 
117 //------------------------------------------------------------------------
118 G4ParameterisationTrdX::~G4ParameterisationTrdX()
119 {
120 }
121 
122 //------------------------------------------------------------------------
123 G4double G4ParameterisationTrdX::GetMaxParameter() const
124 {
125  G4Trd* msol = (G4Trd*)(fmotherSolid);
126  return (msol->GetXHalfLength1()+msol->GetXHalfLength2());
127 }
128 
129 //------------------------------------------------------------------------
130 void
131 G4ParameterisationTrdX::
132 ComputeTransformation( const G4int copyNo,
133  G4VPhysicalVolume *physVol ) const
134 {
135  G4Trd* msol = (G4Trd*)(fmotherSolid );
136  G4double mdx = ( msol->GetXHalfLength1() + msol->GetXHalfLength2() ) / 2.;
137 
138  //----- translation
139  G4ThreeVector origin(0.,0.,0.);
140  G4double posi;
141  if( !bDivInTrap )
142  {
143  posi = -mdx + foffset + (copyNo+0.5)*fwidth;
144  }
145  else
146  {
147  G4double aveHL = (msol->GetXHalfLength1()+msol->GetXHalfLength2())/2.;
148  posi = - aveHL + foffset + (copyNo+0.5)*aveHL/fnDiv*2;
149  }
150  if( faxis == kXAxis )
151  {
152  origin.setX( posi );
153  }
154  else
155  {
156  std::ostringstream message;
157  message << "Only axes along X are allowed ! Axis: " << faxis;
158  G4Exception("G4ParameterisationTrdX::ComputeTransformation()",
159  "GeomDiv0002", FatalException, message);
160  }
161 
162 #ifdef G4DIVDEBUG
163  if( verbose >= 2 )
164  {
165  G4cout << std::setprecision(8)
166  << " G4ParameterisationTrdX::ComputeTransformation() "
167  << copyNo << G4endl
168  << " Position: " << origin << " - Axis: " << faxis << G4endl;
169  }
170 #endif
171 
172  //----- set translation
173  physVol->SetTranslation( origin );
174 }
175 
176 //--------------------------------------------------------------------------
177 void
178 G4ParameterisationTrdX::
179 ComputeDimensions( G4Trd& trd, const G4int, const G4VPhysicalVolume* ) const
180 {
181  G4Trd* msol = (G4Trd*)(fmotherSolid);
182 
183  G4double pDy1 = msol->GetYHalfLength1();
184  G4double pDy2 = msol->GetYHalfLength2();
185  G4double pDz = msol->GetZHalfLength();
186  G4double pDx = fwidth/2. - fhgap;
187 
188  trd.SetAllParameters ( pDx, pDx, pDy1, pDy2, pDz );
189 
190 #ifdef G4DIVDEBUG
191  if( verbose >= 2 )
192  {
193  G4cout << " G4ParameterisationTrdX::ComputeDimensions():"
194  << copyNo << G4endl;
195  trd.DumpInfo();
196  }
197 #endif
198 }
199 
200 G4VSolid* G4ParameterisationTrdX::
201 ComputeSolid(const G4int i, G4VPhysicalVolume* pv)
202 {
203  if( bDivInTrap )
204  {
205  return G4VDivisionParameterisation::ComputeSolid(i, pv);
206  }
207  else
208  {
209  return fmotherSolid;
210  }
211 }
212 
213 
214 //--------------------------------------------------------------------------
215 void
216 G4ParameterisationTrdX::ComputeDimensions( G4Trap& trap, const G4int copyNo,
217  const G4VPhysicalVolume* ) const
218 {
219  G4Trd* msol = (G4Trd*)(fmotherSolid);
220  G4double pDy1 = msol->GetYHalfLength1();
221  G4double pDy2 = msol->GetYHalfLength2();
222  G4double pDz = msol->GetZHalfLength();
223  G4double pDx1 = msol->GetXHalfLength1()/fnDiv;
224  G4double pDx2 = msol->GetXHalfLength2()/fnDiv;
225 
226  G4double cxy1 = -msol->GetXHalfLength1() + foffset
227  + (copyNo+0.5)*pDx1*2;// centre of the side at y=-pDy1
228  G4double cxy2 = -msol->GetXHalfLength2() + foffset
229  + (copyNo+0.5)*pDx2*2;// centre of the side at y=+pDy1
230  G4double alp = std::atan( (cxy2-cxy1)/pDz );
231 
232  trap.SetAllParameters ( pDz,
233  0.,
234  0.,
235  pDy1,
236  pDx1,
237  pDx2,
238  alp,
239  pDy2,
240  pDx1 - fhgap,
241  pDx2 - fhgap * pDx2/pDx1,
242  alp);
243 
244 #ifdef G4DIVDEBUG
245  if( verbose >= 2 )
246  {
247  G4cout << " G4ParameterisationTrdX::ComputeDimensions():"
248  << copyNo << G4endl;
249  trap.DumpInfo();
250  }
251 #endif
252 }
253 
254 //--------------------------------------------------------------------------
255 void G4ParameterisationTrdX::CheckParametersValidity()
256 {
257  G4VDivisionParameterisation::CheckParametersValidity();
258 /*
259  G4Trd* msol = (G4Trd*)(fmotherSolid);
260 
261  G4double mpDx1 = msol->GetXHalfLength1();
262  G4double mpDx2 = msol->GetXHalfLength2();
263  bDivInTrap = false;
264 
265  if( std::fabs(mpDx1 - mpDx2) > kCarTolerance )
266  {
267  std::ostringstream message;
268  message << "Invalid solid specification. NOT supported." << G4endl
269  << "Making a division of a TRD along axis X," << G4endl
270  << "while the X half lengths are not equal," << G4endl
271  << "is not (yet) supported. It will result" << G4endl
272  << "in non-equal division solids.";
273  G4Exception("G4ParameterisationTrdX::CheckParametersValidity()",
274  "GeomDiv0001", FatalException, message);
275  }
276 */
277 }
278 
279 //--------------------------------------------------------------------------
280 void G4ParameterisationTrdX::
281 ComputeTrapParams()
282 {
283 }
284 
285 //--------------------------------------------------------------------------
286 G4ParameterisationTrdY::
287 G4ParameterisationTrdY( EAxis axis, G4int nDiv,
288  G4double width, G4double offset,
289  G4VSolid* msolid, DivisionType divType)
290  : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
291 {
292  CheckParametersValidity();
293  SetType( "DivisionTrdY" );
294 
295  G4Trd* msol = (G4Trd*)(fmotherSolid);
296  if( divType == DivWIDTH )
297  {
298  fnDiv = CalculateNDiv( 2*msol->GetYHalfLength1(),
299  width, offset );
300  }
301  else if( divType == DivNDIV )
302  {
303  fwidth = CalculateWidth( 2*msol->GetYHalfLength1(),
304  nDiv, offset );
305  }
306 
307 #ifdef G4DIVDEBUG
308  if( verbose >= 1 )
309  {
310  G4cout << " G4ParameterisationTrdY no divisions " << fnDiv
311  << " = " << nDiv << G4endl
312  << " Offset " << foffset << " = " << offset << G4endl
313  << " width " << fwidth << " = " << width << G4endl;
314  }
315 #endif
316 }
317 
318 //------------------------------------------------------------------------
319 G4ParameterisationTrdY::~G4ParameterisationTrdY()
320 {
321 }
322 
323 //------------------------------------------------------------------------
324 G4double G4ParameterisationTrdY::GetMaxParameter() const
325 {
326  G4Trd* msol = (G4Trd*)(fmotherSolid);
327  return 2*msol->GetYHalfLength1();
328 }
329 
330 //--------------------------------------------------------------------------
331 void
332 G4ParameterisationTrdY::
333 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
334 {
335  G4Trd* msol = (G4Trd*)(fmotherSolid );
336  G4double mdy = msol->GetYHalfLength1();
337 
338  //----- translation
339  G4ThreeVector origin(0.,0.,0.);
340  G4double posi = -mdy + foffset + (copyNo+0.5)*fwidth;
341 
342  if( faxis == kYAxis )
343  {
344  origin.setY( posi );
345  }
346  else
347  {
348  std::ostringstream message;
349  message << "Only axes along Y are allowed ! Axis: " << faxis;
350  G4Exception("G4ParameterisationTrdY::ComputeTransformation()",
351  "GeomDiv0002", FatalException, message);
352  }
353 
354 #ifdef G4DIVDEBUG
355  if( verbose >= 2 )
356  {
357  G4cout << std::setprecision(8)
358  << " G4ParameterisationTrdY::ComputeTransformation " << copyNo
359  << " pos " << origin << " rot mat " << " axis " << faxis << G4endl;
360  }
361 #endif
362 
363  //----- set translation
364  physVol->SetTranslation( origin );
365 }
366 
367 //--------------------------------------------------------------------------
368 void
369 G4ParameterisationTrdY::
370 ComputeDimensions(G4Trd& trd, const G4int, const G4VPhysicalVolume*) const
371 {
372  //---- The division along Y of a Trd will result a Trd, only
373  //--- if Y at -Z and +Z are equal, else use the G4Trap version
374  G4Trd* msol = (G4Trd*)(fmotherSolid);
375 
376  G4double pDx1 = msol->GetXHalfLength1();
377  G4double pDx2 = msol->GetXHalfLength2();
378  G4double pDz = msol->GetZHalfLength();
379  G4double pDy = fwidth/2. - fhgap;
380 
381  trd.SetAllParameters ( pDx1, pDx2, pDy, pDy, pDz );
382 
383 #ifdef G4DIVDEBUG
384  if( verbose >= 2 )
385  {
386  G4cout << " G4ParameterisationTrdY::ComputeDimensions():" << G4endl;
387  trd.DumpInfo();
388  }
389 #endif
390 }
391 
392 //--------------------------------------------------------------------------
393 void G4ParameterisationTrdY::CheckParametersValidity()
394 {
395  G4VDivisionParameterisation::CheckParametersValidity();
396 
397  G4Trd* msol = (G4Trd*)(fmotherSolid);
398 
399  G4double mpDy1 = msol->GetYHalfLength1();
400  G4double mpDy2 = msol->GetYHalfLength2();
401 
402  if( std::fabs(mpDy1 - mpDy2) > kCarTolerance )
403  {
404  std::ostringstream message;
405  message << "Invalid solid specification. NOT supported." << G4endl
406  << "Making a division of a TRD along axis Y while" << G4endl
407  << "the Y half lengths are not equal is not (yet)" << G4endl
408  << "supported. It will result in non-equal" << G4endl
409  << "division solids.";
410  G4Exception("G4ParameterisationTrdY::CheckParametersValidity()",
411  "GeomDiv0001", FatalException, message);
412  }
413 }
414 
415 //--------------------------------------------------------------------------
416 G4ParameterisationTrdZ::
417 G4ParameterisationTrdZ( EAxis axis, G4int nDiv,
418  G4double width, G4double offset,
419  G4VSolid* msolid, DivisionType divType )
420  : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
421 {
422  CheckParametersValidity();
423  SetType( "DivTrdZ" );
424 
425  G4Trd* msol = (G4Trd*)(fmotherSolid);
426  if( divType == DivWIDTH )
427  {
428  fnDiv = CalculateNDiv( 2*msol->GetZHalfLength(),
429  width, offset );
430  }
431  else if( divType == DivNDIV )
432  {
433  fwidth = CalculateWidth( 2*msol->GetZHalfLength(),
434  nDiv, offset );
435  }
436 
437 #ifdef G4DIVDEBUG
438  if( verbose >= 1 )
439  {
440  G4cout << " G4ParameterisationTrdZ no divisions " << fnDiv
441  << " = " << nDiv << G4endl
442  << " Offset " << foffset << " = " << offset << G4endl
443  << " Width " << fwidth << " = " << width << G4endl;
444  }
445 #endif
446 }
447 
448 //------------------------------------------------------------------------
449 G4ParameterisationTrdZ::~G4ParameterisationTrdZ()
450 {
451 }
452 
453 //------------------------------------------------------------------------
454 G4double G4ParameterisationTrdZ::GetMaxParameter() const
455 {
456  G4Trd* msol = (G4Trd*)(fmotherSolid);
457  return 2*msol->GetZHalfLength();
458 }
459 
460 //--------------------------------------------------------------------------
461 void
462 G4ParameterisationTrdZ::
463 ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
464 {
465  G4Trd* msol = (G4Trd*)(fmotherSolid );
466  G4double mdz = msol->GetZHalfLength();
467 
468  //----- translation
469  G4ThreeVector origin(0.,0.,0.);
470  G4double posi = -mdz + OffsetZ() + (copyNo+0.5)*fwidth;
471  if( faxis == kZAxis )
472  {
473  origin.setZ( posi );
474  }
475  else
476  {
477  std::ostringstream message;
478  message << "Only axes along Z are allowed ! Axis: " << faxis;
479  G4Exception("G4ParameterisationTrdZ::ComputeTransformation()",
480  "GeomDiv0002", FatalException, message);
481  }
482 
483 #ifdef G4DIVDEBUG
484  if( verbose >= 1 )
485  {
486  G4cout << std::setprecision(8) << " G4ParameterisationTrdZ: "
487  << copyNo << G4endl
488  << " Position: " << origin << " - Offset: " << foffset
489  << " - Width: " << fwidth << " Axis " << faxis << G4endl;
490  }
491 #endif
492 
493  //----- set translation
494  physVol->SetTranslation( origin );
495 }
496 
497 //--------------------------------------------------------------------------
498 void
499 G4ParameterisationTrdZ::
500 ComputeDimensions(G4Trd& trd, const G4int copyNo,
501  const G4VPhysicalVolume*) const
502 {
503  //---- The division along Z of a Trd will result a Trd
504  G4Trd* msol = (G4Trd*)(fmotherSolid);
505 
506  G4double pDx1 = msol->GetXHalfLength1();
507  G4double DDx = (msol->GetXHalfLength2() - msol->GetXHalfLength1() );
508  G4double pDy1 = msol->GetYHalfLength1();
509  G4double DDy = (msol->GetYHalfLength2() - msol->GetYHalfLength1() );
510  G4double pDz = fwidth/2. - fhgap;
511  G4double zLength = 2*msol->GetZHalfLength();
512 
513  trd.SetAllParameters( pDx1+DDx*(OffsetZ()+copyNo*fwidth+fhgap)/zLength,
514  pDx1+DDx*(OffsetZ()+(copyNo+1)*fwidth-fhgap)/zLength,
515  pDy1+DDy*(OffsetZ()+copyNo*fwidth+fhgap)/zLength,
516  pDy1+DDy*(OffsetZ()+(copyNo+1)*fwidth-fhgap)/zLength,
517  pDz );
518 
519 #ifdef G4DIVDEBUG
520  if( verbose >= 1 )
521  {
522  G4cout << " G4ParameterisationTrdZ::ComputeDimensions()"
523  << " - Mother TRD " << G4endl;
524  msol->DumpInfo();
525  G4cout << " - Parameterised TRD: "
526  << copyNo << G4endl;
527  trd.DumpInfo();
528  }
529 #endif
530 }