source: trunk/source/geometry/divisions/src/G4ParameterisationTrd.cc@ 1272

Last change on this file since 1272 was 1228, checked in by garnier, 16 years ago

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