BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
BesEmcConstruction.cc
Go to the documentation of this file.
1
2
3//---------------------------------------------------------------------------//
4// BOOST --- BESIII Object_Oreiented Simulation Tool //
5//---------------------------------------------------------------------------//
6// Descpirtion: EMC detector
7// Author: Fu Chengdong
8// Created: Sep 4, 2003
9// Comment:
10//---------------------------------------------------------------------------//
11//
12#include "G4ReflectionFactory.hh"
13#include "G4ThreeVector.hh"
14
15#include "EmcSim/BesCrystalParameterisation.hh"
16#include "EmcSim/BesEmcConstruction.hh"
17#include "EmcSim/BesEmcDetectorMessenger.hh"
18#include "EmcSim/BesEmcEndGeometry.hh"
19#include "EmcSim/BesEmcGeometry.hh"
20#include "EmcSim/BesEmcParameter.hh"
21#include "SimUtil/ReadBoostRoot.hh"
22
23#include "EmcSim/BesEmcSD.hh"
24#include "G4AffineTransform.hh"
25#include "G4Box.hh"
26#include "G4Color.hh"
27#include "G4Cons.hh"
28#include "G4FieldManager.hh"
29#include "G4IntersectionSolid.hh"
30#include "G4LogicalVolume.hh"
31#include "G4Material.hh"
32#include "G4PVParameterised.hh"
33#include "G4PVPlacement.hh"
34#include "G4PVReplica.hh"
35#include "G4Polycone.hh"
36#include "G4Polyhedra.hh"
37#include "G4RunManager.hh"
38#include "G4SDManager.hh"
39#include "G4Subscribers/G4IrregBox.hh"
40#include "G4SubtractionSolid.hh"
41#include "G4Transform3D.hh"
42#include "G4TransportationManager.hh"
43#include "G4Trap.hh"
44#include "G4Tubs.hh"
45#include "G4UniformMagField.hh"
46#include "G4UnionSolid.hh"
47#include "G4VPhysicalVolume.hh"
48#include "G4VisAttributes.hh"
49#include "globals.hh"
50
51#include "G4Geo/EmcG4Geo.h"
52#include "G4ios.hh"
53
54// Ma Qiumei Add
55#include <TError.h>
56#include <TGraphErrors.h>
57
58#define CHECKLV0 false
59#define CHECKLV1 false
60#define CHECKLV2 false
61#define CHECKLV3 false
62#define CHECKLV4 false
63#define CHECKLV5 false
64#define CHECKIrreg false
65
66BesEmcConstruction* BesEmcConstruction::fBesEmcConstruction = 0;
67
69
71 : verboseLevel( 0 )
72 , solidEMC( 0 )
73 , logicEMC( 0 )
74 , physiEMC( 0 )
75 , logicBSCWorld( 0 )
76 , solidBSCPhi( 0 )
77 , logicBSCPhi( 0 )
78 , physiBSCPhi( 0 )
79 , solidBSCTheta( 0 )
80 , logicBSCTheta( 0 )
81 , physiBSCTheta( 0 )
82 , solidBSCCrystal( 0 )
83 , logicBSCCrystal( 0 )
84 , physiBSCCrystal( 0 )
85 , magField( 0 )
86 , detectorMessenger( 0 )
87 , besEMCSD( 0 )
88 , crystalParam( 0 )
89 , logicEnd( 0 )
90 , logicEndPhi( 0 )
91 , logicEndCasing( 0 )
92 , logicEndCrystal( 0 )
93 , logicRear( 0 )
94 , logicRearCasing( 0 )
95 , logicOrgGlass( 0 )
96 , logicPD( 0 )
97 , logicAlPlate( 0 )
98 , logicPreAmpBox( 0 )
99 , logicAirInPABox( 0 )
100 , logicHangingPlate( 0 )
101 , logicOCGirder( 0 )
102 , logicCable( 0 )
103 , logicWaterPipe( 0 )
104 , logicSupportBar( 0 )
105 , logicSupportBar1( 0 )
106 , logicEndRing( 0 )
107 , logicGear( 0 )
108 , logicTaperRing1( 0 )
109 , logicTaperRing2( 0 )
110 , logicTaperRing3( 0 ) {
111 if ( fBesEmcConstruction )
112 {
113 G4cout << "BesEmcConstruction constructed twice." << G4endl;
114 exit( -1 );
115 }
116 fBesEmcConstruction = this;
117 // for debug
118 // G4Exception("BesEmcConstruction::BesEmcConstruction() starting........");
119 startID = 1;
120 phiNbCrystals = 0;
121 thetaNbCrystals = 0;
122 besEMCGeometry = new BesEmcGeometry();
123 emcEnd = new BesEmcEndGeometry();
124}
125
127 if ( detectorMessenger ) delete detectorMessenger;
128 if ( crystalParam ) delete crystalParam;
129 if ( besEMCGeometry ) delete besEMCGeometry;
130 if ( emcEnd ) delete emcEnd;
131
133}
134
135void BesEmcConstruction::Construct( G4LogicalVolume* logicBes ) {
136 gErrorIgnoreLevel = kFatal;
137 besEMCGeometry->ComputeEMCParameters();
138 detectorMessenger = new BesEmcDetectorMessenger( this, besEMCGeometry );
139 emcEnd->ComputeParameters();
140
141 G4SDManager* SDman = G4SDManager::GetSDMpointer();
142 if ( !besEMCSD )
143 {
144 besEMCSD = new BesEmcSD( "CalorSD", this, besEMCGeometry );
145 SDman->AddNewDetector( besEMCSD );
146 }
147
148 // Construction
149 G4cout << "--------ReadBoostRoot::GetEmc()=" << ReadBoostRoot::GetEmc() << G4endl;
150 if ( ReadBoostRoot::GetEmc() == 2 )
151 {
152 logicEMC = EmcG4Geo::Instance()->GetTopVolume();
153
154 if ( logicEMC )
155 {
156 physiEMC = new G4PVPlacement( 0, G4ThreeVector( 0.0, 0.0, 0.0 ), logicEMC, "physicalEMC",
157 logicBes, false, 0, CHECKLV0 );
158 G4cout << "logicEmc: === " << logicEMC << " physiEmc " << physiEMC << G4endl;
159
161 SetVisAndSD();
162 }
163 }
164 else
165 {
166 // for debug
167 // G4Exception("BesEmcConstruction::Construct() starting............");
168 //
169 DefineMaterials();
170 phiNbCrystals = ( *besEMCGeometry ).BSCNbPhi;
171 thetaNbCrystals = ( *besEMCGeometry ).BSCNbTheta * 2;
172
173 G4double da = 0.001 * deg; // delta angle to avoid overlap
174
175 //
176 // BSC
177 //
178 solidBSC = new G4Tubs(
179 "solidBSC", ( *besEMCGeometry ).TaperRingRmin1,
180 ( *besEMCGeometry ).BSCRmax + ( *besEMCGeometry ).SPBarThickness +
181 ( *besEMCGeometry ).SPBarThickness1 + 2.1 * mm, // radius from 942mm to 940 mm
182 ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 +
183 ( *besEMCGeometry ).EndRingDz,
184 0. * deg, 360. * deg );
185
186 solidESC = new G4Cons( "solidESC", ( *emcEnd ).WorldRmin1, ( *emcEnd ).WorldRmax1,
187 ( *emcEnd ).WorldRmin2, ( *emcEnd ).WorldRmax2,
188 ( *emcEnd ).WorldDz / 2, 0. * deg, 360. * deg );
189
190 solidEMC = new G4UnionSolid( "solidEMC0", solidBSC, solidESC, 0,
191 G4ThreeVector( 0, 0, ( *emcEnd ).WorldZPosition ) );
192
193 G4RotationMatrix* rotateESC = new G4RotationMatrix();
194 rotateESC->rotateY( 180. * deg );
195
196 solidEMC = new G4UnionSolid( "solidEMC", solidEMC, solidESC, rotateESC,
197 G4ThreeVector( 0, 0, -( *emcEnd ).WorldZPosition ) );
198
199 logicEMC = new G4LogicalVolume( solidEMC, G4Material::GetMaterial( "Air" ), "logicalEMC" );
200
201 physiEMC = new G4PVPlacement( 0, G4ThreeVector(), logicEMC, "physicalEMC", logicBes, false,
202 0, CHECKLV0 );
203
204 solidBSCWorld =
205 new G4SubtractionSolid( "solidBSCWorld0", solidBSC, solidESC, 0,
206 G4ThreeVector( 0, 0, ( *emcEnd ).WorldZPosition ) );
207
208 solidBSCWorld =
209 new G4SubtractionSolid( "solidBSCWorld", solidBSCWorld, solidESC, rotateESC,
210 G4ThreeVector( 0, 0, -( *emcEnd ).WorldZPosition ) );
211
212 logicBSCWorld = new G4LogicalVolume( solidBSCWorld, G4Material::GetMaterial( "Air" ),
213 "logicalBSCWorld" );
214
215 G4RotationMatrix* rotBSC = new G4RotationMatrix();
216 rotBSC->rotateY( 180. * deg );
217 physiBSCWorld = new G4PVPlacement( rotBSC, G4ThreeVector(), logicBSCWorld,
218 "physicalBSCWorld", logicEMC, false, 0, CHECKLV1 );
219
220 G4RotationMatrix* rotateMatrix[200];
221 G4double oOp, ox, oy, oz;
222 G4double delta = 0 * deg;
223 G4ThreeVector axis = G4ThreeVector( 0, 0, 0 );
224 oOp = ( *besEMCGeometry ).BSCRmin /
225 sin( 0.5 * ( *besEMCGeometry ).BSCPhiDphi + 90 * deg ) *
226 sin( ( *besEMCGeometry ).BSCAngleRotat );
227 G4double ll = ( *besEMCGeometry ).BSCCryLength;
228 G4double rr = ( *besEMCGeometry ).BSCRmin;
229 G4double oj = sqrt( ll * ll + rr * rr -
230 2 * ll * rr * cos( 180. * deg - ( *besEMCGeometry ).BSCAngleRotat ) );
231 G4double oij =
232 90. * deg - ( *besEMCGeometry ).BSCPhiDphi / 2. - ( *besEMCGeometry ).BSCAngleRotat;
233 G4double doj = asin( sin( 180. * deg - ( *besEMCGeometry ).BSCAngleRotat ) / oj * ll );
234 G4double ioj = ( *besEMCGeometry ).BSCPhiDphi / 2. + doj;
235 G4double ij = oj / sin( oij ) * sin( ioj );
236 G4double dOp = rr / sin( 90. * deg - ( *besEMCGeometry ).BSCPhiDphi / 2. ) *
237 sin( 90. * deg + ( *besEMCGeometry ).BSCPhiDphi / 2. -
238 ( *besEMCGeometry ).BSCAngleRotat );
239 G4double cOp = rr / sin( 90. * deg + ( *besEMCGeometry ).BSCPhiDphi / 2. ) *
240 sin( 90. * deg - ( *besEMCGeometry ).BSCPhiDphi / 2. -
241 ( *besEMCGeometry ).BSCAngleRotat );
242 G4double ch = ( dOp + ll ) / cos( ( *besEMCGeometry ).BSCPhiDphi ) - cOp;
243 G4double hi = ( dOp + ll ) * tan( ( *besEMCGeometry ).BSCPhiDphi ) - ij;
244 G4double oh = sqrt( ch * ch + rr * rr -
245 2 * ch * rr * cos( 180 * deg - ( *besEMCGeometry ).BSCAngleRotat ) );
246 G4double hoi = asin( sin( 180 * deg - oij ) / oh * hi );
247 G4double dok = asin( sin( 180 * deg - ( *besEMCGeometry ).BSCAngleRotat ) / oh * ch );
248 if ( verboseLevel > 3 )
249 G4cout << "oj=" << oj / cm << G4endl << "oij=" << oij / deg << G4endl
250 << "doj=" << doj / deg << G4endl << "ioj=" << ioj / deg << G4endl
251 << "ij=" << ij / cm << G4endl << "dOp=" << dOp / cm << G4endl
252 << "cOp=" << cOp / cm << G4endl << "ch=" << ch / cm << G4endl << "hi=" << hi / cm
253 << G4endl << "oh=" << oh / cm << G4endl << "hoi=" << hoi / deg << G4endl
254 << "dok=" << dok / deg << G4endl;
255
256 // Phi Cell
257 G4double cmo = asin( sin( 180 * degree - ( *besEMCGeometry ).BSCAngleRotat ) /
258 ( *besEMCGeometry ).BSCRmax2 * dOp );
259 G4double cm = ( *besEMCGeometry ).BSCRmax2 /
260 sin( 180 * degree - ( *besEMCGeometry ).BSCAngleRotat ) *
261 sin( ( *besEMCGeometry ).BSCAngleRotat - cmo );
262 G4ThreeVector Pout( dOp + cm * cos( ( *besEMCGeometry ).BSCAngleRotat ),
263 -cm * sin( ( *besEMCGeometry ).BSCAngleRotat ),
264 ( *besEMCGeometry ).BSCDz );
265
266 G4double rTaperRingOuter1 =
267 ( *besEMCGeometry ).TaperRingRmin1 + ( *besEMCGeometry ).TaperRingThickness1;
268 G4double rTaperRingOuter2 =
269 ( *besEMCGeometry ).TaperRingRmin2 + ( *besEMCGeometry ).TaperRingDr;
270 G4double zTaperRing1 = ( *besEMCGeometry ).BSCDz +
271 ( *besEMCGeometry ).TaperRingThickness3 -
272 ( *besEMCGeometry ).TaperRingDz;
273 G4double zTaperRing2 = ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3;
274
275 G4RotationMatrix* rotIntersection = new G4RotationMatrix();
276 rotIntersection->rotateZ( -( *besEMCGeometry ).BSCAngleRotat -
277 ( *besEMCGeometry ).BSCPhiDphi / 2. - hoi );
278 G4ThreeVector P( oOp * cos( -90. * deg + ( *besEMCGeometry ).BSCAngleRotat + hoi ),
279 oOp * sin( -90. * deg + ( *besEMCGeometry ).BSCAngleRotat + hoi ), 0 );
280 G4AffineTransform Td( rotIntersection, P );
281 G4ThreeVector md = Td.Inverse().TransformPoint( G4ThreeVector( 0, 0, 0 ) );
282 G4ThreeVector vd = Td.Inverse().TransformPoint( Pout );
283
284 G4double zPlane0[4] = { -( *besEMCGeometry ).BSCDz, -zTaperRing1, zTaperRing1,
285 ( *besEMCGeometry ).BSCDz };
286 G4double rInner0[4] = { vd.perp(), cOp, cOp, vd.perp() };
287 G4double rOuter0[4] = { ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).BSCPhiRmax,
288 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).BSCPhiRmax };
289
290 G4double zPlane1[4] = { -zTaperRing2, -zTaperRing1, zTaperRing1, zTaperRing2 };
291 G4double rInner1[4] = { rTaperRingOuter2, rTaperRingOuter1, rTaperRingOuter1,
292 rTaperRingOuter2 };
293 G4double rOuter1[4] = { ( *besEMCGeometry ).BSCRmax, ( *besEMCGeometry ).BSCRmax,
294 ( *besEMCGeometry ).BSCRmax, ( *besEMCGeometry ).BSCRmax };
295
296 G4VSolid* solidBSCPhi1 =
297 new G4Polycone( "solidBSCPhi1", 0, 360 * deg, 4, zPlane1, rInner1, rOuter1 );
298 G4VSolid* solidBSCPhi0 =
299 new G4Polyhedra( "solidBSCPhi0", 360. * deg - ( *besEMCGeometry ).BSCPhiDphi,
300 ( *besEMCGeometry ).BSCPhiDphi, ( *besEMCGeometry ).BSCNbPhi / 2, 4,
301 zPlane0, rInner0, rOuter0 );
302 G4IntersectionSolid* solidBSCPhi =
303 new G4IntersectionSolid( "solidBSCPhi", solidBSCPhi0, solidBSCPhi1, 0, md );
304
305 logicBSCPhi =
306 new G4LogicalVolume( solidBSCPhi, G4Material::GetMaterial( "Air" ), "logicalBSCPhi" );
307
308 G4int i;
309 for ( G4int j = 0; j < ( *besEMCGeometry ).BSCNbPhi; j++ ) //=============
310 {
311 if ( j < ( *besEMCGeometry ).BSCNbPhi / 2 )
312 { // 0~59
313 i = ( *besEMCGeometry ).BSCNbPhi / 2 - j - 1;
314 }
315 else
316 { // 60~119
317 i = ( *besEMCGeometry ).BSCNbPhi * 3 / 2 - j - 1;
318 }
319 rotateMatrix[i] = new G4RotationMatrix();
320 rotateMatrix[i]->rotateZ( -i * ( *besEMCGeometry ).BSCPhiDphi -
321 ( *besEMCGeometry ).BSCAngleRotat -
322 ( *besEMCGeometry ).BSCPhiDphi / 2. - hoi );
323 rotateMatrix[i]->getAngleAxis( delta, axis );
324 // G4cout << "The axis of crystals in the world system is: "
325 // << delta/deg << "(deg)(delta) "
326 //<< axis << "(Z axis)"<< G4endl;
327 ox = oOp * cos( -90. * deg + ( *besEMCGeometry ).BSCAngleRotat + hoi +
328 i * ( *besEMCGeometry ).BSCPhiDphi );
329 oy = oOp * sin( -90. * deg + ( *besEMCGeometry ).BSCAngleRotat + hoi +
330 i * ( *besEMCGeometry ).BSCPhiDphi );
331 oz = 0 * cm;
332
333 ostringstream strPhi;
334 strPhi << "physicalBSCPhi" << j;
335
336 physiBSCPhi =
337 new G4PVPlacement( rotateMatrix[i], G4ThreeVector( ox, oy, oz ), logicBSCPhi,
338 strPhi.str(), logicBSCWorld, false, j, CHECKLV2 );
339 // G4cout << G4ThreeVector(ox/cm,oy/cm,oz/cm) <<"(cm)" << G4endl
340 // << (-(*besEMCGeometry).BSCAngleRotat+(i-1)*(*besEMCGeometry).BSCPhiDphi)/deg
341 // <<"(degree)" << G4endl;
342 }
343
344 //
345 // Crystals
346 //
347 G4double zHalfLength[50];
348 G4double thetaAxis[50];
349 G4double phiAxis[50];
350 G4double yHalfLength1[50];
351 G4double xHalfLength2[50];
352 G4double xHalfLength1[50];
353 G4double tanAlpha1[50];
354 G4double yHalfLength2[50];
355 G4double xHalfLength4[50];
356 G4double xHalfLength3[50];
357 G4double tanAlpha2[50];
358 G4double xPosition[50];
359 G4double yPosition[50];
360 G4double zPosition[50];
361 G4double thetaPosition[50];
362 for ( i = 0; i < ( *besEMCGeometry ).BSCNbTheta; i++ )
363 {
364 zHalfLength[i] = ( *besEMCGeometry ).zHalfLength[i];
365 thetaAxis[i] = ( *besEMCGeometry ).thetaAxis[i];
366 phiAxis[i] = ( *besEMCGeometry ).phiAxis[i];
367 yHalfLength1[i] = ( *besEMCGeometry ).yHalfLength1[i];
368 xHalfLength2[i] = ( *besEMCGeometry ).xHalfLength2[i];
369 xHalfLength1[i] = ( *besEMCGeometry ).xHalfLength1[i];
370 tanAlpha1[i] = ( *besEMCGeometry ).tanAlpha1[i];
371 yHalfLength2[i] = ( *besEMCGeometry ).yHalfLength2[i];
372 xHalfLength4[i] = ( *besEMCGeometry ).xHalfLength4[i];
373 xHalfLength3[i] = ( *besEMCGeometry ).xHalfLength3[i];
374 tanAlpha2[i] = ( *besEMCGeometry ).tanAlpha2[i];
375 xPosition[i] = ( *besEMCGeometry ).xPosition[i];
376 yPosition[i] = ( *besEMCGeometry ).yPosition[i];
377 zPosition[i] = ( *besEMCGeometry ).zPosition[i];
378 thetaPosition[i] = ( *besEMCGeometry ).thetaPosition[i];
379 if ( verboseLevel > 4 )
380 G4cout << "The sizes of the " << i + 1 << " crystal are:" << G4endl
381 << "zHalfLength =" << zHalfLength[i] / cm << "(cm)," << G4endl
382 << "thetaAxis =" << thetaAxis[i] / deg << "(deg)," << G4endl
383 << "phiAxis =" << phiAxis[i] / deg << "(deg)," << G4endl
384 << "yHalfLength1=" << yHalfLength1[i] / cm << "(cm)," << G4endl
385 << "xHalfLength1=" << xHalfLength1[i] / cm << "(cm)," << G4endl
386 << "xHalfLength2=" << xHalfLength2[i] / cm << "(cm)," << G4endl
387 << "tanAlpha1 =" << tanAlpha1[i] << G4endl
388 << "yHalfLength2=" << yHalfLength2[i] / cm << "(cm)," << G4endl
389 << "xHalfLength3=" << xHalfLength3[i] / cm << "(cm)," << G4endl
390 << "xHalfLength4=" << xHalfLength4[i] / cm << "(cm)," << G4endl
391 << "tanAlpha2 =" << tanAlpha2[i] << "." << G4endl;
392 }
393 besEMCGeometry->ModifyForCasing();
394
395 solidBSCCrystal =
396 new G4Trap( "solidCrystal", 100 * cm, 100 * deg, 100 * deg, 100 * cm, 100 * cm,
397 100 * cm, 100 * deg, 100 * cm, 100 * cm, 100 * cm, 100 * deg );
398
399 logicBSCCrystal =
400 new G4LogicalVolume( solidBSCCrystal, fCrystalMaterial, "logicalCrystal" );
401
402 crystalParam = new BesCrystalParameterisation( startID, thetaNbCrystals,
403 ( *besEMCGeometry ).BSCNbTheta * 2,
404 besEMCGeometry, verboseLevel );
405
406 //---------------------------------------------------------------------------------
407 // rear substance
408 solidRear =
409 new G4Box( "solidRearBox", ( *besEMCGeometry ).rearBoxLength / 2,
410 ( *besEMCGeometry ).rearBoxLength / 2, ( *besEMCGeometry ).rearBoxDz / 2 );
411
412 logicRear =
413 new G4LogicalVolume( solidRear, G4Material::GetMaterial( "Air" ), "logicalRearBox" );
414
415 // organic glass
416 solidOrgGlass = new G4Box( "solidOrganicGlass", ( *besEMCGeometry ).orgGlassLengthX / 2,
417 ( *besEMCGeometry ).orgGlassLengthY / 2,
418 ( *besEMCGeometry ).orgGlassLengthZ / 2 );
419
420 logicOrgGlass = new G4LogicalVolume( solidOrgGlass, organicGlass, "logicalOrganicGlass" );
421
422 physiOrgGlass = new G4PVPlacement(
423 0,
424 G4ThreeVector(
425 0, 0,
426 -( ( *besEMCGeometry ).rearBoxDz - ( *besEMCGeometry ).orgGlassLengthZ ) / 2 ),
427 logicOrgGlass, "physicalOrganicGlass", logicRear, false, 0, CHECKLV4 );
428
429 // casing
430 solidCasingBox = new G4Box( "solidCasingBox", ( *besEMCGeometry ).rearBoxLength / 2,
431 ( *besEMCGeometry ).rearBoxLength / 2,
432 ( *besEMCGeometry ).rearCasingThickness / 2 );
433
434 solidAirHole =
435 new G4Box( "solidAirHole", ( *besEMCGeometry ).orgGlassLengthX / 2,
436 ( *besEMCGeometry ).orgGlassLengthY / 2,
437 ( *besEMCGeometry ).rearBoxDz / 2 ); // any value more than casing thickness
438
439 solidRearCasing = new G4SubtractionSolid( "solidRearCasing", solidCasingBox, solidAirHole,
440 0, G4ThreeVector() );
441
442 logicRearCasing =
443 new G4LogicalVolume( solidRearCasing, rearCasingMaterial, "logicalRearCasing" );
444
445 physiRearCasing = new G4PVPlacement(
446 0,
447 G4ThreeVector(
448 0, 0,
449 -( ( *besEMCGeometry ).rearBoxDz - ( *besEMCGeometry ).rearCasingThickness ) / 2 ),
450 logicRearCasing, "physicalRearCasing", logicRear, false, 0, CHECKLV4 );
451
452 // Al Plate
453 solidAlBox =
454 new G4Box( "solidAlBox", ( *besEMCGeometry ).rearBoxLength / 2,
455 ( *besEMCGeometry ).rearBoxLength / 2, ( *besEMCGeometry ).AlPlateDz / 2 );
456
457 solidAlPlate =
458 new G4SubtractionSolid( "solidAlPlate", solidAlBox, solidAirHole, 0, G4ThreeVector() );
459
460 logicAlPlate = new G4LogicalVolume( solidAlPlate, G4Material::GetMaterial( "Aluminium" ),
461 "logicalAlPlate" );
462
463 physiAlPlate =
464 new G4PVPlacement( 0,
465 G4ThreeVector( 0, 0,
466 -( ( *besEMCGeometry ).rearBoxDz / 2 -
467 ( *besEMCGeometry ).rearCasingThickness -
468 ( *besEMCGeometry ).AlPlateDz / 2 ) ),
469 logicAlPlate, "physicalAlPlate", logicRear, false, 0, CHECKLV4 );
470
471 // photodiode
472 solidPD =
473 new G4Box( "solidPD",
474 ( *besEMCGeometry ).PDLengthX, // two PD
475 ( *besEMCGeometry ).PDLengthY / 2, ( *besEMCGeometry ).PDLengthZ / 2 );
476
477 logicPD =
478 new G4LogicalVolume( solidPD, G4Material::GetMaterial( "M_Silicon" ), "logicalPD" );
479
480 physiPD = new G4PVPlacement( 0,
481 G4ThreeVector( 0, 0,
482 -( ( *besEMCGeometry ).rearBoxDz / 2 -
483 ( *besEMCGeometry ).orgGlassLengthZ -
484 ( *besEMCGeometry ).PDLengthZ / 2 ) ),
485 logicPD, "physicalPD", logicRear, false, 0, CHECKLV4 );
486
487 // preamplifier box
488 solidPreAmpBox =
489 new G4Box( "solidPreAmpBox", ( *besEMCGeometry ).rearBoxLength / 2,
490 ( *besEMCGeometry ).rearBoxLength / 2, ( *besEMCGeometry ).PABoxDz / 2 );
491
492 logicPreAmpBox = new G4LogicalVolume(
493 solidPreAmpBox, G4Material::GetMaterial( "Aluminium" ), "logicalPreAmpBox" );
494
495 physiPreAmpBox = new G4PVPlacement(
496 0,
497 G4ThreeVector( 0, 0,
498 -( ( *besEMCGeometry ).rearBoxDz / 2 -
499 ( *besEMCGeometry ).rearCasingThickness -
500 ( *besEMCGeometry ).AlPlateDz - ( *besEMCGeometry ).PABoxDz / 2 ) ),
501 logicPreAmpBox, "physicalPreAmpBox", logicRear, false, 0, CHECKLV4 );
502
503 // air in preamplifier box
504 solidAirInPABox =
505 new G4Box( "solidAirInPABox",
506 ( *besEMCGeometry ).rearBoxLength / 2 - ( *besEMCGeometry ).PABoxThickness,
507 ( *besEMCGeometry ).rearBoxLength / 2 - ( *besEMCGeometry ).PABoxThickness,
508 ( *besEMCGeometry ).PABoxDz / 2 - ( *besEMCGeometry ).PABoxThickness );
509
510 logicAirInPABox = new G4LogicalVolume( solidAirInPABox, G4Material::GetMaterial( "Air" ),
511 "logicalAirInPABox" );
512
513 physiAirInPABox =
514 new G4PVPlacement( 0, G4ThreeVector(), logicAirInPABox, "physicalAirInPABox",
515 logicPreAmpBox, false, 0, CHECKLV5 );
516
517 // stainless steel for hanging the crystal
518 solidHangingPlate = new G4Box( "solidHangingPlate", ( *besEMCGeometry ).rearBoxLength / 2,
519 ( *besEMCGeometry ).rearBoxLength / 2,
520 ( *besEMCGeometry ).HangingPlateDz / 2 );
521
522 logicHangingPlate =
523 new G4LogicalVolume( solidHangingPlate, stainlessSteel, "logicalHangingPlate" );
524
525 physiHangingPlate = new G4PVPlacement(
526 0,
527 G4ThreeVector( 0, 0,
528 -( ( *besEMCGeometry ).rearBoxDz / 2 -
529 ( *besEMCGeometry ).rearCasingThickness -
530 ( *besEMCGeometry ).AlPlateDz - ( *besEMCGeometry ).PABoxDz -
531 ( *besEMCGeometry ).HangingPlateDz / 2 ) ),
532 logicHangingPlate, "physicalHangingPlate", logicRear, false, 0, CHECKLV4 );
533
534 // water pipe
535 solidWaterPipe = new G4Tubs( "solidWaterPipe", 0, ( *besEMCGeometry ).waterPipeDr,
536 ( *besEMCGeometry ).BSCDz, 0. * deg, 360. * deg );
537
538 logicWaterPipe = new G4LogicalVolume( solidWaterPipe, waterPipe, "logicalWaterPipe" );
539
540 physiWaterPipe = new G4PVPlacement(
541 0,
542 G4ThreeVector( ( *besEMCGeometry ).cablePosX[0] - 2 * ( *besEMCGeometry ).cableDr,
543 ( *besEMCGeometry ).cablePosY[0] - ( *besEMCGeometry ).cableDr -
544 ( *besEMCGeometry ).waterPipeDr,
545 0 ),
546 logicWaterPipe, "physicalWaterPipe", logicBSCPhi, false, 0, CHECKLV3 );
547 //---------------------------------------------------------------------------------
548
549 //
550 // Theta Cell
551 //
552 G4String nameCrystalAndCasing = "CrystalAndCasing";
553
554 G4int id = 0; // ID of crystals after distinguishing left and right
555 for ( i = startID; i <= thetaNbCrystals; i++ ) //================
556 {
557 ostringstream strSolidCasing;
558 strSolidCasing << "solidBSCCasing" << i - 1;
559 ostringstream strVolumeCasing;
560 strVolumeCasing << "logicalBSCCasing" << i - 1;
561 ostringstream strPhysiCasing;
562 strPhysiCasing << "physicalBSCCasing" << i - 1;
563
564 if ( i > ( *besEMCGeometry ).BSCNbTheta )
565 {
566 id = i - ( *besEMCGeometry ).BSCNbTheta - 1;
567 solidBSCTheta =
568 new G4Trap( strSolidCasing.str(), zHalfLength[id], thetaAxis[id], -phiAxis[id],
569 yHalfLength1[id], xHalfLength2[id], xHalfLength1[id], -tanAlpha1[id],
570 yHalfLength2[id], xHalfLength4[id], xHalfLength3[id], -tanAlpha2[id] );
571
572 // G4cout<<"in EmcConstr1: "<<strSolidCasing.str()<<" x1="<<xHalfLength1[id]<<"
573 // y1="<<yHalfLength1[id]<<" theta="<<thetaAxis[id]
574 //<<" phi="<<-phiAxis[id]<<" a1="<<-tanAlpha1[id]<<G4endl;
575
576 logicBSCTheta =
577 new G4LogicalVolume( solidBSCTheta, fCasingMaterial, strVolumeCasing.str() );
578
579 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1] = new G4RotationMatrix();
580 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1]->rotateZ( -90 * deg );
581 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1]->rotateX( -thetaPosition[id] );
582
583 physiBSCTheta = new G4PVPlacement(
584 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1],
585 G4ThreeVector( xPosition[id], yPosition[id], zPosition[id] ), strPhysiCasing.str(),
586 logicBSCTheta, physiBSCPhi, false, i - 1, CHECKLV3 );
587
588 if ( logicBSCTheta )
589 {
590 G4VisAttributes* rightVisAtt = new G4VisAttributes( G4Colour( 1.0, 0., 0. ) );
591 rightVisAtt->SetVisibility( true );
592 logicBSCTheta->SetVisAttributes( rightVisAtt );
593 logicBSCTheta->SetVisAttributes( G4VisAttributes::Invisible );
594 }
595
596 ostringstream strRear;
597 strRear << "physicalRearBox_1_" << i - 1;
598
599 physiRear =
600 new G4PVPlacement( rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1],
601 G4ThreeVector( ( *besEMCGeometry ).rearBoxPosX[id],
602 ( *besEMCGeometry ).rearBoxPosY[id],
603 ( *besEMCGeometry ).rearBoxPosZ[id] ),
604 strRear.str(), logicRear, physiBSCPhi, false, i - 1, CHECKLV3 );
605
606 ostringstream strGirder;
607 strGirder << "solidOpenningCutGirder_1_" << i - 1;
608 solidOCGirder =
609 new G4Cons( strGirder.str(), ( *besEMCGeometry ).OCGirderRmin1[id],
610 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).OCGirderRmin2[id],
611 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).OCGirderDz[id] / 2,
612 360. * deg - ( *besEMCGeometry ).OCGirderAngle / 2,
613 ( *besEMCGeometry ).OCGirderAngle / 2 - da );
614
615 ostringstream strVGirder;
616 strVGirder << "logicalOpenningCutGirder_1_" << i - 1;
617 logicOCGirder = new G4LogicalVolume( solidOCGirder, stainlessSteel, strVGirder.str() );
618 logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
619
620 ostringstream strPGirder;
621 strPGirder << "physicalOpenningCutGirder_1_" << i - 1;
622 physiOCGirder = new G4PVPlacement(
623 0, G4ThreeVector( 0, 0, ( *besEMCGeometry ).OCGirderPosZ[id] ), logicOCGirder,
624 strPGirder.str(), logicBSCPhi, false, 0, CHECKLV3 );
625
626 if ( id < ( *besEMCGeometry ).BSCNbTheta - 1 )
627 {
628 G4double zLength = ( *besEMCGeometry ).OCGirderPosZ[id + 1] -
629 ( *besEMCGeometry ).OCGirderPosZ[id] -
630 ( *besEMCGeometry ).OCGirderDz[id + 1] / 2 -
631 ( *besEMCGeometry ).OCGirderDz[id] / 2;
632 G4double zPositionOCGirder = ( *besEMCGeometry ).OCGirderPosZ[id + 1] -
633 ( *besEMCGeometry ).OCGirderDz[id + 1] / 2 -
634 zLength / 2;
635
636 ostringstream strGirder2;
637 strGirder2 << "solidOpenningCutGirder_2_" << i - 1;
638 solidOCGirder = new G4Cons( strGirder2.str(), ( *besEMCGeometry ).OCGirderRmin2[id],
639 ( *besEMCGeometry ).BSCPhiRmax,
640 ( *besEMCGeometry ).OCGirderRmin1[id + 1],
641 ( *besEMCGeometry ).BSCPhiRmax, zLength / 2,
642 360. * deg - ( *besEMCGeometry ).OCGirderAngle / 2,
643 ( *besEMCGeometry ).OCGirderAngle / 2 - da );
644
645 ostringstream strVGirder2;
646 strVGirder2 << "logicalOpenningCutGirder_2_" << i - 1;
647 logicOCGirder =
648 new G4LogicalVolume( solidOCGirder, stainlessSteel, strVGirder2.str() );
649 logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
650
651 ostringstream strPGirder2;
652 strPGirder2 << "physicalOpenningCutGirder_2_" << i - 1;
653 physiOCGirder =
654 new G4PVPlacement( 0, G4ThreeVector( 0, 0, zPositionOCGirder ), logicOCGirder,
655 strPGirder2.str(), logicBSCPhi, false, 0, CHECKLV3 );
656 }
657
658 ostringstream strBSCCable;
659 strBSCCable << "solidBSCCable_1_" << i - 1;
660 solidCable =
661 new G4Tubs( strBSCCable.str(), 0, ( *besEMCGeometry ).cableDr,
662 ( *besEMCGeometry ).cableLength[id] / 2, 0. * deg, 360. * deg );
663
664 ostringstream strVBSCCable;
665 strVBSCCable << "logicalBSCCable_1_" << i - 1;
666 logicCable = new G4LogicalVolume( solidCable, cable, strVBSCCable.str() );
667
668 ostringstream strPBSCCable;
669 strPBSCCable << "physicalBSCCable_1_" << i - 1;
670 physiCable = new G4PVPlacement( 0,
671 G4ThreeVector( ( *besEMCGeometry ).cablePosX[id],
672 ( *besEMCGeometry ).cablePosY[id],
673 ( *besEMCGeometry ).cablePosZ[id] ),
674 logicCable, strPBSCCable.str(), logicBSCPhi, false, 0,
675 CHECKLV3 );
676 logicCable->SetVisAttributes( G4VisAttributes::Invisible );
677 }
678 else
679 {
680 id = ( *besEMCGeometry ).BSCNbTheta - i;
681 solidBSCTheta =
682 new G4Trap( strSolidCasing.str(), zHalfLength[id], thetaAxis[id], phiAxis[id],
683 yHalfLength1[id], xHalfLength1[id], xHalfLength2[id], tanAlpha1[id],
684 yHalfLength2[id], xHalfLength3[id], xHalfLength4[id], tanAlpha2[id] );
685
686 // G4cout<<"in EmcConstr2: "<<strSolidCasing.str()<<"
687 // x1="<<xHalfLength1[id]<<" y1="<<yHalfLength1[id]<<" theta="<<thetaAxis[id]
688 // <<" phi="<<phiAxis[id]<<" a1="<<tanAlpha1[id]<<G4endl;
689
690 logicBSCTheta =
691 new G4LogicalVolume( solidBSCTheta, fCasingMaterial, strVolumeCasing.str() );
692
693 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1] = new G4RotationMatrix();
694 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1]->rotateZ( -90 * deg );
695 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1]->rotateX( -180 * deg +
696 thetaPosition[id] );
697 physiBSCTheta = new G4PVPlacement(
698 rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1],
699 G4ThreeVector( xPosition[id], yPosition[id], -zPosition[id] ),
700 strPhysiCasing.str(), logicBSCTheta, physiBSCPhi, false, i - 1, CHECKLV3 );
701 if ( logicBSCTheta )
702 {
703 G4VisAttributes* rightVisAtt = new G4VisAttributes( G4Colour( 1.0, 0., 0. ) );
704 rightVisAtt->SetVisibility( true );
705 logicBSCTheta->SetVisAttributes( rightVisAtt );
706 logicBSCTheta->SetVisAttributes( G4VisAttributes::Invisible );
707 }
708
709 ostringstream strRear;
710 strRear << "physicalRearBox_2_" << i - 1;
711
712 physiRear =
713 new G4PVPlacement( rotateMatrix[( *besEMCGeometry ).BSCNbPhi + i - 1],
714 G4ThreeVector( ( *besEMCGeometry ).rearBoxPosX[id],
715 ( *besEMCGeometry ).rearBoxPosY[id],
716 -( *besEMCGeometry ).rearBoxPosZ[id] ),
717 strRear.str(), logicRear, physiBSCPhi, false, i - 1, CHECKLV3 );
718
719 ostringstream strGirder;
720 strGirder << "solidOpenningCutGirder_3_" << i - 1;
721 solidOCGirder =
722 new G4Cons( strGirder.str(), ( *besEMCGeometry ).OCGirderRmin2[id],
723 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).OCGirderRmin1[id],
724 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).OCGirderDz[id] / 2,
725 360. * deg - ( *besEMCGeometry ).OCGirderAngle / 2,
726 ( *besEMCGeometry ).OCGirderAngle / 2 - da );
727
728 ostringstream strVGirder;
729 strVGirder << "logicalOpenningCutGirder_3_" << i - 1;
730 logicOCGirder = new G4LogicalVolume( solidOCGirder, stainlessSteel, strVGirder.str() );
731 logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
732
733 ostringstream strPGirder;
734 strPGirder << "physicalOpenningCutGirder_3_" << i - 1;
735 physiOCGirder = new G4PVPlacement(
736 0, G4ThreeVector( 0, 0, -( *besEMCGeometry ).OCGirderPosZ[id] ), logicOCGirder,
737 strPGirder.str(), logicBSCPhi, false, 0, CHECKLV3 );
738
739 if ( id < ( *besEMCGeometry ).BSCNbTheta - 1 )
740 {
741 G4double zLength = ( *besEMCGeometry ).OCGirderPosZ[id + 1] -
742 ( *besEMCGeometry ).OCGirderPosZ[id] -
743 ( *besEMCGeometry ).OCGirderDz[id + 1] / 2 -
744 ( *besEMCGeometry ).OCGirderDz[id] / 2;
745 G4double zPositionOCGirder = ( *besEMCGeometry ).OCGirderPosZ[id + 1] -
746 ( *besEMCGeometry ).OCGirderDz[id + 1] / 2 -
747 zLength / 2;
748
749 ostringstream strGirder2;
750 strGirder2 << "solidOpenningCutGirder_4_" << i - 1;
751 solidOCGirder = new G4Cons(
752 strGirder2.str(), ( *besEMCGeometry ).OCGirderRmin1[id + 1],
753 ( *besEMCGeometry ).BSCPhiRmax, ( *besEMCGeometry ).OCGirderRmin2[id],
754 ( *besEMCGeometry ).BSCPhiRmax, zLength / 2,
755 360. * deg - ( *besEMCGeometry ).OCGirderAngle / 2,
756 ( *besEMCGeometry ).OCGirderAngle / 2 - da );
757
758 ostringstream strVGirder2;
759 strVGirder2 << "logicalOpenningCutGirder_4_" << i - 1;
760 logicOCGirder =
761 new G4LogicalVolume( solidOCGirder, stainlessSteel, strVGirder2.str() );
762 logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
763
764 ostringstream strPGirder2;
765 strPGirder2 << "physicalOpenningCutGirder_4_" << i - 1;
766 physiOCGirder =
767 new G4PVPlacement( 0, G4ThreeVector( 0, 0, -zPositionOCGirder ), logicOCGirder,
768 strPGirder2.str(), logicBSCPhi, false, 0, CHECKLV3 );
769 }
770
771 ostringstream strBSCCable;
772 strBSCCable << "solidBSCCable_2_" << i - 1;
773 solidCable =
774 new G4Tubs( strBSCCable.str(), 0, ( *besEMCGeometry ).cableDr,
775 ( *besEMCGeometry ).cableLength[id] / 2, 0. * deg, 360. * deg );
776
777 ostringstream strVBSCCable;
778 strVBSCCable << "logicalBSCCable_2_" << i - 1;
779 logicCable = new G4LogicalVolume( solidCable, cable, strVBSCCable.str() );
780
781 ostringstream strPBSCCable;
782 strPBSCCable << "physicalBSCCable_2_" << i - 1;
783 physiCable = new G4PVPlacement( 0,
784 G4ThreeVector( ( *besEMCGeometry ).cablePosX[id],
785 ( *besEMCGeometry ).cablePosY[id],
786 -( *besEMCGeometry ).cablePosZ[id] ),
787 logicCable, strPBSCCable.str(), logicBSCPhi, false, 0,
788 CHECKLV3 );
789 logicCable->SetVisAttributes( G4VisAttributes::Invisible );
790 }
791
792 ostringstream strCrystal;
793 strCrystal << "physicalCrystal" << i - 1;
794 physiBSCCrystal =
795 new G4PVParameterised( strCrystal.str(), logicBSCCrystal, physiBSCTheta, kZAxis,
796 1, // for this method,it must be 1.
797 crystalParam );
798 ( *besEMCGeometry ).physiBSCCrystal[i] = physiBSCCrystal;
799 // G4cout << (*besEMCGeometry).physiBSCCrystal[i] << G4endl;
800 physiBSCCrystal->SetCopyNo( i );
801 if ( CHECKLV4 ) physiBSCCrystal->CheckOverlaps();
802
803 if ( verboseLevel > 4 )
804 G4cout
805 << "BesEmcConstruction*****************************" << G4endl
806 << "point of crystal =" << physiBSCCrystal
807 << G4endl
808 // << "point of mother =" <<physiBSCCrystal->GetMotherPhysical() << G4endl
809 << "point of excepted=" << physiBSCTheta << G4endl;
810 // G4Exception("BesEMCConstruction::Construct() starting............");
811 }
812 //
813 // always return the physical World
814 //
815 if ( verboseLevel > 0 ) PrintEMCParameters();
816 // return physiBSC;
817
818 ConstructSPFrame( logicBSCWorld, besEMCGeometry );
819 ConstructEndGeometry( logicEMC );
820 }
821
822 // Set vis attributes and sensitive detector
823 SetVisAndSD();
824
825 // list geo tree
826 if ( logicEMC && physiEMC && verboseLevel > 4 )
827 {
828 G4cout << "logicEmc " << logicEMC << " physiEmc " << physiEMC << G4endl;
829 G4cout << "list geo tree" << G4endl;
830
831 int NdaughterofEMC = logicEMC->GetNoDaughters();
832
833 for ( int i = 0; i < NdaughterofEMC; i++ )
834 {
835 G4LogicalVolume* daughterofEmc = logicEMC->GetDaughter( i )->GetLogicalVolume();
836 G4cout << i << "/" << NdaughterofEMC << " name: " << daughterofEmc->GetName() << " "
837 << daughterofEmc << " shape: " << daughterofEmc->GetSolid()->GetName() << G4endl;
838 int NdaughterofEmc_2 = daughterofEmc->GetNoDaughters();
839 for ( int j = 0; j < NdaughterofEmc_2; j++ )
840 {
841 G4LogicalVolume* daughterofEmc_2 = daughterofEmc->GetDaughter( j )->GetLogicalVolume();
842 G4cout << " --> " << j << "/" << NdaughterofEmc_2
843 << " name: " << daughterofEmc_2->GetName() << " " << daughterofEmc_2
844 << " shape: " << daughterofEmc_2->GetSolid()->GetName() << G4endl;
845 int NdaughterofEmc_3 = daughterofEmc_2->GetNoDaughters();
846 for ( int k = 0; k < NdaughterofEmc_3; k++ )
847 {
848 G4LogicalVolume* daughterofEmc_3 =
849 daughterofEmc_2->GetDaughter( k )->GetLogicalVolume();
850 G4cout << " --> " << k << "/" << NdaughterofEmc_3
851 << " name: " << daughterofEmc_3->GetName() << " " << daughterofEmc_3
852 << " shape: " << daughterofEmc_3->GetSolid()->GetName() << G4endl;
853 int NdaughterofEmc_4 = daughterofEmc_3->GetNoDaughters();
854 for ( int m = 0; m < NdaughterofEmc_4; m++ )
855 {
856 G4LogicalVolume* daughterofEmc_4 =
857 daughterofEmc_3->GetDaughter( m )->GetLogicalVolume();
858 G4cout << " --> " << m << "/" << NdaughterofEmc_4
859 << " name: " << daughterofEmc_4->GetName() << " " << daughterofEmc_4
860 << " shape: " << daughterofEmc_4->GetSolid()->GetName() << G4endl;
861 if ( daughterofEmc_3->GetSolid()->GetName().contains( "solidBSCCasing" ) )
862 {
863 G4Trap* Crystal = (G4Trap*)daughterofEmc_3->GetSolid();
864 double hz = Crystal->GetZHalfLength();
865 double hx1 = Crystal->GetXHalfLength1();
866 double hx2 = Crystal->GetXHalfLength2();
867 double hx3 = Crystal->GetXHalfLength3();
868 double hx4 = Crystal->GetXHalfLength4();
869 double hy1 = Crystal->GetYHalfLength1();
870 double hy2 = Crystal->GetYHalfLength2();
871 double tanalpha1 = Crystal->GetTanAlpha1();
872 double tanalpha2 = Crystal->GetTanAlpha2();
873 G4cout << " --> " << hx1 << " " << hx2 << " " << hx3 << " "
874 << hx4 << " " << hy1 << " " << hy2 << " " << hz << " " << tanalpha1 << " "
875 << tanalpha2 << G4endl;
876
877 } // if(SolidCrystal)
878 } // 4
879 } // 3
880 } // 2
881 } // 1
882 }
883}
884
885void BesEmcConstruction::ConstructEndGeometry( G4LogicalVolume* logEMC ) { // logicEMC is not
886 // necessary, since
887 // it is a member of
888 // this class
889 if ( logEMC != logicEMC )
890 G4cout << "BesEmcConstruction::ConstructEndGeometry parameter transmit error!" << G4endl;
891 // has gotten in DefineMaterials()
892 // G4Material* fCrystalMaterial = G4Material::GetMaterial("Cesiumiodide");
893 G4VisAttributes* crystalVisAtt = new G4VisAttributes( G4Colour( 0.5, 0, 1.0 ) );
894 crystalVisAtt->SetVisibility( false );
895 G4VisAttributes* endPhiVisAtt = new G4VisAttributes( G4Colour( 0, 1.0, 0 ) );
896 endPhiVisAtt->SetVisibility( false );
897 const G4double zoomConst = 0.995;
898 const G4double da = 0.001 * deg;
899
900 // world volume of endcap
901 // east end
902 solidEnd = new G4Cons( "solidEndWorld", ( *emcEnd ).WorldRmin1, ( *emcEnd ).WorldRmax1,
903 ( *emcEnd ).WorldRmin2, ( *emcEnd ).WorldRmax2,
904 ( *emcEnd ).WorldDz / 2, 0. * deg, 360. * deg );
905 logicEnd = new G4LogicalVolume( solidEnd, G4Material::GetMaterial( "Aluminium" ),
906 "logicalEndWorld", 0, 0, 0 );
907 physiEnd = new G4PVPlacement( 0, // no rotation
908 G4ThreeVector( 0, 0, ( *emcEnd ).WorldZPosition ),
909 logicEnd, // its logical volume
910 "physicalEndWorld0", // its name
911 logicEMC, // its mother volume
912 false, // no boolean operations
913 0, CHECKLV1 ); // no field specific to volume
914 if ( logicEnd ) logicEnd->SetVisAttributes( G4VisAttributes::Invisible );
915
916 // west end
917 G4RotationMatrix* rotateEnd = new G4RotationMatrix();
918 rotateEnd->rotateY( 180. * deg );
919 physiEnd = new G4PVPlacement( rotateEnd, G4ThreeVector( 0, 0, -( *emcEnd ).WorldZPosition ),
920 logicEnd, "physicalEndWorld2", logicEMC, false, 2, CHECKLV1 );
921
922 ////////////////////////////////////////////////////////////////////////
923 // emc endcap sectors (east) //
924 //////////////////////////////////////////////////////////////////////////
925 // 20mm gap //
926 // || //
927 // \ 7 || 6 / //
928 // - 8 \ || / 5 - //
929 // - \ || / - //
930 // _ 9 - \ || / - 4 _ //
931 // - _ - \ || / - _ - //
932 // - _ - \||/ - _ - //
933 // 10 - -||- - 3 //
934 // ----------------||---------------- //
935 // 11 - -||- - 2 //
936 // _ - - /||\ - - _ //
937 // _ - - / || \ - - _ //
938 // - 12 - / || \ - 1 - //
939 // - / || \ - //
940 // - 13 / || \ 0 - //
941 // / 14 || 15 \ //
942 // || //
943 ////////////////////////////////////////////////////////////////////////
944
945 // 1/16 of endcap world,which has some symmetry
946 // sector 0-6,8-14
947 solidEndPhi = new G4Cons( "solidEndPhi0", ( *emcEnd ).SectorRmin1, ( *emcEnd ).SectorRmax1,
948 ( *emcEnd ).SectorRmin2, ( *emcEnd ).SectorRmax2,
949 ( *emcEnd ).SectorDz / 2, 0. * deg, 22.5 * deg - da );
950 logicEndPhi = new G4LogicalVolume( solidEndPhi, G4Material::GetMaterial( "Air" ),
951 "logicalEndPhi0", 0, 0, 0 );
952 for ( G4int i = 0; i < 14; i++ )
953 {
954 if ( ( i != 6 ) && ( i != 7 ) )
955 {
956 G4RotationMatrix* rotatePhi = new G4RotationMatrix();
957 rotatePhi->rotateZ( -i * 22.5 * deg + 67.5 * deg );
958 ostringstream strEndPhi;
959 strEndPhi << "physicalEndPhi" << i;
960 physiEndPhi =
961 new G4PVPlacement( rotatePhi, // 0,logicEndPhi,strEndPhi.str(),logicEnd,false,i);
962 G4ThreeVector( 0, 0, ( *emcEnd ).SectorZPosition ), logicEndPhi,
963 strEndPhi.str(), logicEnd, false, i, CHECKLV2 );
964 }
965 }
966 if ( logicEndPhi ) logicEndPhi->SetVisAttributes( endPhiVisAtt );
967
968 for ( G4int i = 0; i < 35; i++ )
969 {
970 ostringstream strEndCasing;
971 strEndCasing << "solidEndCasing_0_" << i;
972
973 //-************tranform to new coodinate! liangyt 2007.5.7 *******
974 G4ThreeVector newfPnt[8];
975 G4ThreeVector center( 0.0, 0.0, 0.0 );
976 G4ThreeVector rotAngle( 0.0, 0.0, 0.0 );
977
978 TransformToArb8( ( *emcEnd ).fPnt[i], newfPnt, center, rotAngle );
979
980 emcEnd->Zoom( newfPnt, zoomConst ); // change emcEnd.fPnt[i] to newfPnt
981
982 G4RotationMatrix* rotatePhiIrregBox = new G4RotationMatrix();
983 rotatePhiIrregBox->rotateX( rotAngle.x() );
984 rotatePhiIrregBox->rotateY( rotAngle.y() );
985 rotatePhiIrregBox->rotateZ( rotAngle.z() );
986 //-*******************************************************************
987
988 solidEndCasing = new G4IrregBox( strEndCasing.str(), ( *emcEnd ).zoomPoint ); // liangyt
989
990 ostringstream strVEndCasing;
991 strVEndCasing << "logicalEndCasing_0_" << i;
992 logicEndCasing =
993 new G4LogicalVolume( solidEndCasing, fCasingMaterial, strVEndCasing.str() );
994
995 ostringstream strPEndCasing;
996 strPEndCasing << "physicalEndCasing_0_" << i;
997 physiEndCasing =
998 new G4PVPlacement( rotatePhiIrregBox, center, logicEndCasing, strPEndCasing.str(),
999 logicEndPhi, false, i, CHECKIrreg ); // change with rot and pos now!
1000
1001 ostringstream strEndCrystal;
1002 strEndCrystal << "solidEndCrystal_0_" << i;
1003
1004 emcEnd->ModifyForCasing( ( *emcEnd ).zoomPoint, i );
1005 solidEndCrystal = new G4IrregBox( strEndCrystal.str(), ( *emcEnd ).cryPoint );
1006
1007 ostringstream strVEndCrystal;
1008 strVEndCrystal << "logicalEndCrystal_0_" << i;
1009 logicEndCrystal =
1010 new G4LogicalVolume( solidEndCrystal, fCrystalMaterial, strVEndCrystal.str() );
1011
1012 ostringstream strPEndCrystal;
1013 strPEndCrystal << "physicalEndCrystal_0_" << i;
1014 physiEndCrystal =
1015 new G4PVPlacement( 0, G4ThreeVector(), logicEndCrystal, strPEndCrystal.str(),
1016 logicEndCasing, false, i, CHECKIrreg );
1017
1018 logicEndCasing->SetVisAttributes( G4VisAttributes::Invisible );
1019 logicEndCrystal->SetVisAttributes( crystalVisAtt );
1020 logicEndCrystal->SetSensitiveDetector( besEMCSD );
1021 }
1022
1023 // the top area which has 20 mm gap
1024 // sector 6,14
1025 G4double rmin2 =
1026 ( *emcEnd ).WorldRmin1 + ( ( *emcEnd ).WorldRmin2 - ( *emcEnd ).WorldRmin1 ) *
1027 ( *emcEnd ).SectorDz / ( *emcEnd ).WorldDz;
1028 G4double rmax2 =
1029 ( *emcEnd ).WorldRmax1 + ( ( *emcEnd ).WorldRmax2 - ( *emcEnd ).WorldRmax1 ) *
1030 ( *emcEnd ).SectorDz / ( *emcEnd ).WorldDz;
1031 solidEndPhi =
1032 new G4Cons( "solidEndPhi1", ( *emcEnd ).WorldRmin1, ( *emcEnd ).WorldRmax1, rmin2, rmax2,
1033 ( *emcEnd ).SectorDz / 2, 67.5 * deg, 22.5 * deg - da );
1034 logicEndPhi = new G4LogicalVolume( solidEndPhi, G4Material::GetMaterial( "Air" ),
1035 "logicalEndPhi1", 0, 0, 0 );
1036 for ( G4int i = 0; i < 2; i++ )
1037 {
1038 G4RotationMatrix* rotatePhi = new G4RotationMatrix();
1039 rotatePhi->rotateZ( -i * 180. * deg );
1040 ostringstream strEndPhi;
1041 strEndPhi << "physicalEndPhi" << i * 8 + 6;
1042 physiEndPhi = new G4PVPlacement(
1043 rotatePhi, G4ThreeVector( 0, 0, ( *emcEnd ).SectorZPosition ), logicEndPhi,
1044 strEndPhi.str(), logicEnd, false, i * 8 + 6, CHECKLV2 );
1045 }
1046 if ( logicEndPhi ) logicEndPhi->SetVisAttributes( endPhiVisAtt );
1047
1048 for ( G4int i = 0; i < 35; i++ )
1049 {
1050 ostringstream strEndCasing;
1051 strEndCasing << "solidEndCasing_1_" << i;
1052
1053 //-************tranform to new coodinate! liangyt 2007.5.7 *******
1054 G4ThreeVector newfPnt[8];
1055 G4ThreeVector center( 0.0, 0.0, 0.0 );
1056 G4ThreeVector rotAngle( 0.0, 0.0, 0.0 );
1057
1058 TransformToArb8( ( *emcEnd ).fPnt1[i], newfPnt, center, rotAngle );
1059
1060 emcEnd->Zoom( newfPnt, zoomConst ); // change emcEnd.fPnt[i] to newfPnt
1061
1062 G4RotationMatrix* rotatePhiIrregBox = new G4RotationMatrix();
1063 rotatePhiIrregBox->rotateX( rotAngle.x() );
1064 rotatePhiIrregBox->rotateY( rotAngle.y() );
1065 rotatePhiIrregBox->rotateZ( rotAngle.z() );
1066 //-*******************************************************************
1067
1068 solidEndCasing = new G4IrregBox( strEndCasing.str(), ( *emcEnd ).zoomPoint );
1069
1070 ostringstream strVEndCasing;
1071 strVEndCasing << "logicalEndCasing_1_" << i;
1072 logicEndCasing =
1073 new G4LogicalVolume( solidEndCasing, fCasingMaterial, strVEndCasing.str() );
1074
1075 ostringstream strPEndCasing;
1076 strPEndCasing << "physicalEndCasing_1_" << i;
1077 physiEndCasing =
1078 new G4PVPlacement( rotatePhiIrregBox, center, logicEndCasing, strPEndCasing.str(),
1079 logicEndPhi, false, i, CHECKIrreg ); // change with rot and pos now!
1080
1081 ostringstream strEndCrystal;
1082 strEndCrystal << "solidEndCrystal_1_" << i;
1083
1084 emcEnd->ModifyForCasing( ( *emcEnd ).zoomPoint, i );
1085 solidEndCrystal = new G4IrregBox( strEndCrystal.str(), ( *emcEnd ).cryPoint );
1086
1087 ostringstream strVEndCrystal;
1088 strVEndCrystal << "logicalEndCrystal_1_" << i;
1089 logicEndCrystal =
1090 new G4LogicalVolume( solidEndCrystal, fCrystalMaterial, strVEndCrystal.str() );
1091
1092 ostringstream strPEndCrystal;
1093 strPEndCrystal << "physicalEndCrystal_1_" << i;
1094 physiEndCrystal =
1095 new G4PVPlacement( 0, G4ThreeVector(), logicEndCrystal, strPEndCrystal.str(),
1096 logicEndCasing, false, i, CHECKIrreg );
1097
1098 logicEndCasing->SetVisAttributes( G4VisAttributes::Invisible );
1099 logicEndCrystal->SetVisAttributes( crystalVisAtt );
1100 logicEndCrystal->SetSensitiveDetector( besEMCSD );
1101 }
1102
1103 ( *emcEnd ).ReflectX();
1104
1105 // sector 7,15
1106 for ( G4int i = 0; i < 35; i++ )
1107 for ( G4int j = 0; j < 8; j++ ) ( *emcEnd ).fPnt1[i][j].rotateZ( -90. * deg );
1108
1109 solidEndPhi = new G4Cons( "solidEndPhi2", ( *emcEnd ).WorldRmin1, ( *emcEnd ).WorldRmax1,
1110 rmin2, rmax2, ( *emcEnd ).SectorDz / 2, 0 * deg, 22.5 * deg - da );
1111 logicEndPhi = new G4LogicalVolume( solidEndPhi, G4Material::GetMaterial( "Air" ),
1112 "logicalEndPhi2", 0, 0, 0 );
1113 for ( G4int i = 0; i < 2; i++ )
1114 {
1115 G4RotationMatrix* rotatePhi = new G4RotationMatrix();
1116 rotatePhi->rotateZ( -i * 180. * deg - 90. * deg );
1117 ostringstream strEndPhi;
1118 strEndPhi << "physicalEndPhi" << i * 8 + 7;
1119 physiEndPhi = new G4PVPlacement(
1120 rotatePhi, G4ThreeVector( 0, 0, ( *emcEnd ).SectorZPosition ), logicEndPhi,
1121 strEndPhi.str(), logicEnd, false, i * 8 + 7, CHECKLV2 );
1122 }
1123 if ( logicEndPhi ) logicEndPhi->SetVisAttributes( endPhiVisAtt );
1124
1125 for ( G4int i = 0; i < 35; i++ )
1126 {
1127 ostringstream strEndCasing;
1128 strEndCasing << "solidEndCasing_2_" << i;
1129
1130 //-************tranform to new coodinate! liangyt 2007.5.7 *******
1131 G4ThreeVector newfPnt[8];
1132 G4ThreeVector center( 0.0, 0.0, 0.0 );
1133 G4ThreeVector rotAngle( 0.0, 0.0, 0.0 );
1134
1135 TransformToArb8( ( *emcEnd ).fPnt1[i], newfPnt, center, rotAngle );
1136
1137 emcEnd->Zoom( newfPnt, zoomConst ); // change emcEnd.fPnt[i] to newfPnt
1138
1139 G4RotationMatrix* rotatePhiIrregBox = new G4RotationMatrix();
1140 rotatePhiIrregBox->rotateX( rotAngle.x() );
1141 rotatePhiIrregBox->rotateY( rotAngle.y() );
1142 rotatePhiIrregBox->rotateZ( rotAngle.z() );
1143 //-*******************************************************************
1144
1145 solidEndCasing = new G4IrregBox( strEndCasing.str(), ( *emcEnd ).zoomPoint );
1146
1147 ostringstream strVEndCasing;
1148 strVEndCasing << "logicalEndCasing_2_" << i;
1149 logicEndCasing =
1150 new G4LogicalVolume( solidEndCasing, fCasingMaterial, strVEndCasing.str() );
1151
1152 ostringstream strPEndCasing;
1153 strPEndCasing << "physicalEndCasing_2_" << i;
1154 physiEndCasing =
1155 new G4PVPlacement( rotatePhiIrregBox, center, logicEndCasing, strPEndCasing.str(),
1156 logicEndPhi, false, i, CHECKIrreg ); // change with rot and pos now!
1157
1158 ostringstream strEndCrystal;
1159 strEndCrystal << "solidEndCrystal_2_" << i;
1160
1161 emcEnd->ModifyForCasing( ( *emcEnd ).zoomPoint, i );
1162 solidEndCrystal = new G4IrregBox( strEndCrystal.str(), ( *emcEnd ).cryPoint );
1163
1164 ostringstream strVEndCrystal;
1165 strVEndCrystal << "logicalEndCrystal_2_" << i;
1166 logicEndCrystal =
1167 new G4LogicalVolume( solidEndCrystal, fCrystalMaterial, strVEndCrystal.str() );
1168
1169 ostringstream strPEndCrystal;
1170 strPEndCrystal << "physicalEndCrystal_2_" << i;
1171 physiEndCrystal =
1172 new G4PVPlacement( 0, G4ThreeVector(), logicEndCrystal, strPEndCrystal.str(),
1173 logicEndCasing, false, i, CHECKIrreg );
1174
1175 logicEndCasing->SetVisAttributes( G4VisAttributes::Invisible );
1176 logicEndCrystal->SetVisAttributes( crystalVisAtt );
1177 logicEndCrystal->SetSensitiveDetector( besEMCSD );
1178 }
1179}
1180
1182 G4int copyNb;
1183 switch ( num )
1184 {
1185 case 30: copyNb = 5; break;
1186 case 31: copyNb = 6; break;
1187 case 32: copyNb = 14; break;
1188 case 33: copyNb = 15; break;
1189 case 34: copyNb = 16; break;
1190 default: copyNb = num; break;
1191 }
1192 return copyNb;
1193}
1194
1195void BesEmcConstruction::ConstructSPFrame( G4LogicalVolume* logMother,
1196 BesEmcGeometry* emcGeometry ) {
1197 G4double rmax = ( *besEMCGeometry ).BSCRmax + 2. * mm; // radius from 942mm to 940mm
1198 solidSupportBar = new G4Tubs(
1199 "solidSupportBar0", rmax + ( *besEMCGeometry ).SPBarThickness1,
1200 rmax + ( *besEMCGeometry ).SPBarThickness + ( *besEMCGeometry ).SPBarThickness1,
1201 ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 +
1202 ( *besEMCGeometry ).EndRingDz,
1203 0. * deg, 360. * deg );
1204
1205 logicSupportBar =
1206 new G4LogicalVolume( solidSupportBar, stainlessSteel, "logicalSupportBar0" );
1207
1208 physiSupportBar = new G4PVPlacement( 0, G4ThreeVector(), logicSupportBar,
1209 "physicalSupportBar0", logMother, false, 0, CHECKLV2 );
1210
1211 solidSupportBar1 =
1212 new G4Tubs( "solidSupportBar1", rmax, rmax + ( *besEMCGeometry ).SPBarThickness1,
1213 ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3,
1214 ( *besEMCGeometry ).BSCPhiDphi - ( *besEMCGeometry ).SPBarDphi / 2,
1215 ( *besEMCGeometry ).SPBarDphi );
1216
1217 logicSupportBar1 =
1218 new G4LogicalVolume( solidSupportBar1, stainlessSteel, "logicalSupportBar1" );
1219
1220 for ( G4int i = 0; i < ( *besEMCGeometry ).BSCNbPhi / 2; i++ )
1221 {
1222 G4RotationMatrix* rotateSPBar = new G4RotationMatrix();
1223 rotateSPBar->rotateZ( ( *besEMCGeometry ).BSCPhiDphi -
1224 i * 2 * ( *besEMCGeometry ).BSCPhiDphi );
1225 ostringstream strSupportBar1;
1226 strSupportBar1 << "physicalSupportBar1_" << i;
1227 physiSupportBar1 =
1228 new G4PVPlacement( rotateSPBar, G4ThreeVector(), logicSupportBar1,
1229 strSupportBar1.str(), logMother, false, 0, CHECKLV2 );
1230 }
1231
1232 // end ring
1233 solidEndRing =
1234 new G4Tubs( "solidEndRing", ( *besEMCGeometry ).EndRingRmin,
1235 ( *besEMCGeometry ).EndRingRmin + ( *besEMCGeometry ).EndRingDr / 2,
1236 ( *besEMCGeometry ).EndRingDz / 2, 0. * deg, 360. * deg );
1237
1238 solidGear = new G4Tubs(
1239 "solidGear", ( *besEMCGeometry ).EndRingRmin + ( *besEMCGeometry ).EndRingDr / 2,
1240 ( *besEMCGeometry ).EndRingRmin + ( *besEMCGeometry ).EndRingDr,
1241 ( *besEMCGeometry ).EndRingDz / 2, 0. * deg, ( *besEMCGeometry ).BSCPhiDphi );
1242
1243 // taper ring
1244 solidTaperRing1 =
1245 new G4Tubs( "solidTaperRing1", ( *besEMCGeometry ).TaperRingRmin1,
1246 ( *besEMCGeometry ).TaperRingRmin1 + ( *besEMCGeometry ).TaperRingThickness1,
1247 ( *besEMCGeometry ).TaperRingInnerLength / 2, 0. * deg, 360. * deg );
1248
1249 solidTaperRing2 =
1250 new G4Cons( "solidTaperRing2", ( *besEMCGeometry ).TaperRingRmin1,
1251 ( *besEMCGeometry ).TaperRingRmin1 + ( *besEMCGeometry ).TaperRingThickness1,
1252 ( *besEMCGeometry ).TaperRingRmin2,
1253 ( *besEMCGeometry ).TaperRingRmin2 + ( *besEMCGeometry ).TaperRingDr,
1254 ( *besEMCGeometry ).TaperRingDz / 2, 0. * deg, 360. * deg );
1255
1256 solidTaperRing3 =
1257 new G4Cons( "solidTaperRing3", ( *besEMCGeometry ).BSCRmax2,
1258 ( *besEMCGeometry ).BSCRmax2 + ( *besEMCGeometry ).TaperRingOuterLength1,
1259 ( *besEMCGeometry ).TaperRingRmin2 + ( *besEMCGeometry ).TaperRingDr,
1260 ( *besEMCGeometry ).TaperRingRmin2 + ( *besEMCGeometry ).TaperRingDr +
1261 ( *besEMCGeometry ).TaperRingOuterLength,
1262 ( *besEMCGeometry ).TaperRingThickness3 / 2, 0. * deg, 360. * deg );
1263
1264 logicEndRing = new G4LogicalVolume( solidEndRing, stainlessSteel, "logicalEndRing" );
1265 logicGear = new G4LogicalVolume( solidGear, stainlessSteel, "logicalGear" );
1266 logicTaperRing1 =
1267 new G4LogicalVolume( solidTaperRing1, stainlessSteel, "logicalTaperRing1" );
1268 logicTaperRing2 =
1269 new G4LogicalVolume( solidTaperRing2, stainlessSteel, "logicalTaperRing2" );
1270 logicTaperRing3 =
1271 new G4LogicalVolume( solidTaperRing3, stainlessSteel, "logicalTaperRing3" );
1272
1273 for ( G4int i = 0; i < 2; i++ )
1274 {
1275 G4RotationMatrix* rotateSPRing = new G4RotationMatrix();
1276 G4double zEndRing, z1, z2, z3;
1277 if ( i == 0 )
1278 {
1279 zEndRing = ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 +
1280 ( *besEMCGeometry ).EndRingDz / 2;
1281 z1 = ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 -
1282 ( *besEMCGeometry ).TaperRingDz - ( *besEMCGeometry ).TaperRingInnerLength / 2;
1283 z2 = ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 -
1284 ( *besEMCGeometry ).TaperRingDz / 2;
1285 z3 = ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 / 2;
1286 }
1287 else
1288 {
1289 rotateSPRing->rotateY( 180. * deg );
1290 zEndRing = -( ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 +
1291 ( *besEMCGeometry ).EndRingDz / 2 );
1292 z1 = -( ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 -
1293 ( *besEMCGeometry ).TaperRingDz - ( *besEMCGeometry ).TaperRingInnerLength / 2 );
1294 z2 = -( ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 -
1295 ( *besEMCGeometry ).TaperRingDz / 2 );
1296 z3 = -( ( *besEMCGeometry ).BSCDz + ( *besEMCGeometry ).TaperRingThickness3 / 2 );
1297 }
1298
1299 ostringstream strEndRing;
1300 strEndRing << "physicalEndRing_" << i;
1301 physiEndRing =
1302 new G4PVPlacement( rotateSPRing, G4ThreeVector( 0, 0, zEndRing ), logicEndRing,
1303 strEndRing.str(), logMother, false, 0, CHECKLV2 );
1304
1305 for ( G4int j = 0; j < ( *besEMCGeometry ).BSCNbPhi / 2; j++ )
1306 {
1307 G4RotationMatrix* rotateGear = new G4RotationMatrix();
1308 rotateGear->rotateZ( ( *besEMCGeometry ).BSCPhiDphi / 2 -
1309 j * 2 * ( *besEMCGeometry ).BSCPhiDphi );
1310
1311 ostringstream strGear;
1312 strGear << "physicalGear_" << i << "_" << j;
1313 physiGear = new G4PVPlacement( rotateGear, G4ThreeVector( 0, 0, zEndRing ), logicGear,
1314 strGear.str(), logMother, false, 0, CHECKLV2 );
1315 }
1316
1317 ostringstream strTaperRing1;
1318 strTaperRing1 << "physicalTaperRing1_" << i;
1319 physiTaperRing1 =
1320 new G4PVPlacement( rotateSPRing, G4ThreeVector( 0, 0, z1 ), logicTaperRing1,
1321 strTaperRing1.str(), logMother, false, 0, CHECKLV2 );
1322
1323 ostringstream strTaperRing2;
1324 strTaperRing2 << "physicalTaperRing2_" << i;
1325 physiTaperRing2 =
1326 new G4PVPlacement( rotateSPRing, G4ThreeVector( 0, 0, z2 ), logicTaperRing2,
1327 strTaperRing2.str(), logMother, false, 0, CHECKLV2 );
1328
1329 ostringstream strTaperRing3;
1330 strTaperRing3 << "physicalTaperRing3_" << i;
1331 physiTaperRing3 =
1332 new G4PVPlacement( rotateSPRing, G4ThreeVector( 0, 0, z3 ), logicTaperRing3,
1333 strTaperRing3.str(), logMother, false, 0, CHECKLV2 );
1334 }
1335}
1336
1337// Set vis attributes and sensitive detector
1339 //-------------------------------------------------------------
1340 // Barrel
1341 G4VisAttributes* bscVisAtt = new G4VisAttributes( G4Colour( 0.5, 0.5, 0.5 ) );
1342 bscVisAtt->SetVisibility( false );
1343 logicEMC->SetVisAttributes( bscVisAtt );
1344 if ( logicBSCWorld ) logicBSCWorld->SetVisAttributes( G4VisAttributes::Invisible );
1345
1346 if ( logicBSCCrystal )
1347 {
1348 // G4cout<<"find BSCCrystal "<<G4endl;
1349 G4VisAttributes* crystalVisAtt = new G4VisAttributes( G4Colour( 0, 0, 1.0 ) );
1350 crystalVisAtt->SetVisibility( true );
1351 logicBSCCrystal->SetVisAttributes( crystalVisAtt );
1352 logicBSCCrystal->SetSensitiveDetector( besEMCSD );
1353 }
1354
1355 if ( logicBSCPhi )
1356 {
1357 G4VisAttributes* rightVisAtt = new G4VisAttributes( G4Colour( 1.0, 0., 1.0 ) );
1358 rightVisAtt->SetVisibility( false );
1359 logicBSCPhi->SetVisAttributes( rightVisAtt );
1360 }
1361
1362 if ( logicRear ) logicRear->SetVisAttributes( G4VisAttributes::Invisible );
1363 if ( logicOrgGlass ) logicOrgGlass->SetVisAttributes( G4VisAttributes::Invisible );
1364 if ( logicRearCasing ) logicRearCasing->SetVisAttributes( G4VisAttributes::Invisible );
1365 if ( logicAlPlate ) logicAlPlate->SetVisAttributes( G4VisAttributes::Invisible );
1366 if ( logicPD )
1367 {
1368 logicPD->SetVisAttributes( G4VisAttributes::Invisible );
1369 logicPD->SetSensitiveDetector( besEMCSD );
1370 }
1371 if ( logicPreAmpBox ) logicPreAmpBox->SetVisAttributes( G4VisAttributes::Invisible );
1372 if ( logicAirInPABox ) logicAirInPABox->SetVisAttributes( G4VisAttributes::Invisible );
1373 if ( logicHangingPlate ) logicHangingPlate->SetVisAttributes( G4VisAttributes::Invisible );
1374 if ( logicWaterPipe ) logicWaterPipe->SetVisAttributes( G4VisAttributes::Invisible );
1375
1376 //-------------------------------------------------------------
1377 // Support system
1378 G4VisAttributes* ringVisAtt = new G4VisAttributes( G4Colour( 0.5, 0.25, 0. ) );
1379 ringVisAtt->SetVisibility( false );
1380 if ( logicSupportBar ) logicSupportBar->SetVisAttributes( ringVisAtt );
1381 if ( logicSupportBar1 ) logicSupportBar1->SetVisAttributes( ringVisAtt );
1382 if ( logicEndRing ) logicEndRing->SetVisAttributes( ringVisAtt );
1383 if ( logicGear ) logicGear->SetVisAttributes( ringVisAtt );
1384 if ( logicTaperRing1 ) logicTaperRing1->SetVisAttributes( ringVisAtt );
1385 if ( logicTaperRing2 ) logicTaperRing2->SetVisAttributes( ringVisAtt );
1386 if ( logicTaperRing3 ) logicTaperRing3->SetVisAttributes( ringVisAtt );
1387
1388 //-------------------------------------------------------------
1389 // Endcap
1390 G4VisAttributes* endPhiVisAtt = new G4VisAttributes( G4Colour( 0, 1.0, 0 ) );
1391 endPhiVisAtt->SetVisibility( false );
1392 if ( logicEnd ) logicEnd->SetVisAttributes( endPhiVisAtt );
1393}
1394
1395// Get logical volumes from gdml
1397 //-------------------------------------------------------------
1398 // Barrel
1399 logicBSCWorld = FindLogicalVolume( "logicalBSCWorld" );
1400 logicBSCCrystal = FindLogicalVolume( "logicalCrystal" );
1401 logicBSCPhi = FindLogicalVolume( "logicalBSCPhi" );
1402 logicRear = FindLogicalVolume( "logicalRearBox" );
1403 logicOrgGlass = FindLogicalVolume( "logicalOrganicGlass" );
1404 logicRearCasing = FindLogicalVolume( "logicalRearCasing" );
1405 logicAlPlate = FindLogicalVolume( "logicalAlPlate" );
1406 logicPD = FindLogicalVolume( "logicalPD" );
1407 logicPreAmpBox = FindLogicalVolume( "logicalPreAmpBox" );
1408 logicAirInPABox = FindLogicalVolume( "logicalAirInPABox" );
1409 logicHangingPlate = FindLogicalVolume( "logicalHangingPlate" );
1410 logicWaterPipe = FindLogicalVolume( "logicalWaterPipe" );
1411
1412 for ( int i = 0; i < 44; i++ )
1413 {
1414 std::ostringstream osnameBSCCasing;
1415 osnameBSCCasing << "logicalBSCCasing" << i;
1416 logicBSCTheta = FindLogicalVolume( osnameBSCCasing.str() );
1417 if ( logicBSCTheta )
1418 {
1419 G4VisAttributes* rightVisAtt = new G4VisAttributes( G4Colour( 1.0, 0.0, 0.0 ) );
1420 rightVisAtt->SetVisibility( false );
1421 logicBSCTheta->SetVisAttributes( rightVisAtt );
1422 }
1423
1424 std::ostringstream osnameBSCCable1;
1425 osnameBSCCable1 << "logicalBSCCable_1_" << i;
1426 logicCable = FindLogicalVolume( osnameBSCCable1.str() );
1427 if ( logicCable ) logicCable->SetVisAttributes( G4VisAttributes::Invisible );
1428
1429 std::ostringstream osnameBSCCable2;
1430 osnameBSCCable2 << "logicalBSCCable_2_" << i;
1431 logicCable = FindLogicalVolume( osnameBSCCable2.str() );
1432 if ( logicCable ) logicCable->SetVisAttributes( G4VisAttributes::Invisible );
1433
1434 std::ostringstream osnameOCGirder1;
1435 osnameOCGirder1 << "logicalOpenningCutGirder_1_" << i;
1436 logicOCGirder = FindLogicalVolume( osnameOCGirder1.str() );
1437 if ( logicOCGirder ) logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
1438
1439 std::ostringstream osnameOCGirder2;
1440 osnameOCGirder2 << "logicalOpenningCutGirder_2_" << i;
1441 logicOCGirder = FindLogicalVolume( osnameOCGirder2.str() );
1442 if ( logicOCGirder ) logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
1443
1444 std::ostringstream osnameOCGirder3;
1445 osnameOCGirder3 << "logicalOpenningCutGirder_3_" << i;
1446 logicOCGirder = FindLogicalVolume( osnameOCGirder3.str() );
1447 if ( logicOCGirder ) logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
1448
1449 std::ostringstream osnameOCGirder4;
1450 osnameOCGirder4 << "logicalOpenningCutGirder_4_" << i;
1451 logicOCGirder = FindLogicalVolume( osnameOCGirder4.str() );
1452 if ( logicOCGirder ) logicOCGirder->SetVisAttributes( G4VisAttributes::Invisible );
1453 }
1454
1455 //-------------------------------------------------------------
1456 // Support system
1457 logicSupportBar = FindLogicalVolume( "logicalSupportBar0" );
1458 logicSupportBar1 = FindLogicalVolume( "logicalSupportBar1" );
1459 logicEndRing = FindLogicalVolume( "logicalEndRing" );
1460 logicGear = FindLogicalVolume( "logicalGear" );
1461 logicTaperRing1 = FindLogicalVolume( "logicalTaperRing1" );
1462 logicTaperRing2 = FindLogicalVolume( "logicalTaperRing2" );
1463 logicTaperRing3 = FindLogicalVolume( "logicalTaperRing3" );
1464
1465 //-------------------------------------------------------------
1466 // Endcap
1467 logicEnd = FindLogicalVolume( "logicalEndWorld" );
1468
1469 for ( G4int sector = 0; sector < 3; sector++ )
1470 {
1471 std::ostringstream osnameEndPhi;
1472 osnameEndPhi << "logicalEndPhi" << sector;
1473 logicEndPhi = FindLogicalVolume( osnameEndPhi.str() );
1474 if ( logicEndPhi ) { logicEndPhi->SetVisAttributes( G4VisAttributes::Invisible ); }
1475 else { G4cout << "Can't find logicEndPhi!" << G4endl; }
1476
1477 for ( G4int cryNb = 0; cryNb < 35; cryNb++ )
1478 {
1479
1480 std::ostringstream osnameEndCrystal;
1481 osnameEndCrystal << "logicalEndCrystal_" << sector << "_" << cryNb;
1482 logicEndCrystal = FindLogicalVolume( osnameEndCrystal.str() );
1483 if ( logicEndCrystal )
1484 {
1485 logicEndCrystal->SetSensitiveDetector( besEMCSD );
1486 G4VisAttributes* crystalVisAtt = new G4VisAttributes( G4Colour( 0.5, 0, 1.0 ) );
1487 crystalVisAtt->SetVisibility( false );
1488 logicEndCrystal->SetVisAttributes( crystalVisAtt );
1489 }
1490 else { G4cout << "Can't find: " << osnameEndCrystal.str() << G4endl; }
1491
1492 std::ostringstream osnameEndCasing;
1493 osnameEndCasing << "logicalEndCasing_" << sector << "_" << cryNb;
1494 logicEndCasing = FindLogicalVolume( osnameEndCasing.str() );
1495 if ( logicEndCasing ) { logicEndCasing->SetVisAttributes( G4VisAttributes::Invisible ); }
1496 else { G4cout << "Can't find: " << osnameEndCasing.str() << G4endl; }
1497 }
1498 }
1499}
1500
1501void BesEmcConstruction::DefineMaterials() {
1502 G4String name, symbol; // a=mass of a mole;
1503 G4double a, z, density; // z=mean number of protons;
1504 // G4int iz, n; //iz=number of protons in an isotope;
1505 // n=number of nucleons in an isotope;
1506
1507 G4int ncomponents, natoms;
1508 G4double fractionmass;
1509 // G4double abundance, fractionmass;
1510 // G4double temperature, pressure;
1511
1512 // for debug
1513 // G4Exception("BesEmcConstruction::DefineMaterials() starting...........");
1514 //
1515 // define Elements
1516 //
1517 G4Element* H = G4Element::GetElement( "Hydrogen" );
1518 if ( !H )
1519 {
1520 a = 1.01 * g / mole;
1521 H = new G4Element( name = "Hydrogen", symbol = "H", z = 1., a );
1522 }
1523 G4Element* C = G4Element::GetElement( "Carbon" );
1524 if ( !C )
1525 {
1526 a = 12.01 * g / mole;
1527 C = new G4Element( name = "Carbon", symbol = "C", z = 6., a );
1528 }
1529 G4Element* O = G4Element::GetElement( "Oxygen" );
1530 if ( !O )
1531 {
1532 a = 16.00 * g / mole;
1533 O = new G4Element( name = "Oxygen", symbol = "O", z = 8., a );
1534 }
1535
1536 density = 0.344 * g / cm3;
1537 G4Material* Tyvek = new G4Material( name = "M_Polyethylene", density, ncomponents = 2 );
1538 Tyvek->AddElement( C, natoms = 1 );
1539 Tyvek->AddElement( H, natoms = 2 );
1540
1541 density = 1.39 * g / cm3;
1542 G4Material* Mylar =
1543 new G4Material( name = "M_PolyethyleneTerephthlate", density, ncomponents = 3 );
1544 Mylar->AddElement( C, natoms = 5 );
1545 Mylar->AddElement( H, natoms = 4 );
1546 Mylar->AddElement( O, natoms = 2 );
1547
1548 density = 1.18 * g / cm3;
1549 organicGlass = new G4Material( name = "M_OrganicGlass", density, ncomponents = 3 );
1550 organicGlass->AddElement( C, natoms = 5 );
1551 organicGlass->AddElement( H, natoms = 7 );
1552 organicGlass->AddElement( O, natoms = 2 );
1553
1554 G4Material* Fe = new G4Material( name = "M_Iron", z = 26., a = 55.85 * g / mole,
1555 density = 7.87 * g / cm3 );
1556 G4Material* Cr = new G4Material( name = "M_Chromium", z = 24., a = 52.00 * g / mole,
1557 density = 8.72 * g / cm3 );
1558 G4Material* Ni = new G4Material( name = "M_Nickel", z = 28., a = 58.69 * g / mole,
1559 density = 8.72 * g / cm3 );
1560
1561 stainlessSteel =
1562 new G4Material( name = "M_Cr18Ni9", density = 7.85 * g / cm3, ncomponents = 3 );
1563 stainlessSteel->AddMaterial( Fe, fractionmass = 73. * perCent );
1564 stainlessSteel->AddMaterial( Cr, fractionmass = 18. * perCent );
1565 stainlessSteel->AddMaterial( Ni, fractionmass = 9. * perCent );
1566
1567 G4Material* H2O = G4Material::GetMaterial( "Water" );
1568 G4Material* Cu = G4Material::GetMaterial( "Copper" );
1569 G4double dWater = 1. * g / cm3; // density
1570 G4double dCopper = 8.96 * g / cm3;
1571 G4double aWater =
1572 ( ( *besEMCGeometry ).waterPipeDr - ( *besEMCGeometry ).waterPipeThickness ) *
1573 ( ( *besEMCGeometry ).waterPipeDr - ( *besEMCGeometry ).waterPipeThickness ); // area
1574 G4double aCopper =
1575 ( *besEMCGeometry ).waterPipeDr * ( *besEMCGeometry ).waterPipeDr - aWater;
1576 density = ( dWater * aWater + dCopper * aCopper ) / ( aWater + aCopper );
1577
1578 waterPipe = new G4Material( name = "M_WaterPipe", density, ncomponents = 2 );
1579 fractionmass = dWater * aWater / ( dWater * aWater + dCopper * aCopper );
1580 waterPipe->AddMaterial( H2O, fractionmass );
1581 fractionmass = dCopper * aCopper / ( dWater * aWater + dCopper * aCopper );
1582 waterPipe->AddMaterial( Cu, fractionmass );
1583
1584 cable = new G4Material( name = "M_Cable", density = 4. * g / cm3, ncomponents = 1 );
1585 cable->AddMaterial( Cu, 1 );
1586
1587 // for debug
1588 // G4Exception("BesEmcConstruction::DefineMaterials() running one.........");
1589 //
1590 // predigest the casing of crystals to a mixture
1591 //
1592
1593 G4Material* Al = G4Material::GetMaterial( "Aluminium" );
1594 if ( Al == NULL )
1595 {
1596 Al = new G4Material( name = "Aluminium", z = 13., a = 26.98 * g / mole,
1597 density = 2.700 * g / cm3 );
1598 }
1599
1600 G4Material* Si = G4Material::GetMaterial( "M_Silicon" );
1601 if ( Si == NULL )
1602 {
1603 Si = new G4Material( name = "M_Silicon", z = 14., a = 28.0855 * g / mole,
1604 density = 2.33 * g / cm3 );
1605 }
1606
1607 // for debug
1608 G4double totalThickness = ( *besEMCGeometry ).fTyvekThickness +
1609 ( *besEMCGeometry ).fAlThickness +
1610 ( *besEMCGeometry ).fMylarThickness;
1611 density = ( Tyvek->GetDensity() * ( *besEMCGeometry ).fTyvekThickness +
1612 Al->GetDensity() * ( *besEMCGeometry ).fAlThickness +
1613 Mylar->GetDensity() * ( *besEMCGeometry ).fMylarThickness ) /
1614 totalThickness;
1615 G4Material* Casing = new G4Material( name = "M_Casing", density, ncomponents = 3 );
1616 Casing->AddMaterial( Tyvek, fractionmass = Tyvek->GetDensity() / density *
1617 ( *besEMCGeometry ).fTyvekThickness /
1618 totalThickness );
1619 Casing->AddMaterial( Al, fractionmass = Al->GetDensity() / density *
1620 ( *besEMCGeometry ).fAlThickness / totalThickness );
1621 Casing->AddMaterial( Mylar, fractionmass = Mylar->GetDensity() / density *
1622 ( *besEMCGeometry ).fMylarThickness /
1623 totalThickness );
1624 fCasingMaterial = Casing;
1625 rearCasingMaterial = Tyvek;
1626 // for debug
1627 // G4Exception("BesEmcConstruction::DefineMaterials() running two........");
1628 fCrystalMaterial = G4Material::GetMaterial( "Cesiumiodide" );
1629}
1630
1631//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1633 G4cout << "-------------------------------------------------------" << G4endl
1634 << "---> There are " << phiNbCrystals << "(max=" << ( *besEMCGeometry ).BSCNbPhi
1635 << ") crystals along phi direction and " << thetaNbCrystals
1636 << "(max=" << ( *besEMCGeometry ).BSCNbTheta << ") crystals along theta direction."
1637 << G4endl << "The crystals have sizes of " << ( *besEMCGeometry ).BSCCryLength / cm
1638 << "cm(L) and " << ( *besEMCGeometry ).BSCYFront / cm << "cm(Y) with "
1639 << fCrystalMaterial->GetName() << "." << G4endl << "The casing is layer of "
1640 << ( *besEMCGeometry ).fTyvekThickness / mm << "mm tyvek,"
1641 << ( *besEMCGeometry ).fAlThickness / mm << "mm aluminum and"
1642 << ( *besEMCGeometry ).fMylarThickness / mm << "mm mylar." << G4endl
1643 << "-------------------------------------------------------" << G4endl;
1644 G4cout << G4Material::GetMaterial( "PolyethyleneTerephthlate" ) << G4endl
1645 << G4Material::GetMaterial( "Casing" ) << G4endl
1646 << G4Material::GetMaterial( "Polyethylene" ) << G4endl
1647 << "-------------------------------------------------------" << G4endl;
1648}
1649
1650//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1651
1652void BesEmcConstruction::SetCrystalMaterial( G4String materialChoice ) {
1653 // search the material by its name
1654 G4Material* pttoMaterial = G4Material::GetMaterial( materialChoice );
1655 if ( pttoMaterial )
1656 {
1657 fCrystalMaterial = pttoMaterial;
1658 logicBSCCrystal->SetMaterial( pttoMaterial );
1660 }
1661}
1662
1663//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1664
1665void BesEmcConstruction::SetCasingMaterial( G4String materialChoice ) {
1666 // search the material by its name
1667 G4Material* pttoMaterial = G4Material::GetMaterial( materialChoice );
1668 if ( pttoMaterial )
1669 {
1670 fCasingMaterial = pttoMaterial;
1671 logicBSCTheta->SetMaterial( pttoMaterial );
1673 }
1674}
1675
1676//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1677
1679 // change Gap thickness and recompute the calorimeter parameters
1680 ( *besEMCGeometry ).fTyvekThickness = val( 'X' );
1681 ( *besEMCGeometry ).fAlThickness = val( 'Y' );
1682 ( *besEMCGeometry ).fMylarThickness = val( 'Z' );
1683}
1684
1685//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1686
1687void BesEmcConstruction::SetBSCRmin( G4double val ) { ( *besEMCGeometry ).BSCRmin = val; }
1688
1689//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1690
1691void BesEmcConstruction::SetBSCNbPhi( G4int val ) { ( *besEMCGeometry ).BSCNbPhi = val; }
1692
1693//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1694
1695void BesEmcConstruction::SetBSCNbTheta( G4int val ) { ( *besEMCGeometry ).BSCNbTheta = val; }
1696
1697void BesEmcConstruction::SetStartIDTheta( G4int val ) { startID = val; }
1698
1699//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1700
1702 ( *besEMCGeometry ).BSCCryLength = val;
1703}
1704
1705//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1706
1708 ( *besEMCGeometry ).BSCYFront0 = val;
1709}
1710
1711//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1712
1713void BesEmcConstruction::SetBSCYFront( G4double val ) { ( *besEMCGeometry ).BSCYFront = val; }
1714
1715//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1716
1718 ( *besEMCGeometry ).BSCPosition0 = val;
1719}
1720
1721//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1722
1724 ( *besEMCGeometry ).BSCPosition1 = val;
1725}
1726
1727//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1728
1729void BesEmcConstruction::SetMagField( G4double fieldValue ) {
1730 // apply a global uniform magnetic field along Z axis
1731 G4FieldManager* fieldMgr =
1732 G4TransportationManager::GetTransportationManager()->GetFieldManager();
1733
1734 if ( magField ) delete magField; // delete the existing magn field
1735
1736 if ( fieldValue != 0. ) // create a new one if non nul
1737 {
1738 magField = new G4UniformMagField( G4ThreeVector( 0., 0., fieldValue ) );
1739 fieldMgr->SetDetectorField( magField );
1740 fieldMgr->CreateChordFinder( magField );
1741 fmagField = fieldValue;
1742 }
1743 else
1744 {
1745 magField = 0;
1746 fieldMgr->SetDetectorField( magField );
1747 fmagField = 0.;
1748 }
1749}
1750
1751//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1752//???????????????????????????????????????????????
1754 ; // G4RunManager::GetRunManager()->DefineWorldVolume(BesDetectorConstruction::Construct());
1755}
1756
1757//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1758
1759void BesEmcConstruction::TransformToArb8( const G4ThreeVector fPnt[8],
1760 G4ThreeVector newfPnt[8], G4ThreeVector& center,
1761 G4ThreeVector& rotAngle ) {
1762 HepPoint3D point[8];
1763 center = G4ThreeVector( 0.0, 0.0, 0.0 );
1764 for ( int i = 0; i < 8; i++ )
1765 {
1766 point[i] = HepPoint3D( fPnt[i].x(), fPnt[i].y(), fPnt[i].z() );
1767 center += point[i];
1768 }
1769 center /= 8.0;
1770
1771 HepPlane3D bottomPlane( point[4], point[5], point[6] );
1772 HepPoint3D centerProject = bottomPlane.point( center );
1773 Hep3Vector newZ = center - centerProject;
1774
1775 rotAngle = RotAngleFromNewZ( newZ );
1776 G4RotationMatrix* g4Rot = new G4RotationMatrix();
1777 g4Rot->rotateX( rotAngle.x() );
1778 g4Rot->rotateY( rotAngle.y() );
1779
1780 G4AffineTransform* transform = new G4AffineTransform( g4Rot, center );
1781 transform->Invert();
1782 for ( int i = 0; i < 8; i++ ) { newfPnt[i] = transform->TransformPoint( fPnt[i] ); }
1783 delete g4Rot;
1784 delete transform;
1785}
1786
1787void BesEmcConstruction::ThreeVectorTrans( G4ThreeVector fPnt[8], double x[8], double y[8],
1788 double z[8] ) {
1789 for ( int i = 0; i < 8; i++ )
1790 {
1791 x[i] = fPnt[i].x();
1792 y[i] = fPnt[i].y();
1793 z[i] = fPnt[i].z();
1794 }
1795}
1796
1797Hep3Vector BesEmcConstruction::RotAngleFromNewZ( Hep3Vector newZ ) {
1798 newZ.setMag( 1.0 );
1799 Hep3Vector x0( 1, 0, 0 ), y0( 0, 1, 0 ), z0( 0, 0, 1 );
1800 double dx, dy, dz = 0.0;
1801
1802 Hep3Vector a( 0.0, newZ.y(), newZ.z() );
1803 // three rotated angles, rotate dx angles(rad) by X axis,
1804 // then rotate dy by new Y axis, then dz by new Z axis;
1805 // all rotate in left hand system;
1806
1807 // a, project of final newZ on original YZ plane, 0 <= theta < pi;
1808 if ( a.mag() != 0.0 ) a.setMag( 1.0 );
1809 else cout << "newZ on X axis, a=(0,0,0)" << endl;
1810 dx = acos( a.dot( z0 ) );
1811 if ( a.dot( z0 ) == -1.0 ) dx = 0.0;
1812
1813 // b, rotate dx angle of z0 around original X axis(x0),
1814 Hep3Vector b( 0, sin( dx ), cos( dx ) );
1815 dy = acos( b.dot( newZ ) );
1816 if ( newZ.x() > 0.0 ) dy = -dy;
1817
1818 Hep3Vector v( dx, dy, dz );
1819 return v;
1820}
#define CHECKIrreg
#define CHECKLV0
#define CHECKLV4
#define CHECKLV2
#define CHECKLV3
#define CHECKLV1
#define CHECKLV5
double P(RecMdcKalTrack *trk)
Double_t x[10]
**********Class see also m_nmax DOUBLE PRECISION m_amel DOUBLE PRECISION m_x2 DOUBLE PRECISION m_alfinv DOUBLE PRECISION m_Xenph INTEGER m_KeyWtm INTEGER m_idyfs DOUBLE PRECISION m_zini DOUBLE PRECISION m_q2 DOUBLE PRECISION m_Wt_KF DOUBLE PRECISION m_WtCut INTEGER m_KFfin *COMMON c_KarLud $ !Input CMS energy[GeV] $ !CMS energy after beam spread beam strahlung[GeV] $ !Beam energy spread[GeV] $ !z boost due to beam spread $ !electron beam mass *ff pair spectrum $ !minimum v
Definition KarLud.h:35
***************************************************************************************Pseudo Class RRes *****************************************************************************************Parameters and physical constants **Maarten sept ************************************************************************DOUBLE PRECISION xsmu **************************************************************************PARTICLE DATA all others are from PDG *Only resonances with known widths into electron pairs are sept ************************************************************************C Declarations C
Definition RRes.h:29
static BesEmcConstruction * GetBesEmcConstruction()
void ThreeVectorTrans(G4ThreeVector fPnt[8], double x[8], double y[8], double z[8])
void ConstructSPFrame(G4LogicalVolume *, BesEmcGeometry *)
void Construct(G4LogicalVolume *)
void TransformToArb8(const G4ThreeVector fPnt[8], G4ThreeVector newfPnt[8], G4ThreeVector &center, G4ThreeVector &rotAngle)
void ConstructEndGeometry(G4LogicalVolume *)
void SetBSCPosition1(G4double)
void SetCasingThickness(G4ThreeVector)
void SetBSCCrystalLength(G4double)
void SetCrystalMaterial(G4String)
void SetBSCPosition0(G4double)
Hep3Vector RotAngleFromNewZ(Hep3Vector newZ)
void SetCasingMaterial(G4String)
static void Kill()
static EmcG4Geo * Instance()
Get a pointer to the single instance of EmcG4Geo.
Definition EmcG4Geo.cxx:43
G4LogicalVolume * GetTopVolume()
Get the top(world) volume;.
IMPLICIT REAL *A H
Definition myXsection.h:1