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