Geant4  10.00.p01
G4ParameterisationPolycone.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: G4ParameterisationPolycone.cc 75219 2013-10-29 16:12:13Z gcosmo $
28 //
29 // class G4ParameterisationPolycone Implementation file
30 //
31 // 26.05.03 - P.Arce, Initial version
32 // 08.04.04 - I.Hrivnacova, Implemented reflection
33 //---------------------------------------------------------------------
34 
36 
37 #include <iomanip>
38 #include "G4ThreeVector.hh"
39 #include "G4RotationMatrix.hh"
40 #include "G4VPhysicalVolume.hh"
41 #include "G4LogicalVolume.hh"
42 #include "G4ReflectedSolid.hh"
43 
44 //-----------------------------------------------------------------------
47  G4double offset, G4VSolid* msolid,
48  DivisionType divType )
49  : G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
50 {
51  G4Polycone* msol = (G4Polycone*)(msolid);
52  if (msolid->GetEntityType() == "G4ReflectedSolid")
53  {
54  // Get constituent solid
55  G4VSolid* mConstituentSolid
56  = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
57  msol = (G4Polycone*)(mConstituentSolid);
58 
59  // Get parameters
60  G4int nofZplanes = msol->GetOriginalParameters()->Num_z_planes;
61  G4double* zValues = msol->GetOriginalParameters()->Z_values;
62  G4double* rminValues = msol->GetOriginalParameters()->Rmin;
63  G4double* rmaxValues = msol->GetOriginalParameters()->Rmax;
64 
65  // Invert z values
66  G4double* zValuesRefl = new double[nofZplanes];
67  for (G4int i=0; i<nofZplanes; i++) zValuesRefl[i] = - zValues[i];
68 
69  G4Polycone* newSolid
70  = new G4Polycone(msol->GetName(),
71  msol->GetStartPhi(),
72  msol->GetEndPhi() - msol->GetStartPhi(),
73  nofZplanes, zValuesRefl, rminValues, rmaxValues);
74 
75  delete [] zValuesRefl;
76 
77  msol = newSolid;
78  fmotherSolid = newSolid;
79  fReflectedSolid = true;
80  fDeleteSolid = true;
81  }
82 }
83 
84 //---------------------------------------------------------------------
86 {
87 }
88 
89 //---------------------------------------------------------------------
92  G4double width, G4double offset,
93  G4VSolid* msolid, DivisionType divType )
94  : G4VParameterisationPolycone( axis, nDiv, width, offset, msolid, divType )
95 {
97  SetType( "DivisionPolyconeRho" );
98 
100  G4PolyconeHistorical* origparamMother = msol->GetOriginalParameters();
101 
102  if( divType == DivWIDTH )
103  {
104  fnDiv = CalculateNDiv( origparamMother->Rmax[0]
105  - origparamMother->Rmin[0], width, offset );
106  }
107  else if( divType == DivNDIV )
108  {
109  fwidth = CalculateWidth( origparamMother->Rmax[0]
110  - origparamMother->Rmin[0], nDiv, offset );
111  }
112 
113 #ifdef G4DIVDEBUG
114  if( verbose >= -1 )
115  {
116  G4cout << " G4ParameterisationPolyconeRho - # divisions " << fnDiv
117  << " = " << nDiv << G4endl
118  << " Offset " << foffset << " = " << offset << G4endl
119  << " Width " << fwidth << " = " << width << G4endl;
120  }
121 #endif
122 }
123 
124 //---------------------------------------------------------------------
126 {
127 }
128 
129 //---------------------------------------------------------------------
131 {
133 
134  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
135 
137  {
138  std::ostringstream message;
139  message << "In solid " << msol->GetName() << G4endl
140  << "Division along R will be done with a width "
141  << "different for each solid section." << G4endl
142  << "WIDTH will not be used !";
143  G4Exception("G4VParameterisationPolycone::CheckParametersValidity()",
144  "GeomDiv1001", JustWarning, message);
145  }
146  if( foffset != 0. )
147  {
148  std::ostringstream message;
149  message << "In solid " << msol->GetName() << G4endl
150  << "Division along R will be done with a width "
151  << "different for each solid section." << G4endl
152  << "OFFSET will not be used !";
153  G4Exception("G4VParameterisationPolycone::CheckParametersValidity()",
154  "GeomDiv1001", JustWarning, message);
155  }
156 }
157 
158 //------------------------------------------------------------------------
160 {
161  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
162  G4PolyconeHistorical* original_pars = msol->GetOriginalParameters();
163  return original_pars->Rmax[0] - original_pars->Rmin[0];
164 }
165 
166 
167 //---------------------------------------------------------------------
168 void
171 {
172  //----- translation
173  G4ThreeVector origin(0.,0.,0.);
174  //----- set translation
175  physVol->SetTranslation( origin );
176 
177  //----- calculate rotation matrix: unit
178 
179 #ifdef G4DIVDEBUG
180  if( verbose >= 2 )
181  {
182  G4cout << " G4ParameterisationPolyconeRho " << G4endl
183  << " foffset: " << foffset
184  << " - fwidth: " << fwidth << G4endl;
185  }
186 #endif
187 
188  ChangeRotMatrix( physVol );
189 
190 #ifdef G4DIVDEBUG
191  if( verbose >= 2 )
192  {
193  G4cout << std::setprecision(8) << " G4ParameterisationPolyconeRho "
194  << G4endl
195  << " Position: " << origin/mm
196  << " - Width: " << fwidth/deg
197  << " - Axis: " << faxis << G4endl;
198  }
199 #endif
200 }
201 
202 //---------------------------------------------------------------------
203 void
205 ComputeDimensions( G4Polycone& pcone, const G4int copyNo,
206  const G4VPhysicalVolume* ) const
207 {
208  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
209 
210  G4PolyconeHistorical* origparamMother = msol->GetOriginalParameters();
211  G4PolyconeHistorical origparam( *origparamMother );
212  G4int nZplanes = origparamMother->Num_z_planes;
213 
214  G4double width = 0.;
215  for( G4int ii = 0; ii < nZplanes; ii++ )
216  {
217  width = CalculateWidth( origparamMother->Rmax[ii]
218  - origparamMother->Rmin[ii], fnDiv, foffset );
219  origparam.Rmin[ii] = origparamMother->Rmin[ii]+foffset+width*copyNo;
220  origparam.Rmax[ii] = origparamMother->Rmin[ii]+foffset+width*(copyNo+1);
221  }
222 
223  pcone.SetOriginalParameters(&origparam); // copy values & transfer pointers
224  pcone.Reset(); // reset to new solid parameters
225 
226 #ifdef G4DIVDEBUG
227  if( verbose >= -2 )
228  {
229  G4cout << "G4ParameterisationPolyconeRho::ComputeDimensions()" << G4endl
230  << "-- Parametrised pcone copy-number: " << copyNo << G4endl;
231  pcone.DumpInfo();
232  }
233 #endif
234 }
235 
236 //---------------------------------------------------------------------
239  G4double width, G4double offset,
240  G4VSolid* msolid, DivisionType divType )
241  : G4VParameterisationPolycone( axis, nDiv, width, offset, msolid, divType )
242 {
244  SetType( "DivisionPolyconePhi" );
245 
246  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
247  G4double deltaPhi = msol->GetEndPhi() - msol->GetStartPhi();
248 
249  if( divType == DivWIDTH )
250  {
251  fnDiv = CalculateNDiv( deltaPhi, width, offset );
252  }
253  else if( divType == DivNDIV )
254  {
255  fwidth = CalculateWidth( deltaPhi, nDiv, offset );
256  }
257 
258 #ifdef G4DIVDEBUG
259  if( verbose >= 1 )
260  {
261  G4cout << " G4ParameterisationPolyconePhi - # divisions " << fnDiv
262  << " = " << nDiv << G4endl
263  << " Offset " << foffset/deg << " = " << offset/deg << G4endl
264  << " Width " << fwidth/deg << " = " << width/deg << G4endl;
265  }
266 #endif
267 }
268 
269 //---------------------------------------------------------------------
271 {
272 }
273 
274 //------------------------------------------------------------------------
276 {
277  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
278  return msol->GetEndPhi() - msol->GetStartPhi();
279 }
280 
281 //---------------------------------------------------------------------
282 void
284 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
285 {
286  //----- translation
287  G4ThreeVector origin(0.,0.,0.);
288  //----- set translation
289  physVol->SetTranslation( origin );
290 
291  //----- calculate rotation matrix (so that all volumes point to the centre)
292  G4double posi = foffset + copyNo*fwidth;
293 
294 #ifdef G4DIVDEBUG
295  if( verbose >= 2 )
296  {
297  G4cout << " G4ParameterisationPolyconePhi - position: " << posi/deg
298  << G4endl
299  << " copyNo: " << copyNo << " - foffset: " << foffset/deg
300  << " - fwidth: " << fwidth/deg << G4endl;
301  }
302 #endif
303 
304  ChangeRotMatrix( physVol, -posi );
305 
306 #ifdef G4DIVDEBUG
307  if( verbose >= 2 )
308  {
309  G4cout << std::setprecision(8) << " G4ParameterisationPolyconePhi "
310  << copyNo << G4endl
311  << " Position: " << origin << " - Width: " << fwidth
312  << " - Axis: " << faxis << G4endl;
313  }
314 #endif
315 }
316 
317 //---------------------------------------------------------------------
318 void
321  const G4VPhysicalVolume* ) const
322 {
323  G4Polycone* msol = (G4Polycone*)(fmotherSolid);
324 
325  G4PolyconeHistorical* origparamMother = msol->GetOriginalParameters();
326  G4PolyconeHistorical origparam( *origparamMother );
327  origparam.Start_angle = origparamMother->Start_angle;
328  origparam.Opening_angle = fwidth;
329 
330  pcone.SetOriginalParameters(&origparam); // copy values & transfer pointers
331  pcone.Reset(); // reset to new solid parameters
332 
333 #ifdef G4DIVDEBUG
334  if( verbose >= 2 )
335  {
336  G4cout << "G4ParameterisationPolyconePhi::ComputeDimensions():" << G4endl;
337  pcone.DumpInfo();
338  }
339 #endif
340 }
341 
342 //---------------------------------------------------------------------
345  G4double width, G4double offset,
346  G4VSolid* msolid, DivisionType divType)
347  : G4VParameterisationPolycone( axis, nDiv, width, offset, msolid, divType ),
348  fNSegment(0),
349  fOrigParamMother(((G4Polycone*)fmotherSolid)->GetOriginalParameters())
350 {
351 
353  SetType( "DivisionPolyconeZ" );
354 
355  if( divType == DivWIDTH )
356  {
357  fnDiv =
359  - fOrigParamMother->Z_values[0] , width, offset );
360  }
361  else if( divType == DivNDIV )
362  {
363  fwidth =
365  - fOrigParamMother->Z_values[0] , nDiv, offset );
366  }
367 
368 #ifdef G4DIVDEBUG
369  if( verbose >= 1 )
370  {
371  G4cout << " G4ParameterisationPolyconeZ - # divisions " << fnDiv << " = "
372  << nDiv << G4endl
373  << " Offset " << foffset << " = " << offset << G4endl
374  << " Width " << fwidth << " = " << width << G4endl;
375  }
376 #endif
377 }
378 
379 //---------------------------------------------------------------------
381 {
382 }
383 
384 //------------------------------------------------------------------------
386  G4double z1, G4double r1,
387  G4double z2, G4double r2) const
388 {
389  // Linear parameterisation:
390  // r = az + b
391  // a = (r1 - r2)/(z1-z2)
392  // b = r1 - a*z1
393 
394  return (r1-r2)/(z1-z2)*z + ( r1 - (r1-r2)/(z1-z2)*z1 ) ;
395 }
396 
397 //------------------------------------------------------------------------
399 {
400 // Get Rmin in the given z position for the given polycone segment
401 
402  return GetR(z,
403  fOrigParamMother->Z_values[nseg],
404  fOrigParamMother->Rmin[nseg],
405  fOrigParamMother->Z_values[nseg+1],
406  fOrigParamMother->Rmin[nseg+1]);
407 }
408 
409 //------------------------------------------------------------------------
411 {
412 // Get Rmax in the given z position for the given polycone segment
413 
414  return GetR(z,
415  fOrigParamMother->Z_values[nseg],
416  fOrigParamMother->Rmax[nseg],
417  fOrigParamMother->Z_values[nseg+1],
418  fOrigParamMother->Rmax[nseg+1]);
419 }
420 
421 //------------------------------------------------------------------------
423 {
426 }
427 
428 //---------------------------------------------------------------------
430 {
432 
433  // Division will be following the mother polycone segments
434  if( fDivisionType == DivNDIV ) {
435  if( fnDiv > fOrigParamMother->Num_z_planes-1 ) {
436  std::ostringstream error;
437  error << "Configuration not supported." << G4endl
438  << "Division along Z will be done by splitting in the defined"
439  << G4endl
440  << "Z planes, i.e, the number of division would be: "
442  << ", instead of: " << fnDiv << " !";
443  G4Exception("G4ParameterisationPolyconeZ::CheckParametersValidity()",
444  "GeomDiv0001", FatalException, error);
445  }
446  }
447 
448  // Division will be done within one polycone segment
449  // with applying given width and offset
451  // Check if divided region does not span over more
452  // than one z segment
453 
454  G4int isegstart = -1; // number of the segment containing start position
455  G4int isegend = -1; // number of the segment containing end position
456 
457  if ( ! fReflectedSolid ) {
458  // The start/end position of the divided region
459  G4double zstart
461  G4double zend
463 
464  G4int counter = 0;
465  while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) {
466  // first segment
467  if ( zstart >= fOrigParamMother->Z_values[counter] &&
468  zstart < fOrigParamMother->Z_values[counter+1] ) {
469  isegstart = counter;
470  }
471  // last segment
472  if ( zend > fOrigParamMother->Z_values[counter] &&
473  zend <= fOrigParamMother->Z_values[counter+1] ) {
474  isegend = counter;
475  }
476  ++counter;
477  }
478  }
479  else {
480  // The start/end position of the divided region
481  G4double zstart
483  G4double zend
485 
486  G4int counter = 0;
487  while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) {
488  // first segment
489  if ( zstart <= fOrigParamMother->Z_values[counter] &&
490  zstart > fOrigParamMother->Z_values[counter+1] ) {
491  isegstart = counter;
492  }
493  // last segment
494  if ( zend < fOrigParamMother->Z_values[counter] &&
495  zend >= fOrigParamMother->Z_values[counter+1] ) {
496  isegend = counter;
497  }
498  ++counter;
499  }
500  }
501 
502 
503  if ( isegstart != isegend ) {
504  std::ostringstream message;
505  message << "Condiguration not supported." << G4endl
506  << "Division with user defined width." << G4endl
507  << "Solid " << fmotherSolid->GetName() << G4endl
508  << "Divided region is not between two z planes.";
509  G4Exception("G4ParameterisationPolyconeZ::CheckParametersValidity()",
510  "GeomDiv0001", FatalException, message);
511  }
512 
513  fNSegment = isegstart;
514  }
515 }
516 
517 //---------------------------------------------------------------------
518 void
520 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol) const
521 {
522  if ( fDivisionType == DivNDIV ) {
523  // The position of the centre of copyNo-th mother polycone segment
524  G4double posi
525  = ( fOrigParamMother->Z_values[copyNo]
526  + fOrigParamMother->Z_values[copyNo+1])/2;
527  physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
528  }
529 
531  // The position of the centre of copyNo-th division
533 
534  if ( ! fReflectedSolid )
535  posi += foffset + (2*copyNo + 1) * fwidth/2.;
536  else
537  posi -= foffset + (2*copyNo + 1) * fwidth/2.;
538 
539  physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
540  }
541 
542  //----- calculate rotation matrix: unit
543 
544 #ifdef G4DIVDEBUG
545  if( verbose >= 2 )
546  {
547  G4cout << " G4ParameterisationPolyconeZ - position: " << posi << G4endl
548  << " copyNo: " << copyNo << " - foffset: " << foffset/deg
549  << " - fwidth: " << fwidth/deg << G4endl;
550  }
551 #endif
552 
553  ChangeRotMatrix( physVol );
554 
555 #ifdef G4DIVDEBUG
556  if( verbose >= 2 )
557  {
558  G4cout << std::setprecision(8) << " G4ParameterisationPolyconeZ "
559  << copyNo << G4endl
560  << " Position: " << origin << " - Width: " << fwidth
561  << " - Axis: " << faxis << G4endl;
562  }
563 #endif
564 }
565 
566 //---------------------------------------------------------------------
567 void
569 ComputeDimensions( G4Polycone& pcone, const G4int copyNo,
570  const G4VPhysicalVolume* ) const
571 {
572 
573  // Define division solid
574  G4PolyconeHistorical origparam;
575  G4int nz = 2;
576  origparam.Num_z_planes = nz;
579 
580  // Define division solid z sections
581  origparam.Z_values = new G4double[nz];
582  origparam.Rmin = new G4double[nz];
583  origparam.Rmax = new G4double[nz];
584 
585  if ( fDivisionType == DivNDIV ) {
586  // The position of the centre of copyNo-th mother polycone segment
587  G4double posi = (fOrigParamMother->Z_values[copyNo]
588  + fOrigParamMother->Z_values[copyNo+1])/2;
589 
590  origparam.Z_values[0] = fOrigParamMother->Z_values[copyNo] - posi;
591  origparam.Z_values[1] = fOrigParamMother->Z_values[copyNo+1] - posi;
592  origparam.Rmin[0] = fOrigParamMother->Rmin[copyNo];
593  origparam.Rmin[1] = fOrigParamMother->Rmin[copyNo+1];
594  origparam.Rmax[0] = fOrigParamMother->Rmax[copyNo];
595  origparam.Rmax[1] = fOrigParamMother->Rmax[copyNo+1];
596  }
597 
599  if ( ! fReflectedSolid ) {
600  origparam.Z_values[0] = - fwidth/2.;
601  origparam.Z_values[1] = fwidth/2.;
602 
603  // The position of the centre of copyNo-th division
604  G4double posi
605  = fOrigParamMother->Z_values[0] + foffset + (2*copyNo + 1) * fwidth/2.;
606 
607  // The first and last z sides z values
608  G4double zstart = posi - fwidth/2.;
609  G4double zend = posi + fwidth/2.;
610  origparam.Rmin[0] = GetRmin(zstart, fNSegment);
611  origparam.Rmax[0] = GetRmax(zstart, fNSegment);
612  origparam.Rmin[1] = GetRmin(zend, fNSegment);
613  origparam.Rmax[1] = GetRmax(zend, fNSegment);
614  }
615  else {
616  origparam.Z_values[0] = fwidth/2.;
617  origparam.Z_values[1] = - fwidth/2.;
618 
619  // The position of the centre of copyNo-th division
620  G4double posi
621  = fOrigParamMother->Z_values[0] - ( foffset + (2*copyNo + 1) * fwidth/2.);
622 
623  // The first and last z sides z values
624  G4double zstart = posi + fwidth/2.;
625  G4double zend = posi - fwidth/2.;
626  origparam.Rmin[0] = GetRmin(zstart, fNSegment);
627  origparam.Rmax[0] = GetRmax(zstart, fNSegment);
628  origparam.Rmin[1] = GetRmin(zend, fNSegment);
629  origparam.Rmax[1] = GetRmax(zend, fNSegment);
630  }
631 
632  // It can happen due to rounding errors
633  if ( origparam.Rmin[0] < 0.0 ) origparam.Rmin[0] = 0.0;
634  if ( origparam.Rmin[nz-1] < 0.0 ) origparam.Rmin[1] = 0.0;
635  }
636 
637  pcone.SetOriginalParameters(&origparam); // copy values & transfer pointers
638  pcone.Reset(); // reset to new solid parameters
639 
640 #ifdef G4DIVDEBUG
641  if( verbose >= 2 )
642  {
643  G4cout << "G4ParameterisationPolyconeZ::ComputeDimensions()" << G4endl
644  << "-- Parametrised pcone copy-number: " << copyNo << G4endl;
645  pcone.DumpInfo();
646  }
647 #endif
648 }
G4String GetName() const
G4double GetRmax(G4double z, G4int nsegment) const
CLHEP::Hep3Vector G4ThreeVector
void SetType(const G4String &type)
G4double z
Definition: TRTMaterials.hh:39
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
#define width
G4double CalculateWidth(G4double motherDim, G4int nDiv, G4double offset) const
virtual G4GeometryType GetEntityType() const =0
G4double GetEndPhi() const
int G4int
Definition: G4Types.hh:78
G4ParameterisationPolyconeRho(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *motherSolid, DivisionType divType)
void DumpInfo() const
G4GLOB_DLL std::ostream G4cout
void ChangeRotMatrix(G4VPhysicalVolume *physVol, G4double rotZ=0.) const
static const double deg
Definition: G4SIunits.hh:133
void ComputeDimensions(G4Polycone &pcone, const G4int copyNo, const G4VPhysicalVolume *physVol) const
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
void SetTranslation(const G4ThreeVector &v)
void SetOriginalParameters(G4PolyconeHistorical *pars)
void ComputeDimensions(G4Polycone &pcone, const G4int copyNo, const G4VPhysicalVolume *physVol) const
G4VParameterisationPolycone(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *msolid, DivisionType divType)
G4double GetR(G4double z, G4double z1, G4double r1, G4double z2, G4double r2) const
G4PolyconeHistorical * GetOriginalParameters() const
G4double GetStartPhi() const
G4bool Reset()
Definition: G4Polycone.cc:437
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4ParameterisationPolyconePhi(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *motherSolid, DivisionType divType)
EAxis
Definition: geomdefs.hh:54
G4int CalculateNDiv(G4double motherDim, G4double width, G4double offset) const
#define G4endl
Definition: G4ios.hh:61
void ComputeDimensions(G4Polycone &pcone, const G4int copyNo, const G4VPhysicalVolume *physVol) const
static PROLOG_HANDLER error
Definition: xmlrole.cc:108
double G4double
Definition: G4Types.hh:76
G4double GetRmin(G4double z, G4int nsegment) const
static const double mm
Definition: G4SIunits.hh:102
G4ParameterisationPolyconeZ(EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *motherSolid, DivisionType divType)
void ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const