Geant4 Cross Reference |
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 68040 2013-03-13 14:19:04Z 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