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

Last change on this file since 1351 was 1347, checked in by garnier, 15 years ago

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