source: trunk/source/geometry/divisions/src/G4ParameterisationPolyhedra.cc @ 1202

Last change on this file since 1202 was 1010, checked in by garnier, 15 years ago

update
File size: 20.8 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: G4ParameterisationPolyhedra.cc,v 1.16 2007/05/18 07:27:23 gcosmo Exp $
28// GEANT4 tag $Name: geant4-09-02-ref-02 $
29//
30// class G4ParameterisationPolyhedra Implementation file
31//
32// 14.10.03 - P.Arce, Initial version
33// 08.04.04 - I.Hrivnacova, Implemented reflection
34// --------------------------------------------------------------------
35
36#include "G4ParameterisationPolyhedra.hh"
37
38#include <iomanip>
39#include "G4ThreeVector.hh"
40#include "G4GeometryTolerance.hh"
41#include "G4RotationMatrix.hh"
42#include "G4VPhysicalVolume.hh"
43#include "G4LogicalVolume.hh"
44#include "G4ReflectedSolid.hh"
45#include "G4Polyhedra.hh"
46
47//--------------------------------------------------------------------------
48G4VParameterisationPolyhedra::
49G4VParameterisationPolyhedra( EAxis axis, G4int nDiv, G4double width,
50                              G4double offset, G4VSolid* msolid,
51                              DivisionType divType )
52  :  G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
53{
54  G4Polyhedra* msol = (G4Polyhedra*)(msolid);
55  if ((msolid->GetEntityType() != "G4ReflectedSolid") && (msol->IsGeneric()))
56  {
57    G4String message =
58         "Sorry, generic construct for G4Polyhedra NOT supported.\n Solid: "
59         + msol->GetName();
60    G4Exception("G4VParameterisationPolyhedra::G4VParameterisationPolyhedra()",
61                "NotSupported", FatalException, message);
62  }
63  if (msolid->GetEntityType() == "G4ReflectedSolid")
64  {
65    // Get constituent solid 
66    G4VSolid* mConstituentSolid
67       = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
68    msol = (G4Polyhedra*)(mConstituentSolid);
69 
70    // Get parameters
71    G4int   nofSides = msol->GetOriginalParameters()->numSide;
72    G4int   nofZplanes = msol->GetOriginalParameters()->Num_z_planes;
73    G4double* zValues  = msol->GetOriginalParameters()->Z_values;
74    G4double* rminValues  = msol->GetOriginalParameters()->Rmin;
75    G4double* rmaxValues  = msol->GetOriginalParameters()->Rmax;
76
77    // Invert z values,
78    // convert radius parameters
79    G4double* rminValues2 = new G4double[nofZplanes];
80    G4double* rmaxValues2 = new G4double[nofZplanes];
81    G4double* zValuesRefl = new G4double[nofZplanes];
82    for (G4int i=0; i<nofZplanes; i++)
83    {
84      rminValues2[i] = rminValues[i] * ConvertRadiusFactor(*msol);
85      rmaxValues2[i] = rmaxValues[i] * ConvertRadiusFactor(*msol);
86      zValuesRefl[i] = - zValues[i];
87    } 
88   
89    G4Polyhedra* newSolid
90      = new G4Polyhedra(msol->GetName(),
91                        msol->GetStartPhi(), 
92                        msol->GetEndPhi() - msol->GetStartPhi(),
93                        nofSides,
94                        nofZplanes, zValuesRefl, rminValues2, rmaxValues2);
95
96    delete [] rminValues2;       
97    delete [] rmaxValues2;       
98    delete [] zValuesRefl;       
99
100    msol = newSolid;
101    fmotherSolid = newSolid;
102    fReflectedSolid = true;
103    fDeleteSolid = true;
104  }   
105}
106
107//------------------------------------------------------------------------
108G4VParameterisationPolyhedra::~G4VParameterisationPolyhedra()
109{
110}
111
112//--------------------------------------------------------------------------
113G4double
114G4VParameterisationPolyhedra::
115ConvertRadiusFactor(const G4Polyhedra& phedra) const
116{
117  G4double phiTotal = phedra.GetEndPhi() - phedra.GetStartPhi();
118  G4int nofSides = phedra.GetOriginalParameters()->numSide;
119 
120  if ( (phiTotal <=0) || (phiTotal >
121        2*pi+G4GeometryTolerance::GetInstance()->GetAngularTolerance()) )
122    { phiTotal = 2*pi; }
123 
124  return std::cos(0.5*phiTotal/nofSides);
125} 
126
127//--------------------------------------------------------------------------
128G4ParameterisationPolyhedraRho::
129G4ParameterisationPolyhedraRho( EAxis axis, G4int nDiv,
130                               G4double width, G4double offset,
131                               G4VSolid* msolid, DivisionType divType )
132  :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType )
133{
134
135  CheckParametersValidity();
136  SetType( "DivisionPolyhedraRho" );
137
138  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
139  G4PolyhedraHistorical* original_pars = msol->GetOriginalParameters();
140
141  if( divType == DivWIDTH )
142  {
143    fnDiv = CalculateNDiv( original_pars->Rmax[0]
144                         - original_pars->Rmin[0], width, offset );
145  }
146  else if( divType == DivNDIV )
147  {
148    fwidth = CalculateWidth( original_pars->Rmax[0]
149                           - original_pars->Rmin[0], nDiv, offset );
150  }
151
152#ifdef G4DIVDEBUG
153  if( verbose >= 1 )
154  {
155    G4cout << " G4ParameterisationPolyhedraRho - # divisions " << fnDiv
156           << " = " << nDiv << G4endl
157           << " Offset " << foffset << " = " << offset << G4endl
158           << " Width " << fwidth << " = " << width << G4endl;
159  }
160#endif
161}
162
163//------------------------------------------------------------------------
164G4ParameterisationPolyhedraRho::~G4ParameterisationPolyhedraRho()
165{
166}
167
168//---------------------------------------------------------------------
169void G4ParameterisationPolyhedraRho::CheckParametersValidity()
170{
171  G4VDivisionParameterisation::CheckParametersValidity();
172
173  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
174
175  if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
176  {
177    G4cerr << "WARNING - "
178           << "G4ParameterisationPolyhedraRho::CheckParametersValidity()"
179           << G4endl
180           << "          Solid " << msol->GetName() << G4endl
181           << "          Division along R will be done with a width "
182           << "different for each solid section." << G4endl
183           << "          WIDTH will not be used !" << G4endl;
184  }
185  if( foffset != 0. )
186  {
187    G4cerr << "WARNING - "
188           << "G4ParameterisationPolyhedraRho::CheckParametersValidity()"
189           << G4endl
190           << "          Solid " << msol->GetName() << G4endl
191           << "          Division along  R will be done with a width "
192           << "different for each solid section." << G4endl
193           << "          OFFSET will not be used !" << G4endl;
194  }
195}
196
197//------------------------------------------------------------------------
198G4double G4ParameterisationPolyhedraRho::GetMaxParameter() const
199{
200  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
201  G4PolyhedraHistorical* original_pars = msol->GetOriginalParameters();
202  return original_pars->Rmax[0] - original_pars->Rmin[0];
203}
204
205//--------------------------------------------------------------------------
206void
207G4ParameterisationPolyhedraRho::
208ComputeTransformation( const G4int, G4VPhysicalVolume* physVol ) const
209{
210  //----- translation
211  G4ThreeVector origin(0.,0.,0.);
212
213  //----- set translation
214  physVol->SetTranslation( origin );
215
216  //----- calculate rotation matrix: unit
217
218#ifdef G4DIVDEBUG
219  if( verbose >= 2 )
220  {
221    G4cout << " G4ParameterisationPolyhedraRho " << G4endl
222           << " foffset: " << foffset/deg
223           << " - fwidth: " << fwidth/deg << G4endl;
224  }
225#endif
226
227  ChangeRotMatrix( physVol );
228
229#ifdef G4DIVDEBUG
230  if( verbose >= 2 )
231  {
232    G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraRho "
233           << G4endl
234           << " Position: " << origin
235           << " - Width: " << fwidth
236           << " - Axis: " << faxis  << G4endl;
237  }
238#endif
239}
240
241//--------------------------------------------------------------------------
242void
243G4ParameterisationPolyhedraRho::
244ComputeDimensions( G4Polyhedra& phedra, const G4int copyNo,
245                   const G4VPhysicalVolume* ) const
246{
247  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
248
249  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
250  G4PolyhedraHistorical origparam( *origparamMother );
251  G4int nZplanes = origparamMother->Num_z_planes;
252
253  G4double width = 0.;
254  for( G4int ii = 0; ii < nZplanes; ii++ )
255  {
256    width = CalculateWidth( origparamMother->Rmax[ii]
257                          - origparamMother->Rmin[ii], fnDiv, foffset );
258    origparam.Rmin[ii] = origparamMother->Rmin[ii]+foffset+width*copyNo;
259    origparam.Rmax[ii] = origparamMother->Rmin[ii]+foffset+width*(copyNo+1);
260  }
261
262  phedra.SetOriginalParameters(&origparam); // copy values & transfer pointers
263  phedra.Reset();                           // reset to new solid parameters
264
265#ifdef G4DIVDEBUG
266  if( verbose >= -2 )
267  {
268    G4cout << "G4ParameterisationPolyhedraRho::ComputeDimensions()" << G4endl
269           << "-- Parametrised phedra copy-number: " << copyNo << G4endl;
270    phedra.DumpInfo();
271  }
272#endif
273}
274
275//--------------------------------------------------------------------------
276G4ParameterisationPolyhedraPhi::
277G4ParameterisationPolyhedraPhi( EAxis axis, G4int nDiv,
278                               G4double width, G4double offset,
279                               G4VSolid* msolid, DivisionType divType )
280  :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType )
281{ 
282  CheckParametersValidity();
283  SetType( "DivisionPolyhedraPhi" );
284
285  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
286  G4double deltaPhi = msol->GetEndPhi() - msol->GetStartPhi();
287
288  if( divType == DivWIDTH )
289  {
290    fnDiv = CalculateNDiv( deltaPhi, width, offset );
291  }
292  else if( divType == DivNDIV )
293  {
294    fwidth = CalculateWidth( deltaPhi, nDiv, offset );
295  }
296
297#ifdef G4DIVDEBUG
298  if( verbose >= 1 )
299  {
300    G4cout << " G4ParameterisationPolyhedraPhi - # divisions " << fnDiv
301           << " = " << nDiv << G4endl
302           << " Offset " << foffset << " = " << offset << G4endl
303           << " Width " << fwidth << " = " << width << G4endl;
304  }
305#endif
306}
307
308//------------------------------------------------------------------------
309G4ParameterisationPolyhedraPhi::~G4ParameterisationPolyhedraPhi()
310{
311}
312
313//------------------------------------------------------------------------
314G4double G4ParameterisationPolyhedraPhi::GetMaxParameter() const
315{
316  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
317  return msol->GetEndPhi() - msol->GetStartPhi();
318}
319
320//---------------------------------------------------------------------
321void G4ParameterisationPolyhedraPhi::CheckParametersValidity()
322{
323  G4VDivisionParameterisation::CheckParametersValidity();
324
325  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
326
327  if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
328  {
329    G4cerr << "WARNING - "
330           << "G4ParameterisationPolyhedraPhi::CheckParametersValidity()"
331           << G4endl
332           << "          Solid " << msol->GetName() << G4endl
333           << "          Division along PHI will be done splitting "
334           << "in the defined numSide." << G4endl
335           << "          WIDTH will not be used !" << G4endl;
336  }
337  if( foffset != 0. )
338  {
339    G4cerr << "WARNING - "
340           << "G4ParameterisationPolyhedraPhi::CheckParametersValidity()"
341           << G4endl
342           << "          Solid " << msol->GetName() << G4endl
343           << "          Division along PHI will be done splitting "
344           << "in the defined numSide." << G4endl
345           << "          OFFSET will not be used !" << G4endl;
346  }
347
348  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
349
350  if( origparamMother->numSide != fnDiv )
351  { 
352    G4cerr << "ERROR - "
353           << "G4ParameterisationPolyhedraPhi::CheckParametersValidity()"
354           << G4endl
355           << "        Division along PHI will be done splitting in the defined"
356           << G4endl
357           << "        numSide, i.e, the number of division would be :"
358           << "        " << origparamMother->numSide
359           << " instead of " << fnDiv << " !"
360           << G4endl; 
361    G4Exception("G4ParameterisationPolyhedraPhi::CheckParametersValidity()",
362                "IllegalConstruct", FatalException,
363                "Not supported configuration.");
364  }
365}
366
367//--------------------------------------------------------------------------
368void
369G4ParameterisationPolyhedraPhi::
370ComputeTransformation( const G4int copyNo, G4VPhysicalVolume *physVol ) const
371{
372  //----- translation
373  G4ThreeVector origin(0.,0.,0.); 
374  //----- set translation
375  physVol->SetTranslation( origin );
376
377  //----- calculate rotation matrix (so that all volumes point to the centre)
378  G4double posi = copyNo*fwidth;
379
380#ifdef G4DIVDEBUG
381  if( verbose >= 2 )
382  {
383    G4cout << " G4ParameterisationPolyhedraPhi - position: " << posi/deg
384           << G4endl
385           << " copyNo: " << copyNo
386           << " - fwidth: " << fwidth/deg << G4endl;
387  }
388#endif
389
390  ChangeRotMatrix( physVol, -posi );
391
392#ifdef G4DIVDEBUG
393  if( verbose >= 2 )
394  {
395    G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraPhi " << copyNo
396           << G4endl
397           << " Position: " << origin << " - Width: " << fwidth
398           << " - Axis: " << faxis  << G4endl;
399  }
400#endif
401}
402
403//--------------------------------------------------------------------------
404void
405G4ParameterisationPolyhedraPhi::
406ComputeDimensions( G4Polyhedra& phedra, const G4int,
407                   const G4VPhysicalVolume* ) const
408{
409  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
410
411  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
412  G4PolyhedraHistorical origparam( *origparamMother );
413
414  origparam.numSide = 1;
415  origparam.Start_angle = origparamMother->Start_angle;
416  origparam.Opening_angle = fwidth;
417
418  phedra.SetOriginalParameters(&origparam);  // copy values & transfer pointers
419  phedra.Reset();                            // reset to new solid parameters
420
421#ifdef G4DIVDEBUG
422  if( verbose >= 2 )
423  {
424    G4cout << "G4ParameterisationPolyhedraPhi::ComputeDimensions():" << G4endl;
425    phedra.DumpInfo();
426  }
427#endif
428}
429
430//--------------------------------------------------------------------------
431G4ParameterisationPolyhedraZ::
432G4ParameterisationPolyhedraZ( EAxis axis, G4int nDiv,
433                             G4double width, G4double offset,
434                             G4VSolid* msolid, DivisionType divType )
435  :  G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType )
436{ 
437  CheckParametersValidity();
438  SetType( "DivisionPolyhedraZ" );
439
440  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
441  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
442 
443  if( divType == DivWIDTH )
444    {
445    fnDiv =
446      CalculateNDiv( origparamMother->Z_values[origparamMother->Num_z_planes-1]
447                     - origparamMother->Z_values[0] , width, offset );
448  }
449  else if( divType == DivNDIV )
450    {
451    fwidth =
452      CalculateNDiv( origparamMother->Z_values[origparamMother->Num_z_planes-1]
453                     - origparamMother->Z_values[0] , nDiv, offset );
454  }
455 
456#ifdef G4DIVDEBUG
457  if( verbose >= 1 )
458    {
459    G4cout << " G4ParameterisationPolyhedraZ - # divisions " << fnDiv << " = "
460           << nDiv << G4endl
461           << " Offset " << foffset << " = " << offset << G4endl
462           << " Width " << fwidth << " = " << width << G4endl;
463  }
464#endif
465}
466
467//---------------------------------------------------------------------
468G4ParameterisationPolyhedraZ::~G4ParameterisationPolyhedraZ()
469{
470}
471
472//------------------------------------------------------------------------
473G4double G4ParameterisationPolyhedraZ::GetMaxParameter() const
474{
475  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
476  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
477  return std::abs (origparamMother->Z_values[origparamMother->Num_z_planes-1]
478             -origparamMother->Z_values[0]);
479}
480
481//---------------------------------------------------------------------
482void G4ParameterisationPolyhedraZ::CheckParametersValidity()
483{
484  G4VDivisionParameterisation::CheckParametersValidity();
485
486  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
487
488  if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH )
489  {
490    G4cerr << "WARNING - "
491           << "G4ParameterisationPolyhedraZ::CheckParametersValidity()"
492           << G4endl
493           << "          Solid " << msol->GetName() << G4endl
494           << "          Division along Z will be done splitting "
495           << "in the defined z_planes." << G4endl
496           << "          WIDTH will not be used !" << G4endl;
497  }
498
499  if( foffset != 0. )
500  {
501    G4cerr << "WARNING - "
502           << "G4ParameterisationPolyhedraZ::CheckParametersValidity()"
503           << G4endl
504           << "          Solid " << msol->GetName() << G4endl
505           << "          Division along Z will be done splitting "
506           << "in the defined z_planes." << G4endl
507           << "          OFFSET will not be used !" << G4endl;
508  }
509
510  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
511
512  if( origparamMother->Num_z_planes-1 != fnDiv )
513  { 
514    G4cerr << "ERROR - "
515           << "G4ParameterisationPolyhedraZ::CheckParametersValidity()"
516           << G4endl
517           << "        Division along Z will be done splitting in the defined"
518           << G4endl
519           << "        z_planes, i.e, the number of division would be :"
520           << "        " << origparamMother->Num_z_planes-1
521           << " instead of " << fnDiv << " !"
522           << G4endl; 
523    G4Exception("G4ParameterisationPolyhedraZ::CheckParametersValidity()",
524                "IllegalConstruct", FatalException,
525                "Not supported configuration.");
526  }
527}
528
529//---------------------------------------------------------------------
530void
531G4ParameterisationPolyhedraZ::
532ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol) const
533{
534  G4Polyhedra* msol = (G4Polyhedra*)(GetMotherSolid());
535
536  //----- set translation: along Z axis
537  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
538  G4double posi = (origparamMother->Z_values[copyNo]
539                   + origparamMother->Z_values[copyNo+1])/2;
540  G4ThreeVector origin(0.,0.,posi); 
541  physVol->SetTranslation( origin );
542
543  //----- calculate rotation matrix: unit
544
545#ifdef G4DIVDEBUG
546  if( verbose >= 2 )
547  {
548    G4cout << " G4ParameterisationPolyhedraZ - position: " << posi << G4endl
549           << " copyNo: " << copyNo << " - foffset: " << foffset/deg
550           << " - fwidth: " << fwidth/deg << G4endl;
551  }
552#endif
553
554  ChangeRotMatrix( physVol );
555
556#ifdef G4DIVDEBUG
557  if( verbose >= 2 )
558  {
559    G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraZ "
560           << copyNo << G4endl
561           << " Position: " << origin << " - Width: " << fwidth
562           << " - Axis: " << faxis  << G4endl;
563  }
564#endif
565}
566
567//---------------------------------------------------------------------
568void
569G4ParameterisationPolyhedraZ::
570ComputeDimensions( G4Polyhedra& phedra, const G4int copyNo,
571                   const G4VPhysicalVolume* ) const
572{
573  // only for mother number of planes = 2!!
574  //
575  G4Polyhedra* msol = (G4Polyhedra*)(fmotherSolid);
576
577  G4PolyhedraHistorical* origparamMother = msol->GetOriginalParameters();
578  G4PolyhedraHistorical origparam( *origparamMother );
579
580  G4double posi = (origparamMother->Z_values[copyNo]
581                   + origparamMother->Z_values[copyNo+1])/2;
582
583  origparam.Num_z_planes = 2;
584  origparam.Z_values[0] = origparamMother->Z_values[copyNo] - posi;
585  origparam.Z_values[1] = origparamMother->Z_values[copyNo+1] - posi;
586  origparam.Rmin[0] = origparamMother->Rmin[copyNo];
587  origparam.Rmin[1] = origparamMother->Rmin[copyNo+1];
588  origparam.Rmax[0] = origparamMother->Rmax[copyNo];
589  origparam.Rmax[1] = origparamMother->Rmax[copyNo+1];
590
591  phedra.SetOriginalParameters(&origparam);  // copy values & transfer pointers
592  phedra.Reset();                            // reset to new solid parameters
593
594#ifdef G4DIVDEBUG
595  if( verbose >= 2 )
596  {
597    G4cout << "G4ParameterisationPolyhedraZ::ComputeDimensions()" << G4endl
598           << "-- Parametrised phedra copy-number: " << copyNo << G4endl;
599    phedra.DumpInfo();
600  }
601#endif
602}
603
Note: See TracBrowser for help on using the repository browser.