Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/divisions/src/G4ParameterisationTrd.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 ]

Diff markup

Differences between /geometry/divisions/src/G4ParameterisationTrd.cc (Version 11.3.0) and /geometry/divisions/src/G4ParameterisationTrd.cc (Version 10.7.p4)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 // G4ParameterisationTrd[X/Y/Z] implementation     26 // G4ParameterisationTrd[X/Y/Z] implementation
 27 //                                                 27 //
 28 // 26.05.03 - P.Arce, Initial version              28 // 26.05.03 - P.Arce, Initial version
 29 // 08.04.04 - I.Hrivnacova, Implemented reflec     29 // 08.04.04 - I.Hrivnacova, Implemented reflection
 30 // 21.04.10 - M.Asai, Added gaps                   30 // 21.04.10 - M.Asai, Added gaps
 31 // -------------------------------------------     31 // --------------------------------------------------------------------
 32                                                    32 
 33 #include "G4ParameterisationTrd.hh"                33 #include "G4ParameterisationTrd.hh"
 34                                                    34 
 35 #include <iomanip>                                 35 #include <iomanip>
 36 #include "G4ThreeVector.hh"                        36 #include "G4ThreeVector.hh"
 37 #include "G4RotationMatrix.hh"                     37 #include "G4RotationMatrix.hh"
 38 #include "G4VPhysicalVolume.hh"                    38 #include "G4VPhysicalVolume.hh"
 39 #include "G4LogicalVolume.hh"                      39 #include "G4LogicalVolume.hh"
 40 #include "G4ReflectedSolid.hh"                     40 #include "G4ReflectedSolid.hh"
 41 #include "G4Trd.hh"                                41 #include "G4Trd.hh"
 42 #include "G4Trap.hh"                               42 #include "G4Trap.hh"
 43                                                    43 
 44 //--------------------------------------------     44 //--------------------------------------------------------------------------
 45 G4VParameterisationTrd::                           45 G4VParameterisationTrd::
 46 G4VParameterisationTrd( EAxis axis, G4int nDiv     46 G4VParameterisationTrd( EAxis axis, G4int nDiv, G4double width,
 47                         G4double offset, G4VSo     47                         G4double offset, G4VSolid* msolid,
 48                         DivisionType divType )     48                         DivisionType divType )
 49   :  G4VDivisionParameterisation( axis, nDiv,      49   :  G4VDivisionParameterisation( axis, nDiv, width, offset, divType, msolid )
 50 {                                                  50 {
 51   auto msol = (G4Trd*)(msolid);                <<  51   G4Trd* msol = (G4Trd*)(msolid);
 52   if (msolid->GetEntityType() == "G4ReflectedS     52   if (msolid->GetEntityType() == "G4ReflectedSolid")
 53   {                                                53   {
 54     // Get constituent solid                       54     // Get constituent solid  
 55     G4VSolid* mConstituentSolid                    55     G4VSolid* mConstituentSolid 
 56        = ((G4ReflectedSolid*)msolid)->GetConst     56        = ((G4ReflectedSolid*)msolid)->GetConstituentMovedSolid();
 57     msol = (G4Trd*)(mConstituentSolid);            57     msol = (G4Trd*)(mConstituentSolid);
 58                                                    58   
 59     // Create a new solid with inversed parame     59     // Create a new solid with inversed parameters
 60     auto newSolid                              <<  60     G4Trd* newSolid
 61       = new G4Trd(msol->GetName(),                 61       = new G4Trd(msol->GetName(),
 62                   msol->GetXHalfLength2(), mso     62                   msol->GetXHalfLength2(), msol->GetXHalfLength1(),
 63                   msol->GetYHalfLength2(), mso     63                   msol->GetYHalfLength2(), msol->GetYHalfLength1(),
 64                   msol->GetZHalfLength());         64                   msol->GetZHalfLength());
 65     msol = newSolid;                               65     msol = newSolid;
 66     fmotherSolid = newSolid;                       66     fmotherSolid = newSolid;
 67     fReflectedSolid = true;                        67     fReflectedSolid = true;
 68     fDeleteSolid = true;                           68     fDeleteSolid = true;
 69   }                                                69   }    
 70 }                                                  70 }
 71                                                    71 
 72 //--------------------------------------------     72 //------------------------------------------------------------------------
 73 G4VParameterisationTrd::~G4VParameterisationTr <<  73 G4VParameterisationTrd::~G4VParameterisationTrd()
                                                   >>  74 {
                                                   >>  75 }
 74                                                    76 
 75 //--------------------------------------------     77 //------------------------------------------------------------------------
 76 G4ParameterisationTrdX::                           78 G4ParameterisationTrdX::
 77 G4ParameterisationTrdX( EAxis axis, G4int nDiv     79 G4ParameterisationTrdX( EAxis axis, G4int nDiv,
 78                         G4double width, G4doub     80                         G4double width, G4double offset,
 79                         G4VSolid* msolid, Divi     81                         G4VSolid* msolid, DivisionType divType )
 80   :  G4VParameterisationTrd( axis, nDiv, width     82   :  G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
 81 {                                                  83 {
 82   CheckParametersValidity();                       84   CheckParametersValidity();
 83   SetType( "DivisionTrdX" );                       85   SetType( "DivisionTrdX" );
 84                                                    86 
 85   auto msol = (G4Trd*)(fmotherSolid);          <<  87   G4Trd* msol = (G4Trd*)(fmotherSolid);
 86   if( divType == DivWIDTH )                        88   if( divType == DivWIDTH )
 87   {                                                89   {
 88     fnDiv = CalculateNDiv( msol->GetXHalfLengt     90     fnDiv = CalculateNDiv( msol->GetXHalfLength1()+msol->GetXHalfLength2(),
 89                            width, offset );        91                            width, offset );
 90   }                                                92   }
 91   else if( divType == DivNDIV )                    93   else if( divType == DivNDIV )
 92   {                                                94   {
 93     fwidth = CalculateWidth( msol->GetXHalfLen     95     fwidth = CalculateWidth( msol->GetXHalfLength1()+msol->GetXHalfLength2(),
 94                              nDiv, offset );       96                              nDiv, offset );
 95   }                                                97   }
 96                                                    98 
 97 #ifdef G4DIVDEBUG                                  99 #ifdef G4DIVDEBUG
 98   if( verbose >= 1 )                              100   if( verbose >= 1 )
 99   {                                               101   {
100     G4cout << " G4ParameterisationTrdX - ## di    102     G4cout << " G4ParameterisationTrdX - ## divisions " << fnDiv << " = "
101            << nDiv << G4endl                      103            << nDiv << G4endl
102            << " Offset " << foffset << " = " <    104            << " Offset " << foffset << " = " << offset << G4endl
103            << " Width " << fwidth << " = " <<     105            << " Width " << fwidth << " = " << width << G4endl;
104   }                                               106   }
105 #endif                                            107 #endif
106                                                   108 
107   G4double mpDx1 = msol->GetXHalfLength1();       109   G4double mpDx1 = msol->GetXHalfLength1();
108   G4double mpDx2 = msol->GetXHalfLength2();       110   G4double mpDx2 = msol->GetXHalfLength2();
109   if( std::fabs(mpDx1 - mpDx2) > kCarTolerance    111   if( std::fabs(mpDx1 - mpDx2) > kCarTolerance )
110   {                                               112   {
111     bDivInTrap = true;                            113     bDivInTrap = true;
112   }                                               114   }
113 }                                                 115 }
114                                                   116 
115 //--------------------------------------------    117 //------------------------------------------------------------------------
116 G4ParameterisationTrdX::~G4ParameterisationTrd    118 G4ParameterisationTrdX::~G4ParameterisationTrdX()
117 = default;                                     << 119 {
                                                   >> 120 }
118                                                   121 
119 //--------------------------------------------    122 //------------------------------------------------------------------------
120 G4double G4ParameterisationTrdX::GetMaxParamet    123 G4double G4ParameterisationTrdX::GetMaxParameter() const
121 {                                                 124 {
122   auto msol = (G4Trd*)(fmotherSolid);          << 125   G4Trd* msol = (G4Trd*)(fmotherSolid);
123   return (msol->GetXHalfLength1()+msol->GetXHa    126   return (msol->GetXHalfLength1()+msol->GetXHalfLength2());
124 }                                                 127 }
125                                                   128 
126 //--------------------------------------------    129 //------------------------------------------------------------------------
127 void                                              130 void
128 G4ParameterisationTrdX::                          131 G4ParameterisationTrdX::
129 ComputeTransformation( const G4int copyNo,        132 ComputeTransformation( const G4int copyNo,
130                        G4VPhysicalVolume *phys    133                        G4VPhysicalVolume *physVol ) const
131 {                                                 134 {
132   auto msol = (G4Trd*)(fmotherSolid );         << 135   G4Trd* msol = (G4Trd*)(fmotherSolid );
133   G4double mdx = ( msol->GetXHalfLength1() + m    136   G4double mdx = ( msol->GetXHalfLength1() + msol->GetXHalfLength2() ) / 2.;
134   //----- translation                             137   //----- translation 
135   G4ThreeVector origin(0.,0.,0.);                 138   G4ThreeVector origin(0.,0.,0.); 
136   G4double posi;                                  139   G4double posi;
137   posi = -mdx + foffset + (copyNo+0.5)*fwidth;    140   posi = -mdx + foffset + (copyNo+0.5)*fwidth;
138   if( faxis == kXAxis )                           141   if( faxis == kXAxis )
139   {                                               142   {
140     origin.setX( posi );                          143     origin.setX( posi ); 
141   }                                               144   }
142   else                                            145   else
143   {                                               146   { 
144     std::ostringstream message;                   147     std::ostringstream message;
145     message << "Only axes along X are allowed     148     message << "Only axes along X are allowed !  Axis: " << faxis;
146     G4Exception("G4ParameterisationTrdX::Compu    149     G4Exception("G4ParameterisationTrdX::ComputeTransformation()",
147                 "GeomDiv0002", FatalException,    150                 "GeomDiv0002", FatalException, message);
148   }                                               151   }
149                                                   152 
150 #ifdef G4DIVDEBUG                                 153 #ifdef G4DIVDEBUG
151   if( verbose >= 2 )                              154   if( verbose >= 2 )
152   {                                               155   {
153     G4cout << std::setprecision(8)                156     G4cout << std::setprecision(8)
154            << " G4ParameterisationTrdX::Comput    157            << " G4ParameterisationTrdX::ComputeTransformation() "
155            << copyNo << G4endl                    158            << copyNo << G4endl
156            << " Position: " << origin << " - A    159            << " Position: " << origin << " - Axis: " << faxis << G4endl;
157   }                                               160   }
158 #endif                                            161 #endif
159                                                   162 
160    //----- set translation                        163    //----- set translation 
161   physVol->SetTranslation( origin );              164   physVol->SetTranslation( origin );
162 }                                                 165 }
163                                                   166 
164 //--------------------------------------------    167 //--------------------------------------------------------------------------
165 void                                              168 void
166 G4ParameterisationTrdX::                          169 G4ParameterisationTrdX::
167 ComputeDimensions( G4Trd& trd, [[maybe_unused] << 170 ComputeDimensions( G4Trd& trd, const G4int /* copyNo */,
                                                   >> 171                                const G4VPhysicalVolume* ) const
168 {                                                 172 {  
169   auto msol = (G4Trd*)(fmotherSolid);          << 173   G4Trd* msol = (G4Trd*)(fmotherSolid);
170   G4double pDy1 = msol->GetYHalfLength1();        174   G4double pDy1 = msol->GetYHalfLength1();
171   G4double pDy2 = msol->GetYHalfLength2();        175   G4double pDy2 = msol->GetYHalfLength2();
172   G4double pDz = msol->GetZHalfLength();          176   G4double pDz = msol->GetZHalfLength();
173   G4double pDx = fwidth/2. - fhgap;               177   G4double pDx = fwidth/2. - fhgap;
174                                                   178 
175   trd.SetAllParameters ( pDx, pDx, pDy1, pDy2,    179   trd.SetAllParameters ( pDx, pDx, pDy1, pDy2, pDz );
176                                                   180 
177 #ifdef G4DIVDEBUG                                 181 #ifdef G4DIVDEBUG
178   if( verbose >= 2 )                              182   if( verbose >= 2 )
179   {                                               183   {
180     G4cout << " G4ParameterisationTrdX::Comput    184     G4cout << " G4ParameterisationTrdX::ComputeDimensions():"
181            << copyNo << G4endl;                << 185 //           << copyNo
                                                   >> 186            << G4endl;
182     trd.DumpInfo();                               187     trd.DumpInfo();
183   }                                               188   }
184 #endif                                            189 #endif
185 }                                                 190 }
186                                                   191 
187 //--------------------------------------------    192 //--------------------------------------------------------------------------
188 void                                              193 void
189 G4ParameterisationTrdX::ComputeDimensions( G4T    194 G4ParameterisationTrdX::ComputeDimensions( G4Trap& trap, const G4int copyNo,
190                                            con    195                                            const G4VPhysicalVolume* ) const
191 {                                                 196 {
192   auto msol = (G4Trd*)(fmotherSolid);          << 197   G4Trd* msol = (G4Trd*)(fmotherSolid); 
193   G4double pDy1 = msol->GetYHalfLength1();        198   G4double pDy1 = msol->GetYHalfLength1();
194   G4double pDy2 = msol->GetYHalfLength2();        199   G4double pDy2 = msol->GetYHalfLength2();
195   G4double pDz = msol->GetZHalfLength();          200   G4double pDz = msol->GetZHalfLength();
196   //fwidth is at Z=0                              201   //fwidth is at Z=0
197   G4double pDx1 = msol->GetXHalfLength1();        202   G4double pDx1 = msol->GetXHalfLength1(); 
198   G4double pDx2 = msol->GetXHalfLength2();        203   G4double pDx2 = msol->GetXHalfLength2();
199   //  G4double pDxAVE = (pDx1+pDx2)/2.;           204   //  G4double pDxAVE = (pDx1+pDx2)/2.;
200   G4double xChangeRatio = (pDx2-pDx1) / (pDx2+    205   G4double xChangeRatio = (pDx2-pDx1) / (pDx2+pDx1);
201   G4double fWidChange = xChangeRatio*fwidth;      206   G4double fWidChange = xChangeRatio*fwidth;
202   G4double fWid1 = fwidth - fWidChange;           207   G4double fWid1 = fwidth - fWidChange;
203   G4double fWid2 = fwidth + fWidChange;           208   G4double fWid2 = fwidth + fWidChange;
204   G4double fOffsetChange = xChangeRatio*foffse    209   G4double fOffsetChange = xChangeRatio*foffset/2.;
205   G4double fOffset1 = foffset - fOffsetChange;    210   G4double fOffset1 = foffset - fOffsetChange;
206   G4double fOffset2 = foffset + fOffsetChange;    211   G4double fOffset2 = foffset + fOffsetChange;
207   G4double cxy1 = -pDx1+fOffset1 + (copyNo+0.5    212   G4double cxy1 = -pDx1+fOffset1 + (copyNo+0.5)*fWid1;// centre of the side at y=-pDy1
208   G4double cxy2 = -pDx2+fOffset2 + (copyNo+0.5    213   G4double cxy2 = -pDx2+fOffset2 + (copyNo+0.5)*fWid2;// centre of the side at y=+pDy1
209   G4double alp = std::atan( (cxy2-cxy1)/(pDz*2    214   G4double alp = std::atan( (cxy2-cxy1)/(pDz*2.) );
210                                                   215 
211   pDx1 = fwidth/2. - fWidChange/2.;               216   pDx1 = fwidth/2. - fWidChange/2.;
212   pDx2 = fwidth/2. + fWidChange/2.;               217   pDx2 = fwidth/2. + fWidChange/2.;
213                                                   218 
214                                                   219 
215   trap.SetAllParameters ( pDz,                    220   trap.SetAllParameters ( pDz,
216         alp,                                      221         alp,
217         0.,                                       222         0.,
218         pDy1,                                     223         pDy1,
219         pDx1,                                     224         pDx1,
220         pDx1,                                     225         pDx1,
221         0.,                                       226         0.,
222         pDy2,                                     227         pDy2,
223         pDx2 - fhgap,                             228         pDx2 - fhgap,
224         pDx2 - fhgap * pDx2/pDx1,                 229         pDx2 - fhgap * pDx2/pDx1,
225         0.);                                      230         0.);
226                                                   231 
227 #ifdef G4DIVDEBUG                                 232 #ifdef G4DIVDEBUG
228   if( verbose >= 2 )                              233   if( verbose >= 2 )
229   {                                               234   {
230     G4cout << " G4ParameterisationTrdX::Comput    235     G4cout << " G4ParameterisationTrdX::ComputeDimensions():"
231            << copyNo << G4endl;                   236            << copyNo << G4endl;
232     trap.DumpInfo();                              237     trap.DumpInfo();
233   }                                               238   }
234 #endif                                            239 #endif
235 }                                                 240 }
236                                                   241 
237                                                   242 
238 //--------------------------------------------    243 //--------------------------------------------------------------------------
239 G4ParameterisationTrdY::                          244 G4ParameterisationTrdY::
240 G4ParameterisationTrdY( EAxis axis, G4int nDiv    245 G4ParameterisationTrdY( EAxis axis, G4int nDiv,
241                         G4double width, G4doub    246                         G4double width, G4double offset,
242                         G4VSolid* msolid, Divi    247                         G4VSolid* msolid, DivisionType divType)
243   : G4VParameterisationTrd( axis, nDiv, width,    248   : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
244 {                                                 249 {
245   CheckParametersValidity();                      250   CheckParametersValidity();
246   SetType( "DivisionTrdY" );                      251   SetType( "DivisionTrdY" );
247                                                   252 
248   auto msol = (G4Trd*)(fmotherSolid);          << 253   G4Trd* msol = (G4Trd*)(fmotherSolid);
249   if( divType == DivWIDTH )                       254   if( divType == DivWIDTH )
250   {                                               255   {
251     fnDiv = CalculateNDiv( 2*msol->GetYHalfLen    256     fnDiv = CalculateNDiv( 2*msol->GetYHalfLength1(),
252                            width, offset );       257                            width, offset );
253   }                                               258   }
254   else if( divType == DivNDIV )                   259   else if( divType == DivNDIV )
255   {                                               260   {
256     fwidth = CalculateWidth( 2*msol->GetYHalfL    261     fwidth = CalculateWidth( 2*msol->GetYHalfLength1(),
257                              nDiv, offset );      262                              nDiv, offset );
258   }                                               263   }
259                                                   264 
260 #ifdef G4DIVDEBUG                                 265 #ifdef G4DIVDEBUG
261   if( verbose >= 1 )                              266   if( verbose >= 1 )
262   {                                               267   {
263      G4cout << " G4ParameterisationTrdY no div    268      G4cout << " G4ParameterisationTrdY no divisions " << fnDiv
264             << " = " << nDiv << G4endl            269             << " = " << nDiv << G4endl
265             << " Offset " << foffset << " = "     270             << " Offset " << foffset << " = " << offset << G4endl
266             << " width " << fwidth << " = " <<    271             << " width " << fwidth << " = " << width << G4endl;
267   }                                               272   }
268 #endif                                            273 #endif
269 }                                                 274 }
270                                                   275 
271 //--------------------------------------------    276 //------------------------------------------------------------------------
272 G4ParameterisationTrdY::~G4ParameterisationTrd << 277 G4ParameterisationTrdY::~G4ParameterisationTrdY()
                                                   >> 278 {
                                                   >> 279 }
273                                                   280 
274 //--------------------------------------------    281 //------------------------------------------------------------------------
275 G4double G4ParameterisationTrdY::GetMaxParamet    282 G4double G4ParameterisationTrdY::GetMaxParameter() const
276 {                                                 283 {
277   auto msol = (G4Trd*)(fmotherSolid);          << 284   G4Trd* msol = (G4Trd*)(fmotherSolid);
278   return (msol->GetYHalfLength1()+msol->GetYHa    285   return (msol->GetYHalfLength1()+msol->GetYHalfLength2());
279 }                                                 286 }
280                                                   287 
281 //--------------------------------------------    288 //--------------------------------------------------------------------------
282 void                                              289 void
283 G4ParameterisationTrdY::                          290 G4ParameterisationTrdY::
284 ComputeTransformation( const G4int copyNo, G4V    291 ComputeTransformation( const G4int copyNo, G4VPhysicalVolume* physVol ) const
285 {                                                 292 {
286   auto msol = (G4Trd*)(fmotherSolid );         << 293   G4Trd* msol = (G4Trd*)(fmotherSolid );
287   G4double mdy =  ( msol->GetYHalfLength1() +     294   G4double mdy =  ( msol->GetYHalfLength1() + msol->GetYHalfLength2() ) / 2.;
288                                                   295 
289   //----- translation                             296   //----- translation 
290   G4ThreeVector origin(0.,0.,0.);                 297   G4ThreeVector origin(0.,0.,0.); 
291   G4double posi = -mdy + foffset + (copyNo+0.5    298   G4double posi = -mdy + foffset + (copyNo+0.5)*fwidth;
292                                                   299 
293   if( faxis == kYAxis )                           300   if( faxis == kYAxis )
294   {                                               301   {
295     origin.setY( posi );                          302     origin.setY( posi ); 
296   }                                               303   }
297   else                                            304   else
298   {                                               305   { 
299     std::ostringstream message;                   306     std::ostringstream message;
300     message << "Only axes along Y are allowed     307     message << "Only axes along Y are allowed !  Axis: " << faxis;
301     G4Exception("G4ParameterisationTrdY::Compu    308     G4Exception("G4ParameterisationTrdY::ComputeTransformation()",
302                 "GeomDiv0002", FatalException,    309                 "GeomDiv0002", FatalException, message);
303   }                                               310   }
304                                                   311 
305 #ifdef G4DIVDEBUG                                 312 #ifdef G4DIVDEBUG
306   if( verbose >= 2 )                              313   if( verbose >= 2 )
307   {                                               314   {
308     G4cout << std::setprecision(8)                315     G4cout << std::setprecision(8)
309            << " G4ParameterisationTrdY::Comput    316            << " G4ParameterisationTrdY::ComputeTransformation " << copyNo
310            << " pos " << origin << " rot mat "    317            << " pos " << origin << " rot mat " << " axis " << faxis << G4endl;
311   }                                               318   }
312 #endif                                            319 #endif
313                                                   320 
314    //----- set translation                        321    //----- set translation 
315   physVol->SetTranslation( origin );              322   physVol->SetTranslation( origin );
316 }                                                 323 }
317                                                   324 
318 //--------------------------------------------    325 //--------------------------------------------------------------------------
319 void                                              326 void
320 G4ParameterisationTrdY::                          327 G4ParameterisationTrdY::
321 ComputeDimensions(G4Trd& trd, const G4int, con    328 ComputeDimensions(G4Trd& trd, const G4int, const G4VPhysicalVolume*) const
322 {                                                 329 {
323   //---- The division along Y of a Trd will re    330   //---- The division along Y of a Trd will result a Trd, only 
324   //--- if Y at -Z and +Z are equal, else use     331   //--- if Y at -Z and +Z are equal, else use the G4Trap version
325   auto msol = (G4Trd*)(fmotherSolid);          << 332   G4Trd* msol = (G4Trd*)(fmotherSolid);
326                                                   333   
327   G4double pDx1 = msol->GetXHalfLength1();        334   G4double pDx1 = msol->GetXHalfLength1();
328   G4double pDx2 = msol->GetXHalfLength2();        335   G4double pDx2 = msol->GetXHalfLength2();
329   G4double pDz = msol->GetZHalfLength();          336   G4double pDz = msol->GetZHalfLength();
330   G4double pDy = fwidth/2. - fhgap;               337   G4double pDy = fwidth/2. - fhgap;
331                                                   338  
332   trd.SetAllParameters ( pDx1, pDx2, pDy, pDy,    339   trd.SetAllParameters ( pDx1, pDx2, pDy, pDy, pDz );
333                                                   340 
334 #ifdef G4DIVDEBUG                                 341 #ifdef G4DIVDEBUG
335   if( verbose >= 2 )                              342   if( verbose >= 2 )
336   {                                               343   {
337     G4cout << " G4ParameterisationTrdY::Comput    344     G4cout << " G4ParameterisationTrdY::ComputeDimensions():" << G4endl;
338     trd.DumpInfo();                               345     trd.DumpInfo();
339   }                                               346   }
340 #endif                                            347 #endif
341 }                                                 348 }
342                                                   349 
343 //--------------------------------------------    350 //--------------------------------------------------------------------------
344 void                                              351 void
345 G4ParameterisationTrdY::ComputeDimensions( G4T    352 G4ParameterisationTrdY::ComputeDimensions( G4Trap& trap, const G4int copyNo,
346                                            con    353                                            const G4VPhysicalVolume* ) const
347 {                                                 354 {
348   auto msol = (G4Trd*)(fmotherSolid);          << 355   G4Trd* msol = (G4Trd*)(fmotherSolid); 
349   G4double pDx1 = msol->GetXHalfLength1();        356   G4double pDx1 = msol->GetXHalfLength1();
350   G4double pDx2 = msol->GetXHalfLength2();        357   G4double pDx2 = msol->GetXHalfLength2();
351   G4double pDz = msol->GetZHalfLength();          358   G4double pDz = msol->GetZHalfLength();
352   //fwidth is at Z=0                              359   //fwidth is at Z=0
353   G4double pDy1 = msol->GetYHalfLength1();        360   G4double pDy1 = msol->GetYHalfLength1(); 
354   G4double pDy2 = msol->GetYHalfLength2();        361   G4double pDy2 = msol->GetYHalfLength2();
355   //  G4double pDxAVE = (pDy1+pDy2)/2.;           362   //  G4double pDxAVE = (pDy1+pDy2)/2.;
356   G4double yChangeRatio = (pDy2-pDy1) / (pDy2+    363   G4double yChangeRatio = (pDy2-pDy1) / (pDy2+pDy1);
357   G4double fWidChange = yChangeRatio*fwidth;      364   G4double fWidChange = yChangeRatio*fwidth;
358   G4double fWid1 = fwidth - fWidChange;           365   G4double fWid1 = fwidth - fWidChange;
359   G4double fWid2 = fwidth + fWidChange;           366   G4double fWid2 = fwidth + fWidChange;
360   G4double fOffsetChange = yChangeRatio*foffse    367   G4double fOffsetChange = yChangeRatio*foffset/2.;
361   G4double fOffset1 = foffset - fOffsetChange;    368   G4double fOffset1 = foffset - fOffsetChange;
362   G4double fOffset2 = foffset + fOffsetChange;    369   G4double fOffset2 = foffset + fOffsetChange;
363   G4double cyx1 = -pDy1+fOffset1 + (copyNo+0.5    370   G4double cyx1 = -pDy1+fOffset1 + (copyNo+0.5)*fWid1;// centre of the side at y=-pDy1
364   G4double cyx2 = -pDy2+fOffset2 + (copyNo+0.5    371   G4double cyx2 = -pDy2+fOffset2 + (copyNo+0.5)*fWid2;// centre of the side at y=+pDy1
365   G4double alp = std::atan( (cyx2-cyx1)/(pDz*2    372   G4double alp = std::atan( (cyx2-cyx1)/(pDz*2.) );
366                                                   373 
367   pDy1 = fwidth/2. - fWidChange/2.;               374   pDy1 = fwidth/2. - fWidChange/2.;
368   pDy2 = fwidth/2. + fWidChange/2.;               375   pDy2 = fwidth/2. + fWidChange/2.;
369                                                   376 
370                                                   377 
371   trap.SetAllParameters ( pDz,                    378   trap.SetAllParameters ( pDz,
372         alp,                                      379         alp,
373         90*CLHEP::deg,                            380         90*CLHEP::deg,
374         pDy1,                                     381         pDy1,
375         pDx1,                                     382         pDx1,
376         pDx1,                                     383         pDx1,
377         0.,                                       384         0.,
378         pDy2,                                     385         pDy2,
379         pDx2 - fhgap,                             386         pDx2 - fhgap,
380         pDx2 - fhgap * pDx2/pDx1,                 387         pDx2 - fhgap * pDx2/pDx1,
381         0.);                                      388         0.);
382                                                   389 
383 #ifdef G4DIVDEBUG                                 390 #ifdef G4DIVDEBUG
384   if( verbose >= 2 )                              391   if( verbose >= 2 )
385   {                                               392   {
386     G4cout << " G4ParameterisationTrdY::Comput    393     G4cout << " G4ParameterisationTrdY::ComputeDimensions():"
387            << copyNo << G4endl;                   394            << copyNo << G4endl;
388     trap.DumpInfo();                              395     trap.DumpInfo();
389   }                                               396   }
390 #endif                                            397 #endif
391 }                                                 398 }
392                                                   399 
393 //--------------------------------------------    400 //--------------------------------------------------------------------------
394 G4ParameterisationTrdZ::                          401 G4ParameterisationTrdZ::
395 G4ParameterisationTrdZ( EAxis axis, G4int nDiv    402 G4ParameterisationTrdZ( EAxis axis, G4int nDiv,
396                         G4double width, G4doub    403                         G4double width, G4double offset,
397                         G4VSolid* msolid, Divi    404                         G4VSolid* msolid, DivisionType divType )
398   : G4VParameterisationTrd( axis, nDiv, width,    405   : G4VParameterisationTrd( axis, nDiv, width, offset, msolid, divType )
399 {                                                 406 { 
400   CheckParametersValidity();                      407   CheckParametersValidity();
401   SetType( "DivTrdZ" );                           408   SetType( "DivTrdZ" );
402                                                   409 
403   auto msol = (G4Trd*)(fmotherSolid);          << 410   G4Trd* msol = (G4Trd*)(fmotherSolid);
404   if( divType == DivWIDTH )                       411   if( divType == DivWIDTH )
405   {                                               412   {
406     fnDiv = CalculateNDiv( 2*msol->GetZHalfLen    413     fnDiv = CalculateNDiv( 2*msol->GetZHalfLength(),
407                            width, offset );       414                            width, offset );
408   }                                               415   }
409   else if( divType == DivNDIV )                   416   else if( divType == DivNDIV )
410   {                                               417   {
411     fwidth = CalculateWidth( 2*msol->GetZHalfL    418     fwidth = CalculateWidth( 2*msol->GetZHalfLength(),
412                              nDiv, offset );      419                              nDiv, offset );
413   }                                               420   }
414                                                   421 
415 #ifdef G4DIVDEBUG                                 422 #ifdef G4DIVDEBUG
416   if( verbose >= 1 )                              423   if( verbose >= 1 )
417   {                                               424   {
418     G4cout << " G4ParameterisationTrdZ no divi    425     G4cout << " G4ParameterisationTrdZ no divisions " << fnDiv
419            << " = " << nDiv << G4endl             426            << " = " << nDiv << G4endl
420            << " Offset " << foffset << " = " <    427            << " Offset " << foffset << " = " << offset << G4endl
421            << " Width " << fwidth << " = " <<     428            << " Width " << fwidth << " = " << width << G4endl;
422   }                                               429   }
423 #endif                                            430 #endif
424 }                                                 431 }
425                                                   432 
426 //--------------------------------------------    433 //------------------------------------------------------------------------
427 G4ParameterisationTrdZ::~G4ParameterisationTrd << 434 G4ParameterisationTrdZ::~G4ParameterisationTrdZ()
                                                   >> 435 {
                                                   >> 436 }
428                                                   437 
429 //--------------------------------------------    438 //------------------------------------------------------------------------
430 G4double G4ParameterisationTrdZ::GetMaxParamet    439 G4double G4ParameterisationTrdZ::GetMaxParameter() const
431 {                                                 440 {
432   auto msol = (G4Trd*)(fmotherSolid);          << 441   G4Trd* msol = (G4Trd*)(fmotherSolid);
433   return 2*msol->GetZHalfLength();                442   return 2*msol->GetZHalfLength();
434 }                                                 443 }
435                                                   444 
436 //--------------------------------------------    445 //--------------------------------------------------------------------------
437 void                                              446 void
438 G4ParameterisationTrdZ::                          447 G4ParameterisationTrdZ::
439 ComputeTransformation(const G4int copyNo, G4VP    448 ComputeTransformation(const G4int copyNo, G4VPhysicalVolume *physVol) const
440 {                                                 449 {
441   auto msol = (G4Trd*)(fmotherSolid );         << 450   G4Trd* msol = (G4Trd*)(fmotherSolid );
442   G4double mdz = msol->GetZHalfLength();          451   G4double mdz = msol->GetZHalfLength();
443                                                   452 
444   //----- translation                             453   //----- translation 
445   G4ThreeVector origin(0.,0.,0.);                 454   G4ThreeVector origin(0.,0.,0.); 
446   G4double posi = -mdz + OffsetZ() + (copyNo+0    455   G4double posi = -mdz + OffsetZ() + (copyNo+0.5)*fwidth;
447   if( faxis == kZAxis )                           456   if( faxis == kZAxis )
448   {                                               457   {
449     origin.setZ( posi );                          458     origin.setZ( posi ); 
450   }                                               459   }
451   else                                            460   else
452   {                                               461   { 
453     std::ostringstream message;                   462     std::ostringstream message;
454     message << "Only axes along Z are allowed     463     message << "Only axes along Z are allowed !  Axis: " << faxis;
455     G4Exception("G4ParameterisationTrdZ::Compu    464     G4Exception("G4ParameterisationTrdZ::ComputeTransformation()",
456                 "GeomDiv0002", FatalException,    465                 "GeomDiv0002", FatalException, message);
457   }                                               466   }
458                                                   467 
459 #ifdef G4DIVDEBUG                                 468 #ifdef G4DIVDEBUG
460   if( verbose >= 1 )                              469   if( verbose >= 1 )
461   {                                               470   {
462     G4cout << std::setprecision(8) << " G4Para    471     G4cout << std::setprecision(8) << " G4ParameterisationTrdZ: "
463            << copyNo << G4endl                    472            << copyNo << G4endl
464            << " Position: " << origin << " - O    473            << " Position: " << origin << " - Offset: " << foffset 
465            << " - Width: " << fwidth << " Axis    474            << " - Width: " << fwidth << " Axis " << faxis << G4endl;
466   }                                               475   }
467 #endif                                            476 #endif
468                                                   477 
469    //----- set translation                        478    //----- set translation 
470   physVol->SetTranslation( origin );              479   physVol->SetTranslation( origin );
471 }                                                 480 }
472                                                   481 
473 //--------------------------------------------    482 //--------------------------------------------------------------------------
474 void                                              483 void
475 G4ParameterisationTrdZ::                          484 G4ParameterisationTrdZ::
476 ComputeDimensions(G4Trd& trd, const G4int copy    485 ComputeDimensions(G4Trd& trd, const G4int copyNo,
477                   const G4VPhysicalVolume*) co    486                   const G4VPhysicalVolume*) const
478 {                                                 487 {
479   //---- The division along Z of a Trd will re    488   //---- The division along Z of a Trd will result a Trd
480   auto msol = (G4Trd*)(fmotherSolid);          << 489   G4Trd* msol = (G4Trd*)(fmotherSolid);
481                                                   490 
482   G4double pDx1 = msol->GetXHalfLength1();        491   G4double pDx1 = msol->GetXHalfLength1();
483   G4double DDx = (msol->GetXHalfLength2() - ms    492   G4double DDx = (msol->GetXHalfLength2() - msol->GetXHalfLength1() );
484   G4double pDy1 = msol->GetYHalfLength1();        493   G4double pDy1 = msol->GetYHalfLength1();
485   G4double DDy = (msol->GetYHalfLength2() - ms    494   G4double DDy = (msol->GetYHalfLength2() - msol->GetYHalfLength1() );
486   G4double pDz = fwidth/2. - fhgap;               495   G4double pDz = fwidth/2. - fhgap;
487   G4double zLength = 2*msol->GetZHalfLength();    496   G4double zLength = 2*msol->GetZHalfLength();
488                                                   497  
489   trd.SetAllParameters( pDx1+DDx*(OffsetZ()+co    498   trd.SetAllParameters( pDx1+DDx*(OffsetZ()+copyNo*fwidth+fhgap)/zLength,
490                         pDx1+DDx*(OffsetZ()+(c    499                         pDx1+DDx*(OffsetZ()+(copyNo+1)*fwidth-fhgap)/zLength, 
491                         pDy1+DDy*(OffsetZ()+co    500                         pDy1+DDy*(OffsetZ()+copyNo*fwidth+fhgap)/zLength,
492                         pDy1+DDy*(OffsetZ()+(c    501                         pDy1+DDy*(OffsetZ()+(copyNo+1)*fwidth-fhgap)/zLength,
493                         pDz );                    502                         pDz );
494                                                   503 
495 #ifdef G4DIVDEBUG                                 504 #ifdef G4DIVDEBUG
496   if( verbose >= 1 )                              505   if( verbose >= 1 )
497   {                                               506   {
498     G4cout << " G4ParameterisationTrdZ::Comput    507     G4cout << " G4ParameterisationTrdZ::ComputeDimensions()"
499            << " - Mother TRD " << G4endl;         508            << " - Mother TRD " << G4endl;
500     msol->DumpInfo();                             509     msol->DumpInfo();
501     G4cout << " - Parameterised TRD: "            510     G4cout << " - Parameterised TRD: "
502            << copyNo << G4endl;                   511            << copyNo << G4endl;
503     trd.DumpInfo();                               512     trd.DumpInfo();
504   }                                               513   }
505 #endif                                            514 #endif
506 }                                                 515 }
507                                                   516