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