12#include "G4ThreeVector.hh"
13#include "G4ReflectionFactory.hh"
15#include "BesEmcConstruction.hh"
16#include "BesEmcDetectorMessenger.hh"
17#include "BesEmcGeometry.hh"
18#include "BesCrystalParameterisation.hh"
19#include "BesEmcEndGeometry.hh"
20#include "ReadBoostRoot.hh"
21#include "BesEmcParameter.hh"
24#include "G4IrregBox.hh"
26#include "G4Transform3D.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"
40#include "G4UniformMagField.hh"
41#include "G4FieldManager.hh"
42#include "G4TransportationManager.hh"
43#include "G4SDManager.hh"
44#include "G4RunManager.hh"
45#include "G4VisAttributes.hh"
47#include "G4AffineTransform.hh"
50#include "G4Geo/EmcG4Geo.h"
55{
return fBesEmcConstruction;}
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)
71 if(fBesEmcConstruction)
72 { G4Exception(
"BesEmcConstruction constructed twice."); }
73 fBesEmcConstruction=
this;
85 if(detectorMessenger)
delete detectorMessenger;
86 if(crystalParam)
delete crystalParam;
87 if(besEMCGeometry)
delete besEMCGeometry;
88 if(emcEnd)
delete emcEnd;
99 G4SDManager* SDman = G4SDManager::GetSDMpointer();
101 besEMCSD =
new BesEmcSD(
"CalorSD",
this,besEMCGeometry);
102 SDman->AddNewDetector( besEMCSD );
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;
126 phiNbCrystals = (*besEMCGeometry).BSCNbPhi;
127 thetaNbCrystals = (*besEMCGeometry).BSCNbTheta*2;
129 G4double da=0.001*deg;
134 solidBSC =
new G4Tubs(
"solidBSC",
135 (*besEMCGeometry).TaperRingRmin1,
136 (*besEMCGeometry).BSCRmax+(*besEMCGeometry).SPBarThickness+(*besEMCGeometry).SPBarThickness1+2.1*mm,
137 (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz,
141 solidESC =
new G4Cons(
"solidESC",(*emcEnd).WorldRmin1,(*emcEnd).WorldRmax1,
142 (*emcEnd).WorldRmin2,(*emcEnd).WorldRmax2,
143 (*emcEnd).WorldDz/2,0.*deg,360.*deg);
145 solidEMC =
new G4UnionSolid(
"solidEMC0",
149 G4ThreeVector(0,0,(*emcEnd).WorldZPosition));
151 G4RotationMatrix *rotateESC =
new G4RotationMatrix();
152 rotateESC->rotateY(180.*deg);
154 solidEMC =
new G4UnionSolid(
"solidEMC",
158 G4ThreeVector(0,0,-(*emcEnd).WorldZPosition));
160 logicEMC =
new G4LogicalVolume(solidEMC,
161 G4Material::GetMaterial(
"Air"),
164 physiEMC =
new G4PVPlacement(0,
172 solidBSCWorld =
new G4SubtractionSolid(
"solidBSCWorld0",
176 G4ThreeVector(0,0,(*emcEnd).WorldZPosition));
178 solidBSCWorld =
new G4SubtractionSolid(
"solidBSCWorld",
182 G4ThreeVector(0,0,-(*emcEnd).WorldZPosition));
184 logicBSCWorld =
new G4LogicalVolume(solidBSCWorld,
185 G4Material::GetMaterial(
"Air"),
188 G4RotationMatrix *rotBSC =
new G4RotationMatrix();
189 rotBSC->rotateY(180.*deg);
190 physiBSCWorld =
new G4PVPlacement(rotBSC,
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);
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;
235 solidBSCPhiTub =
new G4Tubs(
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,
250 solidBSCPhi1 =
new G4SubtractionSolid(
"solidBSCPhi1",
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,
263 solidBSCPhi =
new G4SubtractionSolid(
"solidBSCPhi",
267 G4ThreeVector(0,0,(*besEMCGeometry).BSCDz1/2-(*besEMCGeometry).BSCDz));
269 logicBSCPhi =
new G4LogicalVolume(solidBSCPhi,
270 G4Material::GetMaterial(
"Air"),
274 for(G4int j=0;j<(*besEMCGeometry).BSCNbPhi;j++)
276 if(j<(*besEMCGeometry).BSCNbPhi/2) {
277 i=(*besEMCGeometry).BSCNbPhi/2-j-1;
279 i=(*besEMCGeometry).BSCNbPhi*3/2-j-1;
281 rotateMatrix[i] =
new G4RotationMatrix();
282 rotateMatrix[i]->rotateZ(-i*(*besEMCGeometry).BSCPhiDphi
283 -(*besEMCGeometry).BSCAngleRotat
284 -(*besEMCGeometry).BSCPhiDphi/2.
286 rotateMatrix[i]->getAngleAxis(delta, axis);
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);
296 ostringstream strPhi;
297 strPhi <<
"physicalBSCPhi" << j;
299 physiBSCPhi =
new G4PVPlacement(rotateMatrix[i],
300 G4ThreeVector(ox,oy,oz),
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++)
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];
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;
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);
366 logicBSCCrystal =
new G4LogicalVolume(solidBSCCrystal,
373 (*besEMCGeometry).BSCNbTheta*2,
379 solidRear =
new G4Box(
"solidRearBox",
380 (*besEMCGeometry).rearBoxLength/2,
381 (*besEMCGeometry).rearBoxLength/2,
382 (*besEMCGeometry).rearBoxDz/2);
384 logicRear =
new G4LogicalVolume(solidRear,
385 G4Material::GetMaterial(
"Air"),
389 solidOrgGlass =
new G4Box(
"solidOrganicGlass",
390 (*besEMCGeometry).orgGlassLengthX/2,
391 (*besEMCGeometry).orgGlassLengthY/2,
392 (*besEMCGeometry).orgGlassLengthZ/2);
394 logicOrgGlass =
new G4LogicalVolume(solidOrgGlass,
396 "logicalOrganicGlass");
398 physiOrgGlass =
new G4PVPlacement(0,
399 G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz-(*besEMCGeometry).orgGlassLengthZ)/2),
401 "physicalOrganicGlass",
407 solidCasingBox =
new G4Box(
"solidCasingBox",
408 (*besEMCGeometry).rearBoxLength/2,
409 (*besEMCGeometry).rearBoxLength/2,
410 (*besEMCGeometry).rearCasingThickness/2);
412 solidAirHole =
new G4Box(
"solidAirHole",
413 (*besEMCGeometry).orgGlassLengthX/2,
414 (*besEMCGeometry).orgGlassLengthY/2,
415 (*besEMCGeometry).rearBoxDz/2);
417 solidRearCasing =
new G4SubtractionSolid(
"solidRearCasing",
423 logicRearCasing =
new G4LogicalVolume(solidRearCasing,
425 "logicalRearCasing");
427 physiRearCasing =
new G4PVPlacement(0,
428 G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz-(*besEMCGeometry).rearCasingThickness)/2),
430 "physicalRearCasing",
436 solidAlBox =
new G4Box(
"solidAlBox",
437 (*besEMCGeometry).rearBoxLength/2,
438 (*besEMCGeometry).rearBoxLength/2,
439 (*besEMCGeometry).AlPlateDz/2);
441 solidAlPlate =
new G4SubtractionSolid(
"solidAlPlate",
447 logicAlPlate =
new G4LogicalVolume(solidAlPlate,
448 G4Material::GetMaterial(
"Aluminium"),
451 physiAlPlate =
new G4PVPlacement(0,
452 G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
453 -(*besEMCGeometry).rearCasingThickness
454 -(*besEMCGeometry).AlPlateDz/2)),
462 solidPD =
new G4Box(
"solidPD",
463 (*besEMCGeometry).PDLengthX,
464 (*besEMCGeometry).PDLengthY/2,
465 (*besEMCGeometry).PDLengthZ/2);
467 logicPD =
new G4LogicalVolume(solidPD,
468 G4Material::GetMaterial(
"M_Silicon"),
471 physiPD =
new G4PVPlacement(0,
472 G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
473 -(*besEMCGeometry).orgGlassLengthZ
474 -(*besEMCGeometry).PDLengthZ/2)),
482 solidPreAmpBox =
new G4Box(
"solidPreAmpBox",
483 (*besEMCGeometry).rearBoxLength/2,
484 (*besEMCGeometry).rearBoxLength/2,
485 (*besEMCGeometry).PABoxDz/2);
487 logicPreAmpBox =
new G4LogicalVolume(solidPreAmpBox,
488 G4Material::GetMaterial(
"Aluminium"),
491 physiPreAmpBox =
new G4PVPlacement(0,
492 G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
493 -(*besEMCGeometry).rearCasingThickness
494 -(*besEMCGeometry).AlPlateDz
495 -(*besEMCGeometry).PABoxDz/2)),
503 solidAirInPABox =
new G4Box(
"solidAirInPABox",
504 (*besEMCGeometry).rearBoxLength/2-(*besEMCGeometry).PABoxThickness,
505 (*besEMCGeometry).rearBoxLength/2-(*besEMCGeometry).PABoxThickness,
506 (*besEMCGeometry).PABoxDz/2-(*besEMCGeometry).PABoxThickness);
508 logicAirInPABox =
new G4LogicalVolume(solidAirInPABox,
509 G4Material::GetMaterial(
"Air"),
510 "logicalAirInPABox");
512 physiAirInPABox =
new G4PVPlacement(0,
515 "physicalAirInPABox",
521 solidHangingPlate =
new G4Box(
"solidHangingPlate",
522 (*besEMCGeometry).rearBoxLength/2,
523 (*besEMCGeometry).rearBoxLength/2,
524 (*besEMCGeometry).HangingPlateDz/2);
526 logicHangingPlate =
new G4LogicalVolume(solidHangingPlate,stainlessSteel,
"logicalHangingPlate");
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)),
535 "physicalHangingPlate",
541 solidWaterPipe =
new G4Tubs(
"solidWaterPipe",
543 (*besEMCGeometry).waterPipeDr,
544 (*besEMCGeometry).BSCDz,
548 logicWaterPipe =
new G4LogicalVolume(solidWaterPipe,waterPipe,
"logicalWaterPipe");
550 physiWaterPipe =
new G4PVPlacement(0,
551 G4ThreeVector((*besEMCGeometry).cablePosX[0]-2*(*besEMCGeometry).cableDr,
552 (*besEMCGeometry).cablePosY[0]-(*besEMCGeometry).cableDr-(*besEMCGeometry).waterPipeDr,
565 G4String nameCrystalAndCasing=
"CrystalAndCasing";
568 for(i=startID;i<=thetaNbCrystals;i++)
570 ostringstream strSolidCasing;
571 strSolidCasing <<
"solidBSCCasing" << i-1;
572 ostringstream strVolumeCasing;
573 strVolumeCasing <<
"logicalBSCCasing" << i-1;
574 ostringstream strPhysiCasing;
575 strPhysiCasing <<
"physicalBSCCasing" << i-1;
577 if(i>(*besEMCGeometry).BSCNbTheta)
579 id=i-(*besEMCGeometry).BSCNbTheta-1;
580 solidBSCTheta =
new G4Trap(strSolidCasing.str(),
596 logicBSCTheta =
new G4LogicalVolume(solidBSCTheta,
598 strVolumeCasing.str());
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]);
607 new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
608 G4ThreeVector(xPosition[
id],
611 strPhysiCasing.str(),
619 G4VisAttributes* rightVisAtt=
new G4VisAttributes(G4Colour(1.0,0.,0.));
620 rightVisAtt->SetVisibility(
true);
621 logicBSCTheta->SetVisAttributes(rightVisAtt);
622 logicBSCTheta->SetVisAttributes(G4VisAttributes::Invisible);
625 ostringstream strRear;
626 strRear <<
"physicalRearBox_1_" << i-1;
628 physiRear =
new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
629 G4ThreeVector((*besEMCGeometry).rearBoxPosX[
id],
630 (*besEMCGeometry).rearBoxPosY[
id],
631 (*besEMCGeometry).rearBoxPosZ[
id]),
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);
649 ostringstream strVGirder;
650 strVGirder <<
"logicalOpenningCutGirder_1_" << i-1;
651 logicOCGirder =
new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder.str());
652 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
654 ostringstream strPGirder;
655 strPGirder <<
"physicalOpenningCutGirder_1_" << i-1;
656 physiOCGirder =
new G4PVPlacement(0,
657 G4ThreeVector(0,0,(*besEMCGeometry).OCGirderPosZ[
id]),
664 if(
id<(*besEMCGeometry).BSCNbTheta-1)
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;
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,
680 360.*deg-(*besEMCGeometry).OCGirderAngle/2,
681 (*besEMCGeometry).OCGirderAngle/2-da);
683 ostringstream strVGirder2;
684 strVGirder2 <<
"logicalOpenningCutGirder_2_" << i-1;
685 logicOCGirder =
new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder2.str());
686 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
688 ostringstream strPGirder2;
689 strPGirder2 <<
"physicalOpenningCutGirder_2_" << i-1;
690 physiOCGirder =
new G4PVPlacement(0,
691 G4ThreeVector(0,0,zPosition),
699 ostringstream strBSCCable;
700 strBSCCable <<
"solidBSCCable_1_" << i-1;
701 solidCable =
new G4Tubs(strBSCCable.str(),
703 (*besEMCGeometry).cableDr,
704 (*besEMCGeometry).cableLength[
id]/2,
708 ostringstream strVBSCCable;
709 strVBSCCable <<
"logicalBSCCable_1_" << i-1;
710 logicCable =
new G4LogicalVolume(solidCable,cable,strVBSCCable.str());
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]),
723 logicCable->SetVisAttributes(G4VisAttributes::Invisible);
727 id=(*besEMCGeometry).BSCNbTheta-i;
728 solidBSCTheta =
new G4Trap(strSolidCasing.str(),
744 logicBSCTheta =
new G4LogicalVolume(solidBSCTheta,
746 strVolumeCasing.str());
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]);
753 new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
754 G4ThreeVector(xPosition[
id],
757 strPhysiCasing.str(),
764 G4VisAttributes* rightVisAtt=
new G4VisAttributes(G4Colour(1.0,0.,0.));
765 rightVisAtt->SetVisibility(
true);
766 logicBSCTheta->SetVisAttributes(rightVisAtt);
767 logicBSCTheta->SetVisAttributes(G4VisAttributes::Invisible);
770 ostringstream strRear;
771 strRear <<
"physicalRearBox_2_" << i-1;
773 physiRear =
new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
774 G4ThreeVector((*besEMCGeometry).rearBoxPosX[
id],
775 (*besEMCGeometry).rearBoxPosY[
id],
776 -(*besEMCGeometry).rearBoxPosZ[
id]),
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);
794 ostringstream strVGirder;
795 strVGirder <<
"logicalOpenningCutGirder_3_" << i-1;
796 logicOCGirder =
new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder.str());
797 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
799 ostringstream strPGirder;
800 strPGirder <<
"physicalOpenningCutGirder_3_" << i-1;
801 physiOCGirder =
new G4PVPlacement(0,
802 G4ThreeVector(0,0,-(*besEMCGeometry).OCGirderPosZ[
id]),
809 if(
id<(*besEMCGeometry).BSCNbTheta-1)
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;
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,
823 360.*deg-(*besEMCGeometry).OCGirderAngle/2,
824 (*besEMCGeometry).OCGirderAngle/2-da);
826 ostringstream strVGirder2;
827 strVGirder2 <<
"logicalOpenningCutGirder_4_" << i-1;
829 =
new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder2.str());
830 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
832 ostringstream strPGirder2;
833 strPGirder2 <<
"physicalOpenningCutGirder_4_" << i-1;
834 physiOCGirder =
new G4PVPlacement(0,
835 G4ThreeVector(0,0,-zPosition),
843 ostringstream strBSCCable;
844 strBSCCable <<
"solidBSCCable_2_" << i-1;
845 solidCable =
new G4Tubs(strBSCCable.str(),
847 (*besEMCGeometry).cableDr,
848 (*besEMCGeometry).cableLength[
id]/2,
852 ostringstream strVBSCCable;
853 strVBSCCable <<
"logicalBSCCable_2_" << i-1;
854 logicCable =
new G4LogicalVolume(solidCable,cable,strVBSCCable.str());
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]),
867 logicCable->SetVisAttributes(G4VisAttributes::Invisible);
871 ostringstream strCrystal;
872 strCrystal <<
"physicalCrystal" << i-1;
873 physiBSCCrystal =
new G4PVParameterised(
880 (*besEMCGeometry).physiBSCCrystal[i]=physiBSCCrystal;
882 physiBSCCrystal->SetCopyNo(i);
886 G4cout <<
"BesEmcConstruction*****************************"<< G4endl
887 <<
"point of crystal =" <<physiBSCCrystal << G4endl
889 <<
"point of excepted=" <<physiBSCTheta << G4endl;
907 if(logicEMC&&physiEMC&&verboseLevel>4){
908 G4cout<<
"logicEmc "<<logicEMC<<
" physiEmc "<<physiEMC<<G4endl;
909 G4cout<<
"list geo tree"<<G4endl;
911 int NdaughterofEMC = logicEMC->GetNoDaughters();
913 for(
int i = 0; i < NdaughterofEMC; i++)
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++)
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++)
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++)
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"))
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;
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;
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,
973 G4ThreeVector(0,0,(*emcEnd).WorldZPosition),
980 logicEnd->SetVisAttributes(G4VisAttributes::Invisible);
984 G4RotationMatrix *rotateEnd =
new G4RotationMatrix();
985 rotateEnd->rotateY(180.*deg);
986 physiEnd =
new G4PVPlacement(rotateEnd,
987 G4ThreeVector(0,0,-(*emcEnd).WorldZPosition),
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++)
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,
1032 G4ThreeVector(0,0,(*emcEnd).SectorZPosition),logicEndPhi,strEndPhi.str(),logicEnd,
false,i);
1036 logicEndPhi->SetVisAttributes(endPhiVisAtt);
1038 for(G4int i=0;i<35;i++)
1040 ostringstream strEndCasing;
1041 strEndCasing <<
"solidEndCasing_0_" << i;
1044 G4ThreeVector newfPnt[8];
1045 G4ThreeVector center(0.0, 0.0, 0.0);
1046 G4ThreeVector rotAngle(0.0, 0.0, 0.0);
1050 emcEnd->
Zoom(newfPnt,zoomConst);
1052 G4RotationMatrix *rotatePhiIrregBox =
new G4RotationMatrix();
1053 rotatePhiIrregBox->rotateX(rotAngle.x());
1054 rotatePhiIrregBox->rotateY(rotAngle.y());
1055 rotatePhiIrregBox->rotateZ(rotAngle.z());
1058 solidEndCasing =
new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);
1060 ostringstream strVEndCasing;
1061 strVEndCasing <<
"logicalEndCasing_0_" << i;
1062 logicEndCasing =
new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
1064 ostringstream strPEndCasing;
1065 strPEndCasing <<
"physicalEndCasing_0_" << i;
1066 physiEndCasing =
new G4PVPlacement(rotatePhiIrregBox,center,
1067 logicEndCasing,strPEndCasing.str(),logicEndPhi,
false,i);
1069 ostringstream strEndCrystal;
1070 strEndCrystal <<
"solidEndCrystal_0_" << i;
1073 solidEndCrystal =
new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
1075 ostringstream strVEndCrystal;
1076 strVEndCrystal <<
"logicalEndCrystal_0_" << i;
1077 logicEndCrystal =
new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
1079 ostringstream strPEndCrystal;
1080 strPEndCrystal <<
"physicalEndCrystal_0_" << i;
1081 physiEndCrystal =
new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,
false,i);
1083 logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
1084 logicEndCrystal->SetVisAttributes(crystalVisAtt);
1085 logicEndCrystal->SetSensitiveDetector(besEMCSD);
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++)
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);
1104 logicEndPhi->SetVisAttributes(endPhiVisAtt);
1106 for(G4int i=0;i<35;i++)
1108 ostringstream strEndCasing;
1109 strEndCasing <<
"solidEndCasing_1_" << i;
1112 G4ThreeVector newfPnt[8];
1113 G4ThreeVector center(0.0, 0.0, 0.0);
1114 G4ThreeVector rotAngle(0.0, 0.0, 0.0);
1118 emcEnd->
Zoom(newfPnt,zoomConst);
1120 G4RotationMatrix *rotatePhiIrregBox =
new G4RotationMatrix();
1121 rotatePhiIrregBox->rotateX(rotAngle.x());
1122 rotatePhiIrregBox->rotateY(rotAngle.y());
1123 rotatePhiIrregBox->rotateZ(rotAngle.z());
1126 solidEndCasing =
new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);
1128 ostringstream strVEndCasing;
1129 strVEndCasing <<
"logicalEndCasing_1_" << i;
1130 logicEndCasing =
new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
1132 ostringstream strPEndCasing;
1133 strPEndCasing <<
"physicalEndCasing_1_" << i;
1134 physiEndCasing =
new G4PVPlacement(rotatePhiIrregBox,center,
1135 logicEndCasing,strPEndCasing.str(),logicEndPhi,
false,i);
1137 ostringstream strEndCrystal;
1138 strEndCrystal <<
"solidEndCrystal_1_" << i;
1141 solidEndCrystal =
new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
1143 ostringstream strVEndCrystal;
1144 strVEndCrystal <<
"logicalEndCrystal_1_" << i;
1145 logicEndCrystal =
new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
1147 ostringstream strPEndCrystal;
1148 strPEndCrystal <<
"physicalEndCrystal_1_" << i;
1149 physiEndCrystal =
new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,
false,i);
1151 logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
1152 logicEndCrystal->SetVisAttributes(crystalVisAtt);
1153 logicEndCrystal->SetSensitiveDetector(besEMCSD);
1156 (*emcEnd).ReflectX();
1159 for(G4int i=0;i<35;i++)
1160 for (G4int j=0;j<8;j++)
1161 (*emcEnd).fPnt1[i][j].rotateZ(-90.*deg);
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++)
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);
1177 logicEndPhi->SetVisAttributes(endPhiVisAtt);
1179 for(G4int i=0;i<35;i++)
1181 ostringstream strEndCasing;
1182 strEndCasing <<
"solidEndCasing_2_" << i;
1185 G4ThreeVector newfPnt[8];
1186 G4ThreeVector center(0.0, 0.0, 0.0);
1187 G4ThreeVector rotAngle(0.0, 0.0, 0.0);
1191 emcEnd->
Zoom(newfPnt,zoomConst);
1193 G4RotationMatrix *rotatePhiIrregBox =
new G4RotationMatrix();
1194 rotatePhiIrregBox->rotateX(rotAngle.x());
1195 rotatePhiIrregBox->rotateY(rotAngle.y());
1196 rotatePhiIrregBox->rotateZ(rotAngle.z());
1199 solidEndCasing =
new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);
1201 ostringstream strVEndCasing;
1202 strVEndCasing <<
"logicalEndCasing_2_" << i;
1203 logicEndCasing =
new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
1205 ostringstream strPEndCasing;
1206 strPEndCasing <<
"physicalEndCasing_2_" << i;
1207 physiEndCasing =
new G4PVPlacement(rotatePhiIrregBox,center,
1208 logicEndCasing,strPEndCasing.str(),logicEndPhi,
false,i);
1210 ostringstream strEndCrystal;
1211 strEndCrystal <<
"solidEndCrystal_2_" << i;
1214 solidEndCrystal =
new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
1216 ostringstream strVEndCrystal;
1217 strVEndCrystal <<
"logicalEndCrystal_2_" << i;
1218 logicEndCrystal =
new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
1220 ostringstream strPEndCrystal;
1221 strPEndCrystal <<
"physicalEndCrystal_2_" << i;
1222 physiEndCrystal =
new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,
false,i);
1224 logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
1225 logicEndCrystal->SetVisAttributes(crystalVisAtt);
1226 logicEndCrystal->SetSensitiveDetector(besEMCSD);
1258 G4double rmax=(*besEMCGeometry).BSCRmax+2.*mm;
1259 solidSupportBar =
new G4Tubs(
"solidSupportBar0",
1260 rmax+(*besEMCGeometry).SPBarThickness1,
1261 rmax+(*besEMCGeometry).SPBarThickness+(*besEMCGeometry).SPBarThickness1,
1262 (*besEMCGeometry).BSCDz
1263 +(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz,
1267 logicSupportBar =
new G4LogicalVolume(solidSupportBar,stainlessSteel,
"logicalSupportBar0");
1269 physiSupportBar =
new G4PVPlacement(0,0,logicSupportBar,
"physicalSupportBar0",logicEMC,
false,0);
1271 solidSupportBar1 =
new G4Tubs(
"solidSupportBar1",
1273 rmax+(*besEMCGeometry).SPBarThickness1,
1274 (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3,
1275 (*besEMCGeometry).BSCPhiDphi-(*besEMCGeometry).SPBarDphi/2,
1276 (*besEMCGeometry).SPBarDphi);
1278 logicSupportBar1 =
new G4LogicalVolume(solidSupportBar1,stainlessSteel,
"logicalSupportBar1");
1280 for(G4int i=0;i<(*besEMCGeometry).BSCNbPhi/2;i++)
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);
1291 solidEndRing =
new G4Tubs(
"solidEndRing",
1292 (*besEMCGeometry).EndRingRmin,
1293 (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr/2,
1294 (*besEMCGeometry).EndRingDz/2,
1298 solidGear =
new G4Tubs(
"solidGear",
1299 (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr/2,
1300 (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr,
1301 (*besEMCGeometry).EndRingDz/2,
1303 (*besEMCGeometry).BSCPhiDphi);
1306 solidTaperRing1 =
new G4Tubs(
"solidTaperRing1",
1307 (*besEMCGeometry).TaperRingRmin1,
1308 (*besEMCGeometry).TaperRingRmin1+(*besEMCGeometry).TaperRingThickness1,
1309 (*besEMCGeometry).TaperRingInnerLength/2,
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,
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,
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");
1337 for(G4int i=0;i<2;i++)
1339 G4RotationMatrix *rotateSPRing =
new G4RotationMatrix();
1340 G4double zEndRing,z1,z2,z3;
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;
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);
1359 ostringstream strEndRing;
1360 strEndRing <<
"physicalEndRing_" << i;
1361 physiEndRing =
new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,zEndRing),
1362 logicEndRing,strEndRing.str(),logicEMC,
false,0);
1364 for(G4int j=0;j<(*besEMCGeometry).BSCNbPhi/2;j++)
1366 G4RotationMatrix *rotateGear =
new G4RotationMatrix();
1367 rotateGear->rotateZ((*besEMCGeometry).BSCPhiDphi/2-j*2*(*besEMCGeometry).BSCPhiDphi);
1369 ostringstream strGear;
1370 strGear <<
"physicalGear_" << i <<
"_" <<j;
1371 physiGear =
new G4PVPlacement(rotateGear,G4ThreeVector(0,0,zEndRing),
1372 logicGear,strGear.str(),logicEMC,
false,0);
1375 ostringstream strTaperRing1;
1376 strTaperRing1 <<
"physicalTaperRing1_" << i;
1377 physiTaperRing1 =
new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z1),
1378 logicTaperRing1,strTaperRing1.str(),logicEMC,
false,0);
1380 ostringstream strTaperRing2;
1381 strTaperRing2 <<
"physicalTaperRing2_" << i;
1382 physiTaperRing2 =
new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z2),
1383 logicTaperRing2,strTaperRing2.str(),logicEMC,
false,0);
1385 ostringstream strTaperRing3;
1386 strTaperRing3 <<
"physicalTaperRing3_" << i;
1387 physiTaperRing3 =
new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z3),
1388 logicTaperRing3,strTaperRing3.str(),logicEMC,
false,0);
1397 G4VisAttributes* bscVisAtt=
new G4VisAttributes(G4Colour(0.5,0.5,0.5));
1398 bscVisAtt->SetVisibility(
false);
1399 logicEMC->SetVisAttributes(bscVisAtt);
1401 logicBSCWorld->SetVisAttributes(G4VisAttributes::Invisible);
1403 if (logicBSCCrystal) {
1405 G4VisAttributes* crystalVisAtt=
new G4VisAttributes(G4Colour(0,0,1.0));
1406 crystalVisAtt->SetVisibility(
true);
1407 logicBSCCrystal->SetVisAttributes(crystalVisAtt);
1408 logicBSCCrystal->SetSensitiveDetector(besEMCSD);
1412 G4VisAttributes* rightVisAtt=
new G4VisAttributes(G4Colour(1.0,0.,1.0));
1413 rightVisAtt->SetVisibility(
false);
1414 logicBSCPhi->SetVisAttributes(rightVisAtt);
1418 logicRear->SetVisAttributes(G4VisAttributes::Invisible);
1420 logicOrgGlass->SetVisAttributes(G4VisAttributes::Invisible);
1422 logicRearCasing->SetVisAttributes(G4VisAttributes::Invisible);
1424 logicAlPlate->SetVisAttributes(G4VisAttributes::Invisible);
1426 logicPD->SetVisAttributes(G4VisAttributes::Invisible);
1427 logicPD->SetSensitiveDetector(besEMCSD);
1430 logicPreAmpBox->SetVisAttributes(G4VisAttributes::Invisible);
1432 logicAirInPABox->SetVisAttributes(G4VisAttributes::Invisible);
1433 if(logicHangingPlate)
1434 logicHangingPlate->SetVisAttributes(G4VisAttributes::Invisible);
1436 logicWaterPipe->SetVisAttributes(G4VisAttributes::Invisible);
1440 G4VisAttributes* ringVisAtt=
new G4VisAttributes(G4Colour(0.5,0.25,0.));
1441 ringVisAtt->SetVisibility(
false);
1443 logicSupportBar->SetVisAttributes(ringVisAtt);
1444 if(logicSupportBar1)
1445 logicSupportBar1->SetVisAttributes(ringVisAtt);
1447 logicEndRing->SetVisAttributes(ringVisAtt);
1449 logicGear->SetVisAttributes(ringVisAtt);
1451 logicTaperRing1->SetVisAttributes(ringVisAtt);
1453 logicTaperRing2->SetVisAttributes(ringVisAtt);
1455 logicTaperRing3->SetVisAttributes(ringVisAtt);
1459 G4VisAttributes* endPhiVisAtt=
new G4VisAttributes(G4Colour(0,1.0,0));
1460 endPhiVisAtt->SetVisibility(
false);
1462 logicEnd->SetVisAttributes(endPhiVisAtt);
1483 for(
int i = 0; i < 44; i++){
1484 std::ostringstream osnameBSCCasing;
1485 osnameBSCCasing <<
"logicalBSCCasing"<<i;
1489 G4VisAttributes* rightVisAtt=
new G4VisAttributes(G4Colour(1.0, 0.0,0.0));
1490 rightVisAtt->SetVisibility(
false);
1491 logicBSCTheta->SetVisAttributes(rightVisAtt);
1494 std::ostringstream osnameBSCCable1;
1495 osnameBSCCable1 <<
"logicalBSCCable_1_"<<i;
1498 logicCable->SetVisAttributes(G4VisAttributes::Invisible);
1500 std::ostringstream osnameBSCCable2;
1501 osnameBSCCable2 <<
"logicalBSCCable_2_"<<i;
1504 logicCable->SetVisAttributes(G4VisAttributes::Invisible);
1506 std::ostringstream osnameOCGirder1;
1507 osnameOCGirder1 <<
"logicalOpenningCutGirder_1_"<<i;
1510 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
1512 std::ostringstream osnameOCGirder2;
1513 osnameOCGirder2 <<
"logicalOpenningCutGirder_2_"<<i;
1516 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
1518 std::ostringstream osnameOCGirder3;
1519 osnameOCGirder3 <<
"logicalOpenningCutGirder_3_"<<i;
1522 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
1524 std::ostringstream osnameOCGirder4;
1525 osnameOCGirder4 <<
"logicalOpenningCutGirder_4_"<<i;
1528 logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
1545 for(G4int sector=0;sector<3;sector++) {
1546 std::ostringstream osnameEndPhi;
1547 osnameEndPhi<<
"logicalEndPhi"<<sector;
1550 logicEndPhi->SetVisAttributes(G4VisAttributes::Invisible);
1552 G4cout<<
"Can't find logicEndPhi!"<<G4endl;
1555 for(G4int cryNb=0;cryNb<35;cryNb++) {
1557 std::ostringstream osnameEndCrystal;
1558 osnameEndCrystal<<
"logicalEndCrystal_"<<sector<<
"_"<<cryNb;
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);
1567 G4cout<<
"Can't find: "<<osnameEndCrystal.str()<<G4endl;
1570 std::ostringstream osnameEndCasing;
1571 osnameEndCasing<<
"logicalEndCasing_"<<sector<<
"_"<<cryNb;
1573 if(logicEndCasing) {
1574 logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
1576 G4cout<<
"Can't find: "<<osnameEndCasing.str()<<G4endl;
1582void BesEmcConstruction::DefineMaterials()
1584 G4String name, symbol;
1585 G4double a, z, density;
1589 G4int ncomponents, natoms;
1590 G4double fractionmass;
1599 G4Element*
H=G4Element::GetElement(
"Hydrogen");
1603 H =
new G4Element(name=
"Hydrogen",symbol=
"H" , z= 1., a);
1605 G4Element*
C=G4Element::GetElement(
"Carbon");
1609 C =
new G4Element(name=
"Carbon" ,symbol=
"C" , z= 6., a);
1611 G4Element* O=G4Element::GetElement(
"Oxygen");
1615 O =
new G4Element(name=
"Oxygen" ,symbol=
"O" , z= 8., a);
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);
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);
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);
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);
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);
1644 G4Material *H2O = G4Material::GetMaterial(
"Water");
1645 G4Material *Cu = G4Material::GetMaterial(
"Copper");
1646 G4double dWater = 1.*g/cm3;
1647 G4double dCopper = 8.96*g/cm3;
1648 G4double aWater = ((*besEMCGeometry).waterPipeDr-(*besEMCGeometry).waterPipeThickness)
1649 *((*besEMCGeometry).waterPipeDr-(*besEMCGeometry).waterPipeThickness);
1650 G4double aCopper = (*besEMCGeometry).waterPipeDr*(*besEMCGeometry).waterPipeDr-aWater;
1651 density = (dWater*aWater+dCopper*aCopper)/(aWater+aCopper);
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);
1659 cable =
new G4Material(name=
"M_Cable", density=4.*g/cm3, ncomponents=1);
1660 cable->AddMaterial(Cu,1);
1668 G4Material* Al=G4Material::GetMaterial(
"Aluminium");
1671 Al =
new G4Material(name=
"Aluminium", z=13., a=26.98*g/mole, density=2.700*g/cm3);
1674 G4Material *Si=G4Material::GetMaterial(
"M_Silicon");
1677 Si =
new G4Material(name=
"M_Silicon", z=14., a=28.0855*g/mole, density=2.33*g/cm3);
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)
1688 G4Material* Casing =
new G4Material(name=
"M_Casing", density, ncomponents=3);
1689 Casing->AddMaterial(
1691 fractionmass=Tyvek->GetDensity()/density
1692 *(*besEMCGeometry).fTyvekThickness
1694 Casing->AddMaterial(
1696 fractionmass=Al->GetDensity()/density
1697 *(*besEMCGeometry).fAlThickness
1699 Casing->AddMaterial(
1701 fractionmass=Mylar->GetDensity()/density
1702 *(*besEMCGeometry).fMylarThickness
1704 fCasingMaterial = Casing;
1705 rearCasingMaterial = Tyvek;
1708 fCrystalMaterial = G4Material::GetMaterial(
"Cesiumiodide");
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;
1741 G4Material* pttoMaterial = G4Material::GetMaterial(materialChoice);
1743 {fCrystalMaterial = pttoMaterial;
1744 logicBSCCrystal->SetMaterial(pttoMaterial);
1754 G4Material* pttoMaterial = G4Material::GetMaterial(materialChoice);
1756 {fCasingMaterial = pttoMaterial;
1757 logicBSCTheta->SetMaterial(pttoMaterial);
1767 (*besEMCGeometry).fTyvekThickness = val(
'X');
1768 (*besEMCGeometry).fAlThickness = val(
'Y');
1769 (*besEMCGeometry).fMylarThickness = val(
'Z');
1776 (*besEMCGeometry).BSCRmin = val;
1783 (*besEMCGeometry).BSCNbPhi = val;
1790 (*besEMCGeometry).BSCNbTheta = val;
1803 (*besEMCGeometry).BSCCryLength = val;
1810 (*besEMCGeometry).BSCYFront0 = val;
1817 (*besEMCGeometry).BSCYFront = val;
1824 (*besEMCGeometry).BSCPosition0 = val;
1831 (*besEMCGeometry).BSCPosition1 = val;
1840 G4FieldManager* fieldMgr
1841 = G4TransportationManager::GetTransportationManager()->GetFieldManager();
1843 if(magField)
delete magField;
1846 { magField =
new G4UniformMagField(G4ThreeVector(0.,0.,fieldValue));
1847 fieldMgr->SetDetectorField(magField);
1848 fieldMgr->CreateChordFinder(magField);
1849 fmagField=fieldValue;
1852 fieldMgr->SetDetectorField(magField);
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() );
1877 HepPlane3D bottomPlane( point[4], point[5], point[6] );
1878 HepPoint3D centerProject = bottomPlane.point( center );
1879 Hep3Vector newZ = center - centerProject;
1882 G4RotationMatrix *g4Rot =
new G4RotationMatrix();
1883 g4Rot->rotateX( rotAngle.x() );
1884 g4Rot->rotateY( rotAngle.y() );
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]);
1899 for (
int i = 0; i < 8; i++) {
1911 Hep3Vector x0(1, 0, 0), y0(0, 1, 0), z0(0, 0, 1);
1912 double dx, dy, dz = 0.0;
1914 Hep3Vector a(0.0, newZ.y(), newZ.z());
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;
1926 Hep3Vector b(0,
sin(dx),
cos(dx));
1927 dy = acos(b.dot(newZ));
1928 if(newZ.x() > 0.0) dy = -dy;
1930 Hep3Vector
v(dx, dy, dz);
double tan(const BesAngle a)
double sin(const BesAngle a)
double cos(const BesAngle a)
**********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
***************************************************************************************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
HepGeom::Point3D< double > HepPoint3D
HepGeom::Plane3D< double > HepPlane3D
void SetBSCRmin(G4double)
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 ¢er, G4ThreeVector &rotAngle)
void SetMagField(G4double)
void ConstructEndGeometry(G4LogicalVolume *)
void SetBSCPosition1(G4double)
void SetCasingThickness(G4ThreeVector)
void SetBSCYFront0(G4double)
void SetBSCCrystalLength(G4double)
void PrintEMCParameters()
void SetCrystalMaterial(G4String)
void SetBSCPosition0(G4double)
Hep3Vector RotAngleFromNewZ(Hep3Vector newZ)
void SetCasingMaterial(G4String)
void SetStartIDTheta(G4int)
G4int ComputeEndCopyNb(G4int)
void SetBSCNbTheta(G4int)
void ModifyForCasing(G4ThreeVector pos[8], G4int CryNb)
void Zoom(const G4ThreeVector pos[8], const G4double factor)
void ComputeEMCParameters()
G4LogicalVolume * FindLogicalVolume(const G4String &vn)
static EmcG4Geo * Instance()
Get a pointer to the single instance of EmcG4Geo.
G4LogicalVolume * GetTopVolume()
Get the top(world) volume;.