Geant4  10.02.p01
G4ParameterisationBox.cc
Go to the documentation of this file.
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: G4ParameterisationBox.cc 68040 2013-03-13 14:19:04Z gcosmo $
28 //
29 // class G4ParameterisationBox Implementation file
30 //
31 // 26.05.03 - P.Arce, Initial version
32 // 08.04.04 - I.Hrivnacova, Implemented reflection
33 // 21.04.10 - M.Asai, Added gaps
34 // --------------------------------------------------------------------
35 
36 #include "G4ParameterisationBox.hh"
37 
38 #include <iomanip>
39 #include "G4ThreeVector.hh"
40 #include "G4Transform3D.hh"
41 #include "G4RotationMatrix.hh"
42 #include "G4VPhysicalVolume.hh"
43 #include "G4ReflectedSolid.hh"
44 #include "G4Box.hh"
45 
46 //--------------------------------------------------------------------------
49  G4double offset, G4VSolid* msolid,
50  DivisionType divType )
51  : G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
52 {
53  G4Box* msol = (G4Box*)(msolid);
54  if (msolid->GetEntityType() == "G4ReflectedSolid")
55  {
56  // Get constituent solid
57  G4VSolid* mConstituentSolid
58  = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
59  msol = (G4Box*)(mConstituentSolid);
60  fmotherSolid = msol;
61  fReflectedSolid = true;
62  }
63 }
64 
65 //--------------------------------------------------------------------------
67 {
68 }
69 
70 //--------------------------------------------------------------------------
73  G4double offset, G4VSolid* msolid,
74  DivisionType divType )
75  : G4VParameterisationBox( axis, nDiv, width, offset, msolid, divType )
76 {
78  SetType( "DivisionBoxX" );
79 
80  G4Box* mbox = (G4Box*)(fmotherSolid);
81  if( divType == DivWIDTH )
82  {
83  fnDiv = CalculateNDiv( 2*mbox->GetXHalfLength(), width, offset );
84  }
85  else if( divType == DivNDIV )
86  {
87  fwidth = CalculateWidth( 2*mbox->GetXHalfLength(), nDiv, offset );
88  }
89 #ifdef G4DIVDEBUG
90  if( verbose >= 1 )
91  {
92  G4cout << " G4ParameterisationBoxX - no divisions "
93  << fnDiv << " = " << nDiv << G4endl
94  << " Offset " << foffset << " = " << offset << G4endl
95  << " Width " << fwidth << " = " << width << G4endl;
96  }
97 #endif
98 }
99 
100 //------------------------------------------------------------------------
102 {
103 }
104 
105 //------------------------------------------------------------------------
107 {
108  G4Box* msol = (G4Box*)(fmotherSolid);
109  return 2*msol->GetXHalfLength();
110 }
111 
112 //------------------------------------------------------------------------
113 void
115 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
116 {
117  G4Box* msol = (G4Box*)(fmotherSolid );
118  G4double mdx = msol->GetXHalfLength( );
119 
120  //----- translation
121  G4ThreeVector origin(0.,0.,0.);
122  G4double posi = -mdx + foffset+(copyNo+0.5)*fwidth;
123 
124  if( faxis == kXAxis )
125  {
126  origin.setX( posi );
127  }
128  else
129  {
130  std::ostringstream message;
131  message << "Only axes along X are allowed ! Axis: " << faxis;
132  G4Exception("G4ParameterisationBoxX::ComputeTransformation()",
133  "GeomDiv0002", FatalException, message);
134  }
135 #ifdef G4DIVDEBUG
136  if( verbose >= 2 )
137  {
138  G4cout << std::setprecision(8) << " G4ParameterisationBoxX: "
139  << copyNo << G4endl
140  << " Position " << origin << " Axis " << faxis << G4endl;
141  }
142 #endif
143  //----- set translation
144  physVol->SetTranslation( origin );
145 }
146 
147 //------------------------------------------------------------------------
148 void
151  const G4VPhysicalVolume* ) const
152 {
153  G4Box* msol = (G4Box*)(fmotherSolid);
154 
155  G4double pDx = fwidth/2. - fhgap;
156  G4double pDy = msol->GetYHalfLength();
157  G4double pDz = msol->GetZHalfLength();
158 
159  box.SetXHalfLength( pDx );
160  box.SetYHalfLength( pDy );
161  box.SetZHalfLength( pDz );
162 
163 #ifdef G4DIVDEBUG
164  if( verbose >= 2 )
165  {
166  G4cout << " G4ParameterisationBoxX::ComputeDimensions()" << G4endl
167  << " pDx: " << pDz << G4endl;
168  box.DumpInfo();
169  }
170 #endif
171 }
172 
173 //------------------------------------------------------------------------
176  G4double offset, G4VSolid* msolid,
177  DivisionType divType)
178  : G4VParameterisationBox( axis, nDiv, width, offset, msolid, divType )
179 {
181  SetType( "DivisionBoxY" );
182 
183  G4Box* mbox = (G4Box*)(fmotherSolid);
184  if( divType == DivWIDTH )
185  {
186  fnDiv = CalculateNDiv( 2*mbox->GetYHalfLength(), width, offset );
187  }
188  else if( divType == DivNDIV )
189  {
190  fwidth = CalculateWidth( 2*mbox->GetYHalfLength(), nDiv, offset );
191  }
192 
193 #ifdef G4DIVDEBUG
194  if( verbose >= 1 )
195  {
196  G4cout << " G4ParameterisationBoxY - no divisions " << fnDiv << " = "
197  << nDiv << ". Offset " << foffset << " = " << offset
198  << ". Width " << fwidth << " = " << width << G4endl;
199  }
200 #endif
201 }
202 
203 //------------------------------------------------------------------------
205 {
206 }
207 
208 //------------------------------------------------------------------------
210 {
211  G4Box* msol = (G4Box*)(fmotherSolid);
212  return 2*msol->GetYHalfLength();
213 }
214 
215 //------------------------------------------------------------------------
216 void
218 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
219 {
220  G4Box* msol = (G4Box*)(fmotherSolid);
221  G4double mdy = msol->GetYHalfLength();
222 
223  //----- translation
224  G4ThreeVector origin(0.,0.,0.);
225  G4double posi = -mdy + foffset + (copyNo+0.5)*fwidth;
226  if( faxis == kYAxis )
227  {
228  origin.setY( posi );
229  }
230  else
231  {
232  std::ostringstream message;
233  message << "Only axes along Y are allowed ! Axis: " << faxis;
234  G4Exception("G4ParameterisationBoxY::ComputeTransformation()",
235  "GeomDiv0002", FatalException, message);
236  }
237 #ifdef G4DIVDEBUG
238  if( verbose >= 2 )
239  {
240  G4cout << std::setprecision(8) << " G4ParameterisationBoxY: "
241  << copyNo << G4endl
242  << " Position " << origin << " Axis " << faxis << G4endl;
243  }
244 #endif
245  //----- set translation
246  physVol->SetTranslation( origin );
247 }
248 
249 //------------------------------------------------------------------------
250 void
253  const G4VPhysicalVolume* ) const
254 {
255  G4Box* msol = (G4Box*)(fmotherSolid);
256 
257  G4double pDx = msol->GetXHalfLength();
258  G4double pDy = fwidth/2. - fhgap;
259  G4double pDz = msol->GetZHalfLength();
260 
261  box.SetXHalfLength( pDx );
262  box.SetYHalfLength( pDy );
263  box.SetZHalfLength( pDz );
264 
265 #ifdef G4DIVDEBUG
266  if( verbose >= 2 )
267  {
268  G4cout << " G4ParameterisationBoxY::ComputeDimensions()" << G4endl
269  << " pDx: " << pDz << G4endl;
270  box.DumpInfo();
271  }
272 #endif
273 }
274 
275 //------------------------------------------------------------------------
278  G4double offset, G4VSolid* msolid,
279  DivisionType divType )
280  : G4VParameterisationBox( axis, nDiv, width, offset, msolid, divType )
281 {
283  SetType( "DivisionBoxZ" );
284 
285  G4Box* mbox = (G4Box*)(fmotherSolid);
286  if( divType == DivWIDTH )
287  {
288  fnDiv = CalculateNDiv( 2*mbox->GetZHalfLength(), width, offset );
289  }
290  else if ( divType == DivNDIV )
291  {
292  fwidth = CalculateWidth( 2*mbox->GetZHalfLength(), nDiv, offset );
293  }
294 #ifdef G4DIVDEBUG
295  if( verbose >= 1 )
296  {
297  G4cout << " G4ParameterisationBoxZ - no divisions " << fnDiv << " = "
298  << nDiv << ". Offset " << foffset << " = " << offset
299  << ". Width " << fwidth << " = " << width << G4endl;
300  }
301 #endif
302 }
303 
304 //------------------------------------------------------------------------
306 {
307 }
308 
309 //------------------------------------------------------------------------
311 {
312  G4Box* msol = (G4Box*)(fmotherSolid);
313  return 2*msol->GetZHalfLength();
314 }
315 
316 //------------------------------------------------------------------------
317 void
319 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
320 {
321  G4Box* msol = (G4Box*)(fmotherSolid );
322  G4double mdz = msol->GetZHalfLength();
323 
324  //----- translation
325  G4ThreeVector origin(0.,0.,0.);
326  G4double posi = -mdz + OffsetZ() + (copyNo+0.5)*fwidth;
327 
328  if( faxis == kZAxis )
329  {
330  origin.setZ( posi );
331  }
332  else
333  {
334  std::ostringstream message;
335  message << "Only axes along Z are allowed ! Axis: " << faxis;
336  G4Exception("G4ParameterisationBoxZ::ComputeTransformation()",
337  "GeomDiv0002", FatalException, message);
338  }
339 #ifdef G4DIVDEBUG
340  if( verbose >= 2 )
341  {
342  G4cout << std::setprecision(8) << " G4ParameterisationBoxZ: "
343  << copyNo << G4endl
344  << " Position " << origin << " Axis " << faxis << G4endl;
345  }
346 #endif
347  //----- set translation
348  physVol->SetTranslation( origin );
349 }
350 
351 //------------------------------------------------------------------------
352 void
355  const G4VPhysicalVolume* ) const
356 {
357  G4Box* msol = (G4Box*)(fmotherSolid);
358 
359  G4double pDx = msol->GetXHalfLength();
360  G4double pDy = msol->GetYHalfLength();
361  G4double pDz = fwidth/2. - fhgap;
362 
363  box.SetXHalfLength( pDx );
364  box.SetYHalfLength( pDy );
365  box.SetZHalfLength( pDz );
366 
367 #ifdef G4DIVDEBUG
368  if( verbose >= 2 )
369  {
370  G4cout << " G4ParameterisationBoxZ::ComputeDimensions()" << G4endl
371  << " pDx: " << pDz << G4endl;
372  box.DumpInfo();
373  }
374 #endif
375 }
376 
G4double GetXHalfLength() const
void SetZHalfLength(G4double dz)
Definition: G4Box.cc:171
G4ParameterisationBoxZ(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
G4ParameterisationBoxY(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
G4ParameterisationBoxX(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
CLHEP::Hep3Vector G4ThreeVector
void SetType(const G4String &type)
Definition: G4Box.hh:64
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
G4double GetMaxParameter() const
#define width
G4double CalculateWidth(G4double motherDim, G4int nDiv, G4double offset) const
virtual G4GeometryType GetEntityType() const =0
int G4int
Definition: G4Types.hh:78
G4double GetZHalfLength() const
void DumpInfo() const
G4GLOB_DLL std::ostream G4cout
G4double GetMaxParameter() const
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
void SetTranslation(const G4ThreeVector &v)
G4double GetYHalfLength() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
void ComputeDimensions(G4Box &box, const G4int copyNo, const G4VPhysicalVolume *physVol) const
G4VParameterisationBox(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
void ComputeDimensions(G4Box &box, const G4int copyNo, const G4VPhysicalVolume *physVol) const
EAxis
Definition: geomdefs.hh:54
void SetYHalfLength(G4double dy)
Definition: G4Box.cc:151
G4int CalculateNDiv(G4double motherDim, G4double width, G4double offset) const
void SetXHalfLength(G4double dx)
Definition: G4Box.cc:131
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4double GetMaxParameter() const
void ComputeDimensions(G4Box &box, const G4int copyNo, const G4VPhysicalVolume *physVol) const