ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ParameterisationCons.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4ParameterisationCons.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 // G4ParameterisationCons[Rho/Phi/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 "G4ParameterisationCons.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 "G4Cons.hh"
42 
43 //--------------------------------------------------------------------------
44 G4VParameterisationCons::
45 G4VParameterisationCons( EAxis axis, G4int nDiv, G4double width,
46  G4double offset, G4VSolid* msolid,
47  DivisionType divType )
48  : G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
49 {
50  G4Cons* msol = (G4Cons*)(msolid);
51  if (msolid->GetEntityType() == "G4ReflectedSolid")
52  {
53  // Get constituent solid
54  G4VSolid* mConstituentSolid
55  = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
56  msol = (G4Cons*)(mConstituentSolid);
57 
58  // Create a new solid with inversed parameters
59  G4Cons* newSolid
60  = new G4Cons(msol->GetName(),
61  msol->GetInnerRadiusPlusZ(), msol->GetOuterRadiusPlusZ(),
62  msol->GetInnerRadiusMinusZ(), msol->GetOuterRadiusMinusZ(),
63  msol->GetZHalfLength(),
64  msol->GetStartPhiAngle(), msol->GetDeltaPhiAngle());
65  msol = newSolid;
66  fmotherSolid = newSolid;
67  fReflectedSolid = true;
68  fDeleteSolid = true;
69  }
70 }
71 
72 //------------------------------------------------------------------------
73 G4VParameterisationCons::~G4VParameterisationCons()
74 {
75 }
76 
77 //--------------------------------------------------------------------------
78 G4ParameterisationConsRho::
79 G4ParameterisationConsRho( EAxis axis, G4int nDiv,
80  G4double width, G4double offset,
81  G4VSolid* msolid, DivisionType divType )
82  : G4VParameterisationCons( axis, nDiv, width, offset, msolid, divType )
83 {
84  CheckParametersValidity();
85  SetType( "DivisionConsRho" );
86 
87  G4Cons* msol = (G4Cons*)(fmotherSolid);
88  if( msol->GetInnerRadiusPlusZ() == 0. )
89  {
90  std::ostringstream message;
91  message << "OuterRadiusMinusZ = 0" << G4endl
92  << "Width is calculated as that of OuterRadiusMinusZ !";
93  G4Exception("G4ParameterisationConsRho::G4ParameterisationConsRho()",
94  "GeomDiv1001", JustWarning, message);
95  }
96 
97  if( divType == DivWIDTH )
98  {
99  fnDiv = CalculateNDiv( msol->GetOuterRadiusMinusZ()
100  - msol->GetInnerRadiusMinusZ(), width, offset );
101  }
102  else if( divType == DivNDIV )
103  {
104  G4Cons* mconsol = (G4Cons*)(msolid);
105  fwidth = CalculateWidth( mconsol->GetOuterRadiusMinusZ()
106  - mconsol->GetInnerRadiusMinusZ(), nDiv, offset );
107  }
108 
109 #ifdef G4DIVDEBUG
110  if( verbose >= 1 )
111  {
112  G4cout << " G4ParameterisationConsRho - no divisions " << fnDiv << " = "
113  << nDiv << G4endl
114  << " Offset " << foffset << " = " << offset
115  << " - Width " << fwidth << " = " << width << G4endl;
116  }
117 #endif
118 }
119 
120 //--------------------------------------------------------------------------
121 G4ParameterisationConsRho::~G4ParameterisationConsRho()
122 {
123 }
124 
125 //------------------------------------------------------------------------
126 G4double G4ParameterisationConsRho::GetMaxParameter() const
127 {
128  G4Cons* msol = (G4Cons*)(fmotherSolid);
129  return msol->GetOuterRadiusMinusZ() - msol->GetInnerRadiusMinusZ();
130 }
131 
132 //--------------------------------------------------------------------------
133 void
134 G4ParameterisationConsRho::
135 ComputeTransformation( const G4int, G4VPhysicalVolume *physVol ) const
136 {
137  //----- translation
138  G4ThreeVector origin(0.,0.,0.);
139  //----- set translation
140  physVol->SetTranslation( origin );
141 
142  //----- calculate rotation matrix: unit
143 
144 #ifdef G4DIVDEBUG
145  if( verbose >= 2 )
146  {
147  G4cout << " G4ParameterisationConsRho " << G4endl
148  << " Offset: " << foffset
149  << " - Width: " << fwidth << G4endl;
150  }
151 #endif
152 
153  ChangeRotMatrix( physVol );
154 
155 #ifdef G4DIVDEBUG
156  if( verbose >= 2 )
157  {
158  G4cout << std::setprecision(8) << " G4ParameterisationConsRho" << G4endl
159  << " Position: " << origin << " - Width: " << fwidth
160  << " - Axis: " << faxis << G4endl;
161  }
162 #endif
163 }
164 
165 //--------------------------------------------------------------------------
166 void
167 G4ParameterisationConsRho::
168 ComputeDimensions( G4Cons& cons, const G4int copyNo,
169  const G4VPhysicalVolume* ) const
170 {
171  G4Cons* msol = (G4Cons*)(fmotherSolid);
172 
173  G4double pRMin1 = msol->GetInnerRadiusMinusZ() + foffset + fwidth*copyNo;
174  G4double pRMax1 = msol->GetInnerRadiusMinusZ() + foffset + fwidth*(copyNo+1);
175 
176  //width at Z Plus
177  //- G4double fwidthPlus =
178  // fwidth * ( msol->GetOuterRadiusPlusZ()/ msol->GetInnerRadiusPlusZ())
179  //- / ( msol->GetOuterRadiusMinusZ() - msol->GetInnerRadiusMinusZ());
180  G4double fwidthPlus = CalculateWidth( msol->GetOuterRadiusPlusZ()
181  - msol->GetInnerRadiusPlusZ(), fnDiv, foffset );
182  G4double pRMin2 = msol->GetInnerRadiusPlusZ()
183  + foffset + fwidthPlus * copyNo;
184  G4double pRMax2 = msol->GetInnerRadiusPlusZ()
185  + foffset + fwidthPlus * (copyNo+1);
186  G4double pDz = msol->GetZHalfLength();
187 
188  G4double d_half_gap = fhgap * pRMax2 / pRMax1;
189  //- already rotated double pSR = foffset + copyNo*fwidth;
190  G4double pSPhi = msol->GetStartPhiAngle();
191  G4double pDPhi = msol->GetDeltaPhiAngle();;
192 
193  cons.SetInnerRadiusMinusZ( pRMin1 + fhgap );
194  cons.SetOuterRadiusMinusZ( pRMax1 - fhgap );
195  cons.SetInnerRadiusPlusZ( pRMin2 + d_half_gap );
196  cons.SetOuterRadiusPlusZ( pRMax2 - d_half_gap );
197  cons.SetZHalfLength( pDz );
198  cons.SetStartPhiAngle( pSPhi, false );
199  cons.SetDeltaPhiAngle( pDPhi );
200 
201 #ifdef G4DIVDEBUG
202  if( verbose >= 2 )
203  {
204  G4cout << " G4ParameterisationConsRho::ComputeDimensions()" << G4endl
205  << " pRMin: " << pRMin1 << " - pRMax: " << pRMax1 << G4endl;
206  if( verbose >= 4 ) cons.DumpInfo();
207  }
208 #endif
209 }
210 
211 //--------------------------------------------------------------------------
212 G4ParameterisationConsPhi::
213 G4ParameterisationConsPhi( EAxis axis, G4int nDiv,
214  G4double width, G4double offset,
215  G4VSolid* msolid, DivisionType divType )
216  : G4VParameterisationCons( axis, nDiv, width, offset, msolid, divType )
217 {
218  CheckParametersValidity();
219  SetType( "DivisionConsPhi" );
220 
221  G4Cons* msol = (G4Cons*)(fmotherSolid);
222  if( divType == DivWIDTH )
223  {
224  fnDiv = CalculateNDiv( msol->GetDeltaPhiAngle(), width, offset );
225  }
226  else if( divType == DivNDIV )
227  {
228  fwidth = CalculateWidth( msol->GetDeltaPhiAngle(), nDiv, offset );
229  }
230 
231 #ifdef G4DIVDEBUG
232  if( verbose >= 1 )
233  {
234  G4cout << " G4ParameterisationConsPhi no divisions " << fnDiv << " = "
235  << nDiv << G4endl
236  << " Offset " << foffset << " = " << offset << G4endl
237  << " Width " << fwidth << " = " << width << G4endl;
238  }
239 #endif
240 }
241 
242 //--------------------------------------------------------------------------
243 G4ParameterisationConsPhi::~G4ParameterisationConsPhi()
244 {
245 }
246 
247 //------------------------------------------------------------------------
248 G4double G4ParameterisationConsPhi::GetMaxParameter() const
249 {
250  G4Cons* msol = (G4Cons*)(fmotherSolid);
251  return msol->GetDeltaPhiAngle();
252 }
253 
254 //--------------------------------------------------------------------------
255 void
256 G4ParameterisationConsPhi::
257 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
258 {
259  //----- translation
260  G4ThreeVector origin(0.,0.,0.);
261  //----- set translation
262  physVol->SetTranslation( origin );
263 
264  //----- calculate rotation matrix (so that all volumes point to the centre)
265  G4double posi = foffset + copyNo*fwidth;
266 
267 #ifdef G4DIVDEBUG
268  if( verbose >= 2 )
269  {
270  G4cout << " G4ParameterisationConsPhi - position: " << posi/deg << G4endl
271  << " Origin: " << origin << " copyNo: " << copyNo
272  << " - foffset: " << foffset/deg
273  << " - fwidth: " << fwidth/deg << G4endl
274  << " - Axis: " << faxis << G4endl;
275  }
276 #endif
277 
278  ChangeRotMatrix( physVol, -posi );
279 }
280 
281 //--------------------------------------------------------------------------
282 void
283 G4ParameterisationConsPhi::
284 ComputeDimensions( G4Cons& cons, const G4int,
285  const G4VPhysicalVolume* ) const
286 {
287  G4Cons* msol = (G4Cons*)(fmotherSolid);
288 
289  G4double pRMin1 = msol->GetInnerRadiusMinusZ();
290  G4double pRMax1 = msol->GetOuterRadiusMinusZ();
291  G4double pRMin2 = msol->GetInnerRadiusPlusZ();
292  G4double pRMax2 = msol->GetOuterRadiusPlusZ();
293  G4double pDz = msol->GetZHalfLength();
294 
295  //- already rotated double pSPhi = foffset + copyNo*fwidth;
296  G4double pSPhi = foffset + msol->GetStartPhiAngle() + fhgap;
297  G4double pDPhi = fwidth - 2.*fhgap;
298 
299  cons.SetInnerRadiusMinusZ( pRMin1 );
300  cons.SetOuterRadiusMinusZ( pRMax1 );
301  cons.SetInnerRadiusPlusZ( pRMin2 );
302  cons.SetOuterRadiusPlusZ( pRMax2 );
303  cons.SetZHalfLength( pDz );
304  cons.SetStartPhiAngle( pSPhi, false );
305  cons.SetDeltaPhiAngle( pDPhi );
306 
307 #ifdef G4DIVDEBUG
308  if( verbose >= 2 )
309  {
310  G4cout << " G4ParameterisationConsPhi::ComputeDimensions" << G4endl
311  << " pSPhi: " << pSPhi << " - pDPhi: " << pDPhi << G4endl;
312  if( verbose >= 4 ) cons.DumpInfo();
313  }
314 #endif
315 }
316 
317 //--------------------------------------------------------------------------
318 G4ParameterisationConsZ::
319 G4ParameterisationConsZ( EAxis axis, G4int nDiv,
320  G4double width, G4double offset,
321  G4VSolid* msolid, DivisionType divType )
322  : G4VParameterisationCons( axis, nDiv, width, offset, msolid, divType )
323 {
324  CheckParametersValidity();
325  SetType( "DivisionConsZ" );
326 
327  G4Cons* msol = (G4Cons*)(fmotherSolid);
328  if( divType == DivWIDTH )
329  {
330  fnDiv = CalculateNDiv( 2*msol->GetZHalfLength(), width, offset );
331  }
332  else if( divType == DivNDIV )
333  {
334  fwidth = CalculateWidth( 2*msol->GetZHalfLength(), nDiv, offset );
335  }
336 
337 #ifdef G4DIVDEBUG
338  if( verbose >= 1 )
339  {
340  G4cout << " G4ParameterisationConsZ: # divisions " << fnDiv << " = "
341  << nDiv << G4endl
342  << " Offset " << foffset << " = " << offset << G4endl
343  << " Width " << fwidth << " = " << width << G4endl
344  << " - Axis: " << faxis << G4endl;
345  }
346 #endif
347 }
348 
349 //--------------------------------------------------------------------------
350 G4ParameterisationConsZ::~G4ParameterisationConsZ()
351 {
352 }
353 
354 //------------------------------------------------------------------------
355 G4double G4ParameterisationConsZ::GetMaxParameter() const
356 {
357  G4Cons* msol = (G4Cons*)(fmotherSolid);
358  return 2*msol->GetZHalfLength();
359 }
360 
361 //--------------------------------------------------------------------------
362 void
363 G4ParameterisationConsZ::
364 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
365 {
366  //----- set translation: along Z axis
367  G4Cons* motherCons = (G4Cons*)(GetMotherSolid());
368  G4double posi = - motherCons->GetZHalfLength() + OffsetZ()
369  + fwidth/2 + copyNo*fwidth;
370  G4ThreeVector origin(0.,0.,posi);
371  physVol->SetTranslation( origin );
372 
373  //----- calculate rotation matrix: unit
374 
375 #ifdef G4DIVDEBUG
376  if( verbose >= 2 )
377  {
378  G4cout << " G4ParameterisationConsZ::ComputeTransformation()" << G4endl
379  << " Origin: " << origin << " - copyNo: " << copyNo << G4endl
380  << " foffset: " << foffset << " - fwidth: " << fwidth
381  << G4endl;
382  }
383 #endif
384 
385  ChangeRotMatrix( physVol );
386 }
387 
388 
389 //--------------------------------------------------------------------------
390 void
391 G4ParameterisationConsZ::
392 ComputeDimensions( G4Cons& cons, const G4int copyNo,
393  const G4VPhysicalVolume* ) const
394 {
395  G4Cons* msol = (G4Cons*)(fmotherSolid);
396 
397  G4double mHalfLength = msol->GetZHalfLength() - fhgap;
398  G4double aRInner = (msol->GetInnerRadiusPlusZ()
399  - msol->GetInnerRadiusMinusZ()) / (2*mHalfLength);
400  G4double bRInner = (msol->GetInnerRadiusPlusZ()
401  + msol->GetInnerRadiusMinusZ()) / 2;
402  G4double aROuter = (msol->GetOuterRadiusPlusZ()
403  - msol->GetOuterRadiusMinusZ()) / (2*mHalfLength);
404  G4double bROuter = (msol->GetOuterRadiusPlusZ()
405  + msol->GetOuterRadiusMinusZ()) / 2;
406  G4double xMinusZ = -mHalfLength + OffsetZ() + fwidth*copyNo + fhgap;
407  G4double xPlusZ = -mHalfLength + OffsetZ() + fwidth*(copyNo+1) - fhgap;
408  cons.SetInnerRadiusMinusZ( aRInner * xMinusZ + bRInner );
409  cons.SetOuterRadiusMinusZ( aROuter * xMinusZ + bROuter );
410  cons.SetInnerRadiusPlusZ( aRInner * xPlusZ + bRInner );
411  cons.SetOuterRadiusPlusZ( aROuter * xPlusZ + bROuter );
412 
413  G4double pDz = fwidth / 2. - fhgap;
414  G4double pSPhi = msol->GetStartPhiAngle();
415  G4double pDPhi = msol->GetDeltaPhiAngle();
416 
417  cons.SetZHalfLength( pDz );
418  cons.SetStartPhiAngle( pSPhi, false );
419  cons.SetDeltaPhiAngle( pDPhi );
420 
421 #ifdef G4DIVDEBUG
422  if( verbose >= 2 )
423  {
424  G4cout << " G4ParameterisationConsZ::ComputeDimensions()" << G4endl
425  << " pDz: " << pDz << G4endl;
426  if( verbose >= 4 ) cons.DumpInfo();
427  }
428 #endif
429 
430 }