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 9.6.p3)


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