Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/divisions/src/G4ParameterisationPolyhedra.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // G4ParameterisationPolyhedra[Rho/Phi/Z] implementation
 27 //
 28 // 14.10.03 - P.Arce, Initial version
 29 // 08.04.04 - I.Hrivnacova, Implemented reflection
 30 // --------------------------------------------------------------------
 31 
 32 #include "G4ParameterisationPolyhedra.hh"
 33 
 34 #include <iomanip>
 35 #include "G4PhysicalConstants.hh"
 36 #include "G4ThreeVector.hh"
 37 #include "G4GeometryTolerance.hh"
 38 #include "G4RotationMatrix.hh"
 39 #include "G4VPhysicalVolume.hh"
 40 #include "G4LogicalVolume.hh"
 41 #include "G4ReflectedSolid.hh"
 42 #include "G4Polyhedra.hh"
 43 
 44 //--------------------------------------------------------------------------
 45 G4VParameterisationPolyhedra::
 46 G4VParameterisationPolyhedra( EAxis axis, G4int nDiv, G4double width,
 47                               G4double offset, G4VSolid* msolid,
 48                               DivisionType divType )
 49   :  G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
 50 {
 51   std::ostringstream message;
 52   /* #ifdef G4MULTITHREADED
 53   message << "Divisions for G4Polyhedra currently NOT supported in MT-mode."
 54           << G4endl
 55           << "Sorry! Solid: " << msolid->GetName();
 56   G4Exception("G4VParameterisationPolyhedra::G4VParameterisationPolyhedra()",
 57               "GeomDiv0001", FatalException, message);
 58         #endif */
 59 
 60   auto msol = (G4Polyhedra*)(msolid);
 61   if ((msolid->GetEntityType() != "G4ReflectedSolid") && (msol->IsGeneric()))
 62   {
 63     message << "Generic construct for G4Polyhedra NOT supported." << G4endl
 64             << "Sorry! Solid: " << msol->GetName();
 65     G4Exception("G4VParameterisationPolyhedra::G4VParameterisationPolyhedra()",
 66                 "GeomDiv0001", FatalException, message);
 67   }
 68   if (msolid->GetEntityType() == "G4ReflectedSolid")
 69   {
 70     // Get constituent solid
 71     //  
 72     G4VSolid* mConstituentSolid 
 73        = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
 74     msol = (G4Polyhedra*)(mConstituentSolid);
 75   
 76     // Get parameters
 77     //
 78     G4int nofSides = msol->GetOriginalParameters()->numSide;
 79     G4int nofZplanes = msol->GetOriginalParameters()->Num_z_planes;
 80     G4double* zValues = msol->GetOriginalParameters()->Z_values;
 81     G4double* rminValues = msol->GetOriginalParameters()->Rmin;
 82     G4double* rmaxValues = msol->GetOriginalParameters()->Rmax;
 83 
 84     // Invert z values, convert radius parameters
 85     //
 86     auto rminValues2 = new G4double[nofZplanes];
 87     auto rmaxValues2 = new G4double[nofZplanes];
 88     auto zValuesRefl = new G4double[nofZplanes];
 89     for (G4int i=0; i<nofZplanes; ++i)
 90     {
 91       rminValues2[i] = rminValues[i] * ConvertRadiusFactor(*msol);
 92       rmaxValues2[i] = rmaxValues[i] * ConvertRadiusFactor(*msol);
 93       zValuesRefl[i] = - zValues[i];
 94     }  
 95     
 96     auto newSolid
 97       = new G4Polyhedra(msol->GetName(),
 98                         msol->GetStartPhi(), 
 99                         msol->GetEndPhi() - msol->GetStartPhi(),
100                         nofSides,
101                         nofZplanes, zValuesRefl, rminValues2, rmaxValues2);
102 
103     delete [] rminValues2;       
104     delete [] rmaxValues2;       
105     delete [] zValuesRefl;       
106 
107     msol = newSolid;
108     fmotherSolid = newSolid;
109     fReflectedSolid = true;
110     fDeleteSolid = true;
111   }    
112 }
113 
114 //------------------------------------------------------------------------
115 G4VParameterisationPolyhedra::~G4VParameterisationPolyhedra() = default;
116 
117 //--------------------------------------------------------------------------
118 G4double 
119 G4VParameterisationPolyhedra::
120 ConvertRadiusFactor(const G4Polyhedra& phedra) const
121 {
122   G4double phiTotal = phedra.GetEndPhi() - phedra.GetStartPhi();
123   G4int nofSides = phedra.GetOriginalParameters()->numSide;
124   
125   if ( (phiTotal <=0) || (phiTotal >
126         2*pi+G4GeometryTolerance::GetInstance()->GetAngularTolerance()) )
127     { phiTotal = 2*pi; }
128   
129   return std::cos(0.5*phiTotal/nofSides);
130 }  
131 
132 //--------------------------------------------------------------------------
133 G4ParameterisationPolyhedraRho::
134 G4ParameterisationPolyhedraRho( EAxis axis, G4int nDiv,
135                                G4double width, G4double offset,
136                                G4VSolid* msolid, DivisionType divType )
137   :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType )
138 {
139 
140   CheckParametersValidity();
141   SetType( "DivisionPolyhedraRho" );
142 
143   auto msol = (G4Polyhedra*)(fmotherSolid);
144   G4PolyhedraHistorical* original_pars = msol->GetOriginalParameters();
145 
146   if( divType == DivWIDTH )
147   {
148     fnDiv = CalculateNDiv( original_pars->Rmax[0]
149                          - original_pars->Rmin[0], width, offset );
150   }
151   else if( divType == DivNDIV )
152   {
153     fwidth = CalculateWidth( original_pars->Rmax[0]
154                            - original_pars->Rmin[0], nDiv, offset );
155   }
156 
157 #ifdef G4DIVDEBUG
158   if( verbose >= 1 )
159   {
160     G4cout << " G4ParameterisationPolyhedraRho - # divisions " << fnDiv
161            << " = " << nDiv << G4endl
162            << " Offset " << foffset << " = " << offset << G4endl
163            << " Width " << fwidth << " = " << width << G4endl;
164   }
165 #endif
166 }
167 
168 //------------------------------------------------------------------------
169 G4ParameterisationPolyhedraRho::~G4ParameterisationPolyhedraRho() = default;
170 
171 //---------------------------------------------------------------------
172 void G4ParameterisationPolyhedraRho::CheckParametersValidity()
173 {
174   G4VDivisionParameterisation::CheckParametersValidity();
175 
176   auto msol = (G4Polyhedra*)(fmotherSolid);
177 
178   if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
179   {
180     std::ostringstream message;
181     message << "In solid " << msol->GetName() << G4endl
182             << "Division along R will be done with a width "
183             << "different for each solid section." << G4endl
184             << "WIDTH will not be used !";
185     G4Exception("G4ParameterisationPolyhedraRho::CheckParametersValidity()",
186                 "GeomDiv1001", JustWarning, message);
187   }
188   if( foffset != 0. )
189   {
190     std::ostringstream message;
191     message << "In solid " << msol->GetName() << G4endl
192             << "Division along  R will be done with a width "
193             << "different for each solid section." << G4endl
194             << "OFFSET will not be used !";
195     G4Exception("G4ParameterisationPolyhedraRho::CheckParametersValidity()",
196                 "GeomDiv1001", JustWarning, message);
197   }
198 }
199 
200 //------------------------------------------------------------------------
201 G4double G4ParameterisationPolyhedraRho::GetMaxParameter() const
202 {
203   auto msol = (G4Polyhedra*)(fmotherSolid);
204   G4PolyhedraHistorical* original_pars = msol->GetOriginalParameters();
205   return original_pars->Rmax[0] - original_pars->Rmin[0];
206 }
207 
208 //--------------------------------------------------------------------------
209 void
210 G4ParameterisationPolyhedraRho::
211 ComputeTransformation( const G4int, G4VPhysicalVolume* physVol ) const
212 {
213   //----- translation 
214   G4ThreeVector origin(0.,0.,0.);
215 
216   //----- set translation 
217   physVol->SetTranslation( origin );
218 
219   //----- calculate rotation matrix: unit
220 
221 #ifdef G4DIVDEBUG
222   if( verbose >= 2 )
223   {
224     G4cout << " G4ParameterisationPolyhedraRho " << G4endl
225            << " foffset: " << foffset/CLHEP::deg
226            << " - fwidth: " << fwidth/CLHEP::deg << G4endl;
227   }
228 #endif
229 
230   ChangeRotMatrix( physVol );
231 
232 #ifdef G4DIVDEBUG
233   if( verbose >= 2 )
234   {
235     G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraRho "
236            << G4endl
237            << " Position: " << origin
238            << " - Width: " << fwidth
239            << " - Axis: " << faxis  << G4endl;
240   }
241 #endif
242 }
243 
244 //--------------------------------------------------------------------------
245 void
246 G4ParameterisationPolyhedraRho::
247 ComputeDimensions( G4Polyhedra& phedra, const G4int copyNo,
248                    const G4VPhysicalVolume* ) const
249 {
250   auto msol = (G4Polyhedra*)(fmotherSolid);
251 
252   G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
253   G4PolyhedraHistorical origparam( *origparamMother );
254   G4int nZplanes = origparamMother->Num_z_planes;
255 
256   G4double width = 0.;
257   for( G4int ii = 0; ii < nZplanes; ++ii )
258   {
259     width = CalculateWidth( origparamMother->Rmax[ii]
260                           - origparamMother->Rmin[ii], fnDiv, foffset );
261     origparam.Rmin[ii] = origparamMother->Rmin[ii]+foffset+width*copyNo;
262     origparam.Rmax[ii] = origparamMother->Rmin[ii]+foffset+width*(copyNo+1);
263   }
264 
265   phedra.SetOriginalParameters(&origparam); // copy values & transfer pointers
266   phedra.Reset();                           // reset to new solid parameters
267 
268 #ifdef G4DIVDEBUG
269   if( verbose >= -2 )
270   {
271     G4cout << "G4ParameterisationPolyhedraRho::ComputeDimensions()" << G4endl
272            << "-- Parametrised phedra copy-number: " << copyNo << G4endl;
273     phedra.DumpInfo();
274   }
275 #endif
276 }
277 
278 //--------------------------------------------------------------------------
279 G4ParameterisationPolyhedraPhi::
280 G4ParameterisationPolyhedraPhi( EAxis axis, G4int nDiv,
281                                G4double width, G4double offset,
282                                G4VSolid* msolid, DivisionType divType )
283   :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType )
284 { 
285   CheckParametersValidity();
286   SetType( "DivisionPolyhedraPhi" );
287 
288   auto msol = (G4Polyhedra*)(fmotherSolid);
289   G4double deltaPhi = msol->GetEndPhi() - msol->GetStartPhi();
290 
291   if( divType == DivWIDTH )
292   {
293     fnDiv = msol->GetNumSide();
294   }
295 
296   fwidth = CalculateWidth( deltaPhi, fnDiv, 0.0 );
297 
298 #ifdef G4DIVDEBUG
299   if( verbose >= 1 )
300   {
301     G4cout << " G4ParameterisationPolyhedraPhi - # divisions " << fnDiv
302            << " = " << nDiv << G4endl
303            << " Offset " << foffset << " = " << offset << G4endl
304            << " Width " << fwidth << " = " << width << G4endl;
305   }
306 #endif
307 }
308 
309 //------------------------------------------------------------------------
310 G4ParameterisationPolyhedraPhi::~G4ParameterisationPolyhedraPhi() = default;
311 
312 //------------------------------------------------------------------------
313 G4double G4ParameterisationPolyhedraPhi::GetMaxParameter() const
314 {
315   auto msol = (G4Polyhedra*)(fmotherSolid);
316   return msol->GetEndPhi() - msol->GetStartPhi();
317 }
318 
319 //---------------------------------------------------------------------
320 void G4ParameterisationPolyhedraPhi::CheckParametersValidity()
321 {
322   G4VDivisionParameterisation::CheckParametersValidity();
323 
324   auto msol = (G4Polyhedra*)(fmotherSolid);
325 
326   if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
327   {
328     std::ostringstream message;
329     message << "In solid " << msol->GetName() << G4endl
330             << " Division along PHI will be done splitting "
331             << "in the defined numSide." << G4endl
332             << "WIDTH will not be used !";
333     G4Exception("G4ParameterisationPolyhedraPhi::CheckParametersValidity()",
334                 "GeomDiv1001", JustWarning, message);
335   }
336   if( foffset != 0. )
337   {
338     std::ostringstream message;
339     message << "In solid " << msol->GetName() << G4endl
340             << "Division along PHI will be done splitting "
341             << "in the defined numSide." << G4endl
342             << "OFFSET will not be used !";
343     G4Exception("G4ParameterisationPolyhedraPhi::CheckParametersValidity()",
344                 "GeomDiv1001", JustWarning, message);
345   }
346 
347   G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
348 
349   if( origparamMother->numSide != fnDiv &&  fDivisionType != DivWIDTH)
350   { 
351     std::ostringstream message;
352     message << "Configuration not supported." << G4endl
353             << "Division along PHI will be done splitting in the defined"
354             << G4endl
355             << "numSide, i.e, the number of division would be :"
356             << origparamMother->numSide << " instead of " << fnDiv << " !"; 
357     G4Exception("G4ParameterisationPolyhedraPhi::CheckParametersValidity()",
358                 "GeomDiv0001", FatalException, message);
359   }
360 }
361 
362 //--------------------------------------------------------------------------
363 void
364 G4ParameterisationPolyhedraPhi::
365 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
366 {
367   //----- translation 
368   G4ThreeVector origin(0.,0.,0.); 
369   //----- set translation 
370   physVol->SetTranslation( origin );
371 
372   //----- calculate rotation matrix (so that all volumes point to the centre)
373   G4double posi = copyNo*fwidth;
374 
375 #ifdef G4DIVDEBUG
376   if( verbose >= 2 )
377   {
378     G4cout << " G4ParameterisationPolyhedraPhi - position: " << posi/CLHEP::deg
379            << G4endl
380            << " copyNo: " << copyNo
381            << " - fwidth: " << fwidth/CLHEP::deg << G4endl;
382   }
383 #endif
384 
385   ChangeRotMatrix( physVol, -posi );
386 
387 #ifdef G4DIVDEBUG
388   if( verbose >= 2 )
389   {
390     G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraPhi "
391            << copyNo << G4endl
392            << " Position: " << origin << " - Width: " << fwidth
393            << " - Axis: " << faxis  << G4endl;
394   }
395 #endif
396 }
397 
398 //--------------------------------------------------------------------------
399 void
400 G4ParameterisationPolyhedraPhi::
401 ComputeDimensions( G4Polyhedra& phedra, const G4int,
402                    const G4VPhysicalVolume* ) const
403 {
404   auto msol = (G4Polyhedra*)(fmotherSolid);
405 
406   G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
407   G4PolyhedraHistorical origparam( *origparamMother );
408 
409   origparam.numSide = 1;
410   origparam.Start_angle = origparamMother->Start_angle;
411   origparam.Opening_angle = fwidth;
412 
413   phedra.SetOriginalParameters(&origparam);  // copy values & transfer pointers
414   phedra.Reset();                            // reset to new solid parameters
415 
416 #ifdef G4DIVDEBUG
417   if( verbose >= 2 )
418   {
419     G4cout << "G4ParameterisationPolyhedraPhi::ComputeDimensions():" << G4endl;
420     phedra.DumpInfo();
421   }
422 #endif
423 }
424 
425 //--------------------------------------------------------------------------
426 G4ParameterisationPolyhedraZ::
427 G4ParameterisationPolyhedraZ( EAxis axis, G4int nDiv,
428                              G4double width, G4double offset,
429                              G4VSolid* msolid, DivisionType divType )
430   :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType ),
431      fOrigParamMother(((G4Polyhedra*)fmotherSolid)->GetOriginalParameters())
432 { 
433   CheckParametersValidity();
434   SetType( "DivisionPolyhedraZ" );
435 
436   if( divType == DivWIDTH )
437   {
438     fnDiv =
439       CalculateNDiv(fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1]
440                     - fOrigParamMother->Z_values[0] , width, offset);
441   }
442   else if( divType == DivNDIV )
443   {
444     fwidth =
445       CalculateNDiv(fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1]
446                     - fOrigParamMother->Z_values[0] , nDiv, offset);
447   }
448   
449 #ifdef G4DIVDEBUG
450   if( verbose >= 1 )
451   {
452     G4cout << " G4ParameterisationPolyhedraZ - # divisions " << fnDiv << " = "
453            << nDiv << G4endl
454            << " Offset " << foffset << " = " << offset << G4endl
455            << " Width " << fwidth << " = " << width << G4endl;
456   }
457 #endif
458 }
459 
460 //---------------------------------------------------------------------
461 G4ParameterisationPolyhedraZ::~G4ParameterisationPolyhedraZ() = default;
462 
463 //------------------------------------------------------------------------
464 G4double G4ParameterisationPolyhedraZ::GetR(G4double z, 
465                                            G4double z1, G4double r1,  
466                                            G4double z2, G4double r2) const
467 {
468   // Linear parameterisation: 
469   // r = az + b
470   // a = (r1 - r2)/(z1-z2)
471   // b = r1 - a*z1
472 
473   return (r1-r2)/(z1-z2)*z + ( r1 - (r1-r2)/(z1-z2)*z1 ) ;
474 }  
475                                            
476 //------------------------------------------------------------------------
477 G4double G4ParameterisationPolyhedraZ::GetRmin(G4double z, G4int nseg) const
478 {
479 // Get Rmin in the given z position for the given polyhedra segment 
480 
481   return GetR(z, 
482               fOrigParamMother->Z_values[nseg], 
483               fOrigParamMother->Rmin[nseg],
484               fOrigParamMother->Z_values[nseg+1], 
485               fOrigParamMother->Rmin[nseg+1]);
486 }  
487                                            
488 //------------------------------------------------------------------------
489 G4double G4ParameterisationPolyhedraZ::GetRmax(G4double z, G4int nseg) const
490 {
491 // Get Rmax in the given z position for the given polyhedra segment 
492 
493   return GetR(z, 
494               fOrigParamMother->Z_values[nseg], 
495               fOrigParamMother->Rmax[nseg],
496               fOrigParamMother->Z_values[nseg+1], 
497               fOrigParamMother->Rmax[nseg+1]);
498 }  
499                                            
500 //------------------------------------------------------------------------
501 G4double G4ParameterisationPolyhedraZ::GetMaxParameter() const
502 {
503   return std::abs(fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1]
504                  -fOrigParamMother->Z_values[0]);
505 }
506 
507 //---------------------------------------------------------------------
508 void G4ParameterisationPolyhedraZ::CheckParametersValidity()
509 {
510   G4VDivisionParameterisation::CheckParametersValidity();
511 
512   // Division will be following the mother polyhedra segments
513   //
514   if( fDivisionType == DivNDIV )
515   {
516     if( fOrigParamMother->Num_z_planes-1 != fnDiv )
517     { 
518       std::ostringstream message;
519       message << "Configuration not supported." << G4endl
520               << "Division along Z will be done splitting in the defined"
521               << G4endl
522               << "Z planes, i.e, the number of division would be :"
523               << fOrigParamMother->Num_z_planes-1 << " instead of "
524               << fnDiv << " !"; 
525       G4Exception("G4ParameterisationPolyhedraZ::CheckParametersValidity()",
526                   "GeomDiv0001", FatalException, message);
527     }
528   }  
529 
530   // Division will be done within one polyhedra segment
531   // with applying given width and offset
532   //
533   if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
534   {
535     // Check if divided region does not span over more
536     // than one z segment
537 
538     G4int isegstart = -1;  // number of the segment containing start position
539     G4int isegend = -1;    // number of the segment containing end position
540 
541     if ( !fReflectedSolid )
542     {
543       // The start/end position of the divided region
544       //
545       G4double zstart = fOrigParamMother->Z_values[0] + foffset;
546       G4double zend = fOrigParamMother->Z_values[0]
547                     + foffset + fnDiv*fwidth;
548    
549       G4int counter = 0;
550       while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 )
551       {
552         // first segment
553         if ( zstart >= fOrigParamMother->Z_values[counter]  &&
554              zstart  < fOrigParamMother->Z_values[counter+1] )
555         {
556            isegstart = counter;
557         }     
558         // last segment
559         if ( zend  > fOrigParamMother->Z_values[counter] &&
560              zend <= fOrigParamMother->Z_values[counter+1] )
561         {
562            isegend = counter;
563         }   
564         ++counter;   
565       }  // Loop checking, 06.08.2015, G.Cosmo
566     }
567     else
568     {
569       // The start/end position of the divided region
570       //
571       G4double zstart = fOrigParamMother->Z_values[0] - foffset;
572       G4double zend = fOrigParamMother->Z_values[0]
573                     - (foffset + fnDiv* fwidth);
574    
575       G4int counter = 0;
576       while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 )
577       {
578         // first segment
579         if ( zstart <= fOrigParamMother->Z_values[counter]  &&
580              zstart  > fOrigParamMother->Z_values[counter+1] )
581         {
582            isegstart = counter;
583         }     
584         // last segment
585         if ( zend  < fOrigParamMother->Z_values[counter] &&
586              zend >= fOrigParamMother->Z_values[counter+1] )
587         {
588            isegend = counter;
589         }   
590         ++counter;   
591       }  // Loop checking, 06.08.2015, G.Cosmo
592     }
593   
594     if ( isegstart != isegend )
595     {
596       std::ostringstream message;
597       message << "Configuration not supported." << G4endl
598               << "Division with user defined width." << G4endl
599               << "Solid " << fmotherSolid->GetName() << G4endl
600               << "Divided region is not between two Z planes."; 
601       G4Exception("G4ParameterisationPolyhedraZ::CheckParametersValidity()",
602                   "GeomDiv0001", FatalException, message);
603     }
604   
605     fNSegment = isegstart;
606   }  
607 }
608 
609 //---------------------------------------------------------------------
610 void
611 G4ParameterisationPolyhedraZ::
612 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol) const
613 {
614   G4double posi;
615   if ( fDivisionType == DivNDIV )
616   {
617     // The position of the centre of copyNo-th mother polycone segment
618     
619     posi = ( fOrigParamMother->Z_values[copyNo]
620                     + fOrigParamMother->Z_values[copyNo+1])/2;
621     physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
622   }
623   
624   if ( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
625   {
626     // The position of the centre of copyNo-th division
627 
628     posi = fOrigParamMother->Z_values[0];
629     
630     if ( !fReflectedSolid )
631       posi += foffset + (2*copyNo + 1) * fwidth/2.;
632     else
633       posi -= foffset + (2*copyNo + 1) * fwidth/2.;
634     
635     physVol->SetTranslation( G4ThreeVector(0, 0, posi) );
636   }   
637 
638   //----- calculate rotation matrix: unit
639 
640 #ifdef G4DIVDEBUG
641   if( verbose >= 2 )
642   {
643     G4cout << " G4ParameterisationPolyhedraZ - position: " << posi << G4endl
644            << " copyNo: " << copyNo << " - foffset: " << foffset/CLHEP::deg
645            << " - fwidth: " << fwidth/CLHEP::deg << G4endl;
646   }
647 #endif
648 
649   ChangeRotMatrix( physVol );
650 
651 #ifdef G4DIVDEBUG
652   if( verbose >= 2 )
653   {
654     G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraZ "
655            << copyNo << G4endl
656            << " Position: (0,0,0) - Width: " << fwidth
657            << " - Axis: " << faxis  << G4endl;
658   }
659 #endif
660 }
661 
662 //---------------------------------------------------------------------
663 void
664 G4ParameterisationPolyhedraZ::
665 ComputeDimensions( G4Polyhedra& phedra, const G4int copyNo,
666                    const G4VPhysicalVolume* ) const
667 {
668   // Define division solid
669   //
670   G4PolyhedraHistorical origparam;
671   G4int nz = 2; 
672   origparam.Num_z_planes = nz;
673   origparam.numSide = fOrigParamMother->numSide;
674   origparam.Start_angle = fOrigParamMother->Start_angle;
675   origparam.Opening_angle = fOrigParamMother->Opening_angle;
676 
677   // Define division solid z sections
678   //
679   origparam.Z_values = new G4double[nz];
680   origparam.Rmin = new G4double[nz];
681   origparam.Rmax = new G4double[nz];
682   origparam.Z_values[0] = - fwidth/2.;
683   origparam.Z_values[1] = fwidth/2.;
684 
685   if ( fDivisionType == DivNDIV )
686   {
687     // The position of the centre of copyNo-th mother polycone segment
688     //
689     G4double posi = ( fOrigParamMother->Z_values[copyNo]
690                     + fOrigParamMother->Z_values[copyNo+1])/2;
691 
692     origparam.Z_values[0] = fOrigParamMother->Z_values[copyNo] - posi;
693     origparam.Z_values[1] = fOrigParamMother->Z_values[copyNo+1] - posi;
694     origparam.Rmin[0] = fOrigParamMother->Rmin[copyNo];
695     origparam.Rmin[1] = fOrigParamMother->Rmin[copyNo+1];
696     origparam.Rmax[0] = fOrigParamMother->Rmax[copyNo];
697     origparam.Rmax[1] = fOrigParamMother->Rmax[copyNo+1];
698   }  
699 
700   if ( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
701   {
702     if ( !fReflectedSolid )
703     {
704       origparam.Z_values[0] = -fwidth/2.;
705       origparam.Z_values[1] = fwidth/2.;
706 
707       // The position of the centre of copyNo-th division
708       //
709       G4double posi = fOrigParamMother->Z_values[0]
710                     + foffset + (2*copyNo + 1) * fwidth/2.;
711     
712       // The first and last z sides z values
713       G4double zstart = posi - fwidth/2.;
714       G4double zend = posi + fwidth/2.;
715       origparam.Rmin[0] = GetRmin(zstart, fNSegment); 
716       origparam.Rmax[0] = GetRmax(zstart, fNSegment);  
717       origparam.Rmin[1] = GetRmin(zend, fNSegment); 
718       origparam.Rmax[1] = GetRmax(zend, fNSegment); 
719     }
720     else
721     {   
722       origparam.Z_values[0] = fwidth/2.;
723       origparam.Z_values[1] = -fwidth/2.;
724 
725       // The position of the centre of copyNo-th division
726       //
727       G4double posi = fOrigParamMother->Z_values[0]
728                     - ( foffset + (2*copyNo + 1) * fwidth/2.);
729     
730       // The first and last z sides z values
731       //
732       G4double zstart = posi + fwidth/2.;
733       G4double zend = posi - fwidth/2.;
734       origparam.Rmin[0] = GetRmin(zstart, fNSegment); 
735       origparam.Rmax[0] = GetRmax(zstart, fNSegment);  
736       origparam.Rmin[1] = GetRmin(zend, fNSegment); 
737       origparam.Rmax[1] = GetRmax(zend, fNSegment); 
738     } 
739 
740     // It can happen due to rounding errors
741     //
742     if ( origparam.Rmin[0]    < 0.0 ) origparam.Rmin[0] = 0.0;
743     if ( origparam.Rmin[nz-1] < 0.0 ) origparam.Rmin[1] = 0.0;
744   }  
745 
746   phedra.SetOriginalParameters(&origparam);  // copy values & transfer pointers
747   phedra.Reset();                            // reset to new solid parameters
748 
749 #ifdef G4DIVDEBUG
750   if( verbose >= 2 )
751   {
752     G4cout << "G4ParameterisationPolyhedraZ::ComputeDimensions()" << G4endl
753            << "-- Parametrised phedra copy-number: " << copyNo << G4endl;
754     phedra.DumpInfo();
755   }
756 #endif
757 }
758