/home/bes3soft/bes3soft/Boss/7.0.2/dist/7.0.2/Simulation/BOOST/EmcSim/EmcSim-00-00-46/src/BesEmcConstruction.cc

Go to the documentation of this file.
00001 
00002 
00003 //---------------------------------------------------------------------------//
00004 //      BOOST --- BESIII Object_Oreiented Simulation Tool                    //
00005 //---------------------------------------------------------------------------//
00006 //Descpirtion: EMC detector 
00007 //Author: Fu Chengdong
00008 //Created: Sep 4, 2003
00009 //Comment:
00010 //---------------------------------------------------------------------------//
00011 //
00012 #include "G4ThreeVector.hh"
00013 #include "G4ReflectionFactory.hh"
00014 
00015 #include "BesEmcConstruction.hh"
00016 #include "BesEmcDetectorMessenger.hh"
00017 #include "BesEmcGeometry.hh"
00018 #include "BesCrystalParameterisation.hh"
00019 #include "BesEmcEndGeometry.hh"
00020 #include "ReadBoostRoot.hh"
00021 #include "BesEmcParameter.hh"
00022 
00023 #include "BesEmcSD.hh"
00024 #include "G4IrregBox.hh"
00025 #include "G4Box.hh"
00026 #include "G4Transform3D.hh"
00027 #include "G4Tubs.hh"
00028 #include "G4Cons.hh"
00029 #include "G4Trap.hh"
00030 #include "G4UnionSolid.hh"
00031 #include "G4SubtractionSolid.hh"
00032 #include "G4Polyhedra.hh"
00033 #include "G4LogicalVolume.hh"
00034 #include "G4VPhysicalVolume.hh"
00035 #include "G4Material.hh"
00036 #include "G4PVPlacement.hh"
00037 #include "G4PVParameterised.hh"
00038 #include "G4PVReplica.hh"
00039 #include "globals.hh"
00040 #include "G4UniformMagField.hh"
00041 #include "G4FieldManager.hh"
00042 #include "G4TransportationManager.hh"
00043 #include "G4SDManager.hh"
00044 #include "G4RunManager.hh"
00045 #include "G4VisAttributes.hh"
00046 #include "G4Color.hh"
00047 #include "G4AffineTransform.hh"
00048 
00049 #include "G4ios.hh"
00050 #include "G4Geo/EmcG4Geo.h"
00051 
00052 BesEmcConstruction* BesEmcConstruction::fBesEmcConstruction=0;
00053 
00054 BesEmcConstruction* BesEmcConstruction::GetBesEmcConstruction()
00055 { return fBesEmcConstruction;}
00056 
00057 BesEmcConstruction::BesEmcConstruction()
00058   :verboseLevel(0),
00059    solidEMC(0),logicEMC(0),physiEMC(0),logicBSCWorld(0),
00060    solidBSCPhi(0),logicBSCPhi(0),physiBSCPhi(0),
00061    solidBSCTheta(0),logicBSCTheta(0),physiBSCTheta(0),
00062    solidBSCCrystal(0),logicBSCCrystal(0),physiBSCCrystal(0),
00063    magField(0),detectorMessenger(0),besEMCSD(0),crystalParam(0),
00064    logicEnd(0),logicEndPhi(0),logicEndCasing(0),logicEndCrystal(0),
00065    logicRear(0),logicRearCasing(0),logicOrgGlass(0),logicPD(0),
00066    logicAlPlate(0),logicPreAmpBox(0),logicAirInPABox(0),
00067    logicHangingPlate(0),logicOCGirder(0),logicCable(0),logicWaterPipe(0),
00068    logicSupportBar(0),logicSupportBar1(0),logicEndRing(0),logicGear(0),
00069    logicTaperRing1(0),logicTaperRing2(0),logicTaperRing3(0)
00070 {
00071   if(fBesEmcConstruction)
00072     { G4Exception("BesEmcConstruction constructed twice."); }
00073   fBesEmcConstruction=this;
00074   //for debug
00075   //  G4Exception("BesEmcConstruction::BesEmcConstruction() starting........");
00076   startID           = 1;
00077   phiNbCrystals     = 0;
00078   thetaNbCrystals   = 0;
00079   besEMCGeometry = new BesEmcGeometry();
00080   emcEnd = new BesEmcEndGeometry();
00081 }
00082 
00083 BesEmcConstruction::~BesEmcConstruction()
00084 {
00085   if(detectorMessenger) delete detectorMessenger;
00086   if(crystalParam) delete crystalParam;
00087   if(besEMCGeometry) delete besEMCGeometry;
00088   if(emcEnd) delete emcEnd;
00089   
00090   BesEmcParameter::Kill();
00091 }
00092 
00093 void BesEmcConstruction::Construct(G4LogicalVolume* logicBes)
00094 {
00095   besEMCGeometry->ComputeEMCParameters();
00096   detectorMessenger = new BesEmcDetectorMessenger(this,besEMCGeometry);
00097   emcEnd->ComputeParameters();
00098 
00099   G4SDManager* SDman = G4SDManager::GetSDMpointer();
00100   if (!besEMCSD) {
00101     besEMCSD = new BesEmcSD("CalorSD",this,besEMCGeometry);
00102     SDman->AddNewDetector( besEMCSD );
00103   }
00104 
00105   // Construction
00106   G4cout<<"--------ReadBoostRoot::GetEmc()="<<ReadBoostRoot::GetEmc()<<G4endl;
00107   if(ReadBoostRoot::GetEmc()==2)
00108   {
00109     logicEMC = EmcG4Geo::Instance()->GetTopVolume();
00110 
00111     if(logicEMC){
00112       physiEMC = new G4PVPlacement(0,
00113           G4ThreeVector(0.0 ,0.0 ,0.0),
00114           logicEMC, "physicalEMC",logicBes, false, 0);
00115       G4cout<<"logicEmc:  ===  "<<logicEMC<<"  physiEmc "<<physiEMC<<G4endl;
00116 
00117       GetLogicalVolume();
00118       SetVisAndSD();
00119     }
00120   }
00121   else {
00122     //for debug
00123     //  G4Exception("BesEmcConstruction::Construct() starting............");
00124     //
00125     DefineMaterials();
00126     phiNbCrystals     = (*besEMCGeometry).BSCNbPhi;
00127     thetaNbCrystals   = (*besEMCGeometry).BSCNbTheta*2;
00128 
00129     G4double da=0.001*deg;  //delta angle to avoid overlap
00130     
00131     //
00132     //BSC
00133     //
00134     solidBSC = new G4Tubs("solidBSC",
00135       (*besEMCGeometry).TaperRingRmin1,
00136       (*besEMCGeometry).BSCRmax+(*besEMCGeometry).SPBarThickness+(*besEMCGeometry).SPBarThickness1+2.1*mm,   //radius from 942mm to 940 mm
00137       (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz,
00138                         0.*deg,
00139                         360.*deg);
00140 
00141     solidESC = new G4Cons("solidESC",(*emcEnd).WorldRmin1,(*emcEnd).WorldRmax1,
00142         (*emcEnd).WorldRmin2,(*emcEnd).WorldRmax2,
00143         (*emcEnd).WorldDz/2,0.*deg,360.*deg);
00144 
00145     solidEMC = new G4UnionSolid("solidEMC0",
00146         solidBSC,
00147         solidESC,
00148         0,
00149         G4ThreeVector(0,0,(*emcEnd).WorldZPosition));
00150     
00151     G4RotationMatrix *rotateESC = new G4RotationMatrix();
00152     rotateESC->rotateY(180.*deg);
00153     
00154     solidEMC = new G4UnionSolid("solidEMC",
00155         solidEMC,
00156         solidESC,
00157         rotateESC,
00158         G4ThreeVector(0,0,-(*emcEnd).WorldZPosition));
00159     
00160     logicEMC = new G4LogicalVolume(solidEMC,
00161                          G4Material::GetMaterial("Air"),
00162                          "logicalEMC");
00163     
00164     physiEMC = new G4PVPlacement(0,
00165         0,
00166         logicEMC,
00167         "physicalEMC",
00168         logicBes,
00169         false,
00170         0);
00171 
00172     solidBSCWorld = new G4SubtractionSolid("solidBSCWorld0",
00173           solidBSC,
00174           solidESC,
00175           0,
00176           G4ThreeVector(0,0,(*emcEnd).WorldZPosition));
00177 
00178     solidBSCWorld = new G4SubtractionSolid("solidBSCWorld",
00179         solidBSCWorld,
00180         solidESC,
00181         rotateESC,
00182         G4ThreeVector(0,0,-(*emcEnd).WorldZPosition));
00183 
00184     logicBSCWorld = new G4LogicalVolume(solidBSCWorld,
00185         G4Material::GetMaterial("Air"),
00186         "logicalBSCWorld");
00187 
00188     G4RotationMatrix *rotBSC = new G4RotationMatrix();
00189     rotBSC->rotateY(180.*deg);
00190     physiBSCWorld = new G4PVPlacement(rotBSC,
00191         0,
00192         logicBSCWorld,
00193         "physicalBSCWorld",
00194         logicEMC,
00195         false,
00196         0);
00197     
00198     G4RotationMatrix *rotateMatrix[200];
00199     G4double oOp,ox,oy,oz; 
00200     G4double delta = 0*deg;
00201     G4ThreeVector axis = G4ThreeVector(0,0,0);
00202     oOp=(*besEMCGeometry).BSCRmin/sin(0.5*(*besEMCGeometry).BSCPhiDphi+90*deg)
00203       *sin((*besEMCGeometry).BSCAngleRotat);
00204     G4double ll=(*besEMCGeometry).BSCCryLength;
00205     G4double rr=(*besEMCGeometry).BSCRmin;
00206     G4double oj=sqrt(ll*ll+rr*rr-2*ll*rr*cos(180.*deg-(*besEMCGeometry).BSCAngleRotat));
00207     G4double oij=90.*deg-(*besEMCGeometry).BSCPhiDphi/2.-(*besEMCGeometry).BSCAngleRotat;
00208     G4double doj=asin(sin(180.*deg-(*besEMCGeometry).BSCAngleRotat)/oj*ll);
00209     G4double ioj=(*besEMCGeometry).BSCPhiDphi/2.+doj;
00210     G4double ij=oj/sin(oij)*sin(ioj);
00211     G4double dOp=rr/sin(90.*deg-(*besEMCGeometry).BSCPhiDphi/2.)
00212       *sin(90.*deg+(*besEMCGeometry).BSCPhiDphi/2.-(*besEMCGeometry).BSCAngleRotat);
00213     G4double cOp=rr/sin(90.*deg+(*besEMCGeometry).BSCPhiDphi/2.)
00214       *sin(90.*deg-(*besEMCGeometry).BSCPhiDphi/2.-(*besEMCGeometry).BSCAngleRotat);
00215     G4double ch=(dOp+ll)/cos((*besEMCGeometry).BSCPhiDphi)-cOp;
00216     G4double hi=(dOp+ll)*tan((*besEMCGeometry).BSCPhiDphi)-ij;
00217     G4double oh=sqrt(ch*ch+rr*rr-2*ch*rr*cos(180*deg-(*besEMCGeometry).BSCAngleRotat));
00218     G4double hoi=asin(sin(180*deg-oij)/oh*hi);
00219     G4double dok=asin(sin(180*deg-(*besEMCGeometry).BSCAngleRotat)/oh*ch);
00220     if(verboseLevel>3)
00221       G4cout << "oj=" <<oj/cm<<G4endl
00222         << "oij="<<oij/deg<<G4endl
00223         << "doj="<<doj/deg<<G4endl
00224         << "ioj="<<ioj/deg<<G4endl
00225         << "ij="<<ij/cm<<G4endl
00226         << "dOp="<<dOp/cm<<G4endl
00227         << "cOp="<<cOp/cm<<G4endl
00228         << "ch="<<ch/cm<<G4endl
00229         << "hi="<<hi/cm<<G4endl
00230         << "oh="<<oh/cm<<G4endl
00231         << "hoi="<<hoi/deg<<G4endl
00232         << "dok="<<dok/deg<<G4endl;
00233 
00234     // Phi Cell
00235     solidBSCPhiTub = new G4Tubs(
00236         "solidBSCPhiTub",
00237         dOp,
00238         (*besEMCGeometry).BSCPhiRmax,
00239         (*besEMCGeometry).BSCDz,
00240         360.*deg-(*besEMCGeometry).BSCPhiDphi,
00241         (*besEMCGeometry).BSCPhiDphi);
00242     solidConsPhi = new G4Cons("consPhi1",
00243         (*besEMCGeometry).BSCRmin1,
00244         (*besEMCGeometry).BSCRmax1,
00245         (*besEMCGeometry).BSCRmin2,
00246         (*besEMCGeometry).BSCRmax2,
00247         (*besEMCGeometry).BSCDz1/2,
00248         0.*deg,
00249         360.*deg);
00250     solidBSCPhi1 = new G4SubtractionSolid("solidBSCPhi1",
00251         solidBSCPhiTub,
00252         solidConsPhi,
00253         0,
00254         G4ThreeVector(0,0,(*besEMCGeometry).BSCDz-(*besEMCGeometry).BSCDz1/2));
00255     solidConsPhi = new G4Cons("consPhi2",
00256         (*besEMCGeometry).BSCRmin2,
00257         (*besEMCGeometry).BSCRmax2,
00258         (*besEMCGeometry).BSCRmin1,
00259         (*besEMCGeometry).BSCRmax1,
00260         (*besEMCGeometry).BSCDz1/2,
00261         0.*deg,
00262         360.*deg);
00263     solidBSCPhi = new G4SubtractionSolid("solidBSCPhi",
00264         solidBSCPhi1,
00265         solidConsPhi,
00266         0,
00267         G4ThreeVector(0,0,(*besEMCGeometry).BSCDz1/2-(*besEMCGeometry).BSCDz));
00268 
00269     logicBSCPhi = new G4LogicalVolume(solidBSCPhi,
00270         G4Material::GetMaterial("Air"),
00271         "logicalBSCPhi");
00272 
00273     G4int i;
00274     for(G4int j=0;j<(*besEMCGeometry).BSCNbPhi;j++)  //=============
00275     {
00276       if(j<(*besEMCGeometry).BSCNbPhi/2) {  //0~59
00277         i=(*besEMCGeometry).BSCNbPhi/2-j-1;
00278       } else {  //60~119
00279         i=(*besEMCGeometry).BSCNbPhi*3/2-j-1;
00280       }
00281       rotateMatrix[i] = new G4RotationMatrix();
00282       rotateMatrix[i]->rotateZ(-i*(*besEMCGeometry).BSCPhiDphi
00283           -(*besEMCGeometry).BSCAngleRotat
00284           -(*besEMCGeometry).BSCPhiDphi/2.
00285           -hoi);
00286       rotateMatrix[i]->getAngleAxis(delta, axis);
00287       //G4cout << "The axis of crystals in the world system is: "
00288       //   << delta/deg << "(deg)(delta) "
00289       //<< axis  << "(Z axis)"<< G4endl;
00290       ox=oOp*cos(-90.*deg+(*besEMCGeometry).BSCAngleRotat+hoi
00291           +i*(*besEMCGeometry).BSCPhiDphi);
00292       oy=oOp*sin(-90.*deg+(*besEMCGeometry).BSCAngleRotat+hoi
00293           +i*(*besEMCGeometry).BSCPhiDphi);
00294       oz=0*cm;
00295 
00296       ostringstream strPhi;
00297       strPhi << "physicalBSCPhi"  << j;
00298 
00299       physiBSCPhi = new G4PVPlacement(rotateMatrix[i],
00300           G4ThreeVector(ox,oy,oz),
00301           logicBSCPhi,
00302           strPhi.str(),
00303           logicBSCWorld,
00304           false,
00305           j);
00306       //G4cout << G4ThreeVector(ox/cm,oy/cm,oz/cm) <<"(cm)" << G4endl
00307       //   << (-(*besEMCGeometry).BSCAngleRotat+(i-1)*(*besEMCGeometry).BSCPhiDphi)/deg <<"(degree)" << G4endl;
00308     }
00309 
00310     //
00311     //Crystals
00312     //
00313     G4double zHalfLength[50];
00314     G4double thetaAxis[50];
00315     G4double phiAxis[50];
00316     G4double yHalfLength1[50];
00317     G4double xHalfLength2[50];
00318     G4double xHalfLength1[50];
00319     G4double tanAlpha1[50];
00320     G4double yHalfLength2[50];
00321     G4double xHalfLength4[50];
00322     G4double xHalfLength3[50];
00323     G4double tanAlpha2[50];
00324     G4double xPosition[50];
00325     G4double yPosition[50];
00326     G4double zPosition[50];
00327     G4double thetaPosition[50];
00328     for(i=0;i<(*besEMCGeometry).BSCNbTheta;i++)
00329     {
00330       zHalfLength[i]  = (*besEMCGeometry).zHalfLength[i];
00331       thetaAxis[i]    = (*besEMCGeometry).thetaAxis[i];
00332       phiAxis[i]      = (*besEMCGeometry).phiAxis[i];
00333       yHalfLength1[i] = (*besEMCGeometry).yHalfLength1[i];
00334       xHalfLength2[i] = (*besEMCGeometry).xHalfLength2[i];
00335       xHalfLength1[i] = (*besEMCGeometry).xHalfLength1[i];
00336       tanAlpha1[i]    = (*besEMCGeometry).tanAlpha1[i];
00337       yHalfLength2[i] = (*besEMCGeometry).yHalfLength2[i];
00338       xHalfLength4[i] = (*besEMCGeometry).xHalfLength4[i];
00339       xHalfLength3[i] = (*besEMCGeometry).xHalfLength3[i];
00340       tanAlpha2[i]    = (*besEMCGeometry).tanAlpha2[i];
00341       xPosition[i]    = (*besEMCGeometry).xPosition[i];
00342       yPosition[i]    = (*besEMCGeometry).yPosition[i];
00343       zPosition[i]    = (*besEMCGeometry).zPosition[i];
00344       thetaPosition[i]= (*besEMCGeometry).thetaPosition[i];
00345       if(verboseLevel>4)        
00346         G4cout << "The sizes of the "<<i+1<<" crystal are:" << G4endl
00347           <<"zHalfLength ="<<zHalfLength[i]/cm<< "(cm)," << G4endl
00348           << "thetaAxis  ="<<thetaAxis[i]/deg << "(deg),"<< G4endl
00349           << "phiAxis    ="<< phiAxis[i]/deg  << "(deg),"<< G4endl
00350           << "yHalfLength1="<<yHalfLength1[i]/cm<<"(cm),"<< G4endl
00351           << "xHalfLength1="<<xHalfLength1[i]/cm<<"(cm),"<< G4endl
00352           << "xHalfLength2="<<xHalfLength2[i]/cm<<"(cm),"<< G4endl
00353           << "tanAlpha1   ="<< tanAlpha1[i]              << G4endl
00354           << "yHalfLength2="<<yHalfLength2[i]/cm<<"(cm),"<< G4endl
00355           << "xHalfLength3="<<xHalfLength3[i]/cm<<"(cm),"<< G4endl
00356           << "xHalfLength4="<<xHalfLength4[i]/cm<<"(cm),"<< G4endl
00357           << "tanAlpha2   =" << tanAlpha2[i]    << "."   << G4endl;
00358     }
00359     besEMCGeometry->ModifyForCasing();
00360 
00361     solidBSCCrystal = new G4Trap("solidCrystal",
00362         100*cm, 100*deg, 100*deg,
00363         100*cm, 100*cm,  100*cm,  100*deg,
00364         100*cm, 100*cm,  100*cm,  100*deg);
00365 
00366     logicBSCCrystal = new G4LogicalVolume(solidBSCCrystal,
00367         fCrystalMaterial,
00368         "logicalCrystal");
00369 
00370     crystalParam = new BesCrystalParameterisation
00371       (startID,
00372        thetaNbCrystals,
00373        (*besEMCGeometry).BSCNbTheta*2,
00374        besEMCGeometry,
00375        verboseLevel);
00376 
00377     //---------------------------------------------------------------------------------
00378     //rear substance
00379     solidRear = new G4Box("solidRearBox",
00380         (*besEMCGeometry).rearBoxLength/2,
00381         (*besEMCGeometry).rearBoxLength/2,
00382         (*besEMCGeometry).rearBoxDz/2);
00383 
00384     logicRear = new G4LogicalVolume(solidRear,
00385         G4Material::GetMaterial("Air"),
00386         "logicalRearBox");
00387 
00388     //organic glass
00389     solidOrgGlass = new G4Box("solidOrganicGlass",
00390         (*besEMCGeometry).orgGlassLengthX/2,
00391         (*besEMCGeometry).orgGlassLengthY/2,
00392         (*besEMCGeometry).orgGlassLengthZ/2);
00393 
00394     logicOrgGlass = new G4LogicalVolume(solidOrgGlass,
00395         organicGlass,
00396         "logicalOrganicGlass");
00397 
00398     physiOrgGlass = new G4PVPlacement(0,
00399         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz-(*besEMCGeometry).orgGlassLengthZ)/2),
00400         logicOrgGlass,
00401         "physicalOrganicGlass",
00402         logicRear,
00403         false,
00404         0);
00405 
00406     //casing
00407     solidCasingBox = new G4Box("solidCasingBox",
00408         (*besEMCGeometry).rearBoxLength/2,
00409         (*besEMCGeometry).rearBoxLength/2,
00410         (*besEMCGeometry).rearCasingThickness/2);
00411 
00412     solidAirHole = new G4Box("solidAirHole",
00413         (*besEMCGeometry).orgGlassLengthX/2,
00414         (*besEMCGeometry).orgGlassLengthY/2,
00415         (*besEMCGeometry).rearBoxDz/2);      //any value more than casing thickness
00416 
00417     solidRearCasing = new G4SubtractionSolid("solidRearCasing",
00418         solidCasingBox,
00419         solidAirHole,
00420         0,
00421         0);
00422 
00423     logicRearCasing = new G4LogicalVolume(solidRearCasing,
00424         rearCasingMaterial,
00425         "logicalRearCasing");
00426 
00427     physiRearCasing = new G4PVPlacement(0,
00428         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz-(*besEMCGeometry).rearCasingThickness)/2),
00429         logicRearCasing,
00430         "physicalRearCasing",
00431         logicRear,
00432         false,
00433         0);
00434 
00435     //Al Plate
00436     solidAlBox = new G4Box("solidAlBox",
00437         (*besEMCGeometry).rearBoxLength/2,
00438         (*besEMCGeometry).rearBoxLength/2,
00439         (*besEMCGeometry).AlPlateDz/2);
00440 
00441     solidAlPlate = new G4SubtractionSolid("solidAlPlate",
00442         solidAlBox,
00443         solidAirHole,
00444         0,
00445         0);
00446 
00447     logicAlPlate = new G4LogicalVolume(solidAlPlate,
00448         G4Material::GetMaterial("Aluminium"),
00449         "logicalAlPlate");
00450 
00451     physiAlPlate = new G4PVPlacement(0,
00452         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
00453             -(*besEMCGeometry).rearCasingThickness
00454             -(*besEMCGeometry).AlPlateDz/2)),
00455         logicAlPlate,
00456         "physicalAlPlate",
00457         logicRear,
00458         false,
00459         0);
00460 
00461     //photodiode
00462     solidPD = new G4Box("solidPD",
00463         (*besEMCGeometry).PDLengthX,     //two PD
00464         (*besEMCGeometry).PDLengthY/2,
00465         (*besEMCGeometry).PDLengthZ/2);
00466 
00467     logicPD = new G4LogicalVolume(solidPD,
00468         G4Material::GetMaterial("M_Silicon"),
00469         "logicalPD");
00470 
00471     physiPD = new G4PVPlacement(0,
00472         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
00473             -(*besEMCGeometry).orgGlassLengthZ
00474             -(*besEMCGeometry).PDLengthZ/2)),
00475         logicPD,
00476         "physicalPD",
00477         logicRear,
00478         false,
00479         0);
00480 
00481     //preamplifier box
00482     solidPreAmpBox = new G4Box("solidPreAmpBox",
00483         (*besEMCGeometry).rearBoxLength/2,
00484         (*besEMCGeometry).rearBoxLength/2,
00485         (*besEMCGeometry).PABoxDz/2);
00486 
00487     logicPreAmpBox = new G4LogicalVolume(solidPreAmpBox,
00488         G4Material::GetMaterial("Aluminium"),
00489         "logicalPreAmpBox");
00490 
00491     physiPreAmpBox = new G4PVPlacement(0,
00492         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
00493             -(*besEMCGeometry).rearCasingThickness
00494             -(*besEMCGeometry).AlPlateDz
00495             -(*besEMCGeometry).PABoxDz/2)),
00496         logicPreAmpBox,
00497         "physicalPreAmpBox",
00498         logicRear,
00499         false,
00500         0);
00501 
00502     //air in preamplifier box
00503     solidAirInPABox = new G4Box("solidAirInPABox",
00504         (*besEMCGeometry).rearBoxLength/2-(*besEMCGeometry).PABoxThickness,
00505         (*besEMCGeometry).rearBoxLength/2-(*besEMCGeometry).PABoxThickness,
00506         (*besEMCGeometry).PABoxDz/2-(*besEMCGeometry).PABoxThickness);
00507 
00508     logicAirInPABox = new G4LogicalVolume(solidAirInPABox,
00509         G4Material::GetMaterial("Air"),
00510         "logicalAirInPABox");
00511 
00512     physiAirInPABox = new G4PVPlacement(0,
00513         0,
00514         logicAirInPABox,
00515         "physicalAirInPABox",
00516         logicPreAmpBox,
00517         false,
00518         0);
00519 
00520     //stainless steel for hanging the crystal
00521     solidHangingPlate = new G4Box("solidHangingPlate",
00522         (*besEMCGeometry).rearBoxLength/2,
00523         (*besEMCGeometry).rearBoxLength/2,
00524         (*besEMCGeometry).HangingPlateDz/2);
00525 
00526     logicHangingPlate = new G4LogicalVolume(solidHangingPlate,stainlessSteel,"logicalHangingPlate");
00527 
00528     physiHangingPlate = new G4PVPlacement(0,
00529         G4ThreeVector(0,0,-((*besEMCGeometry).rearBoxDz/2
00530             -(*besEMCGeometry).rearCasingThickness
00531             -(*besEMCGeometry).AlPlateDz
00532             -(*besEMCGeometry).PABoxDz
00533             -(*besEMCGeometry).HangingPlateDz/2)),
00534         logicHangingPlate,
00535         "physicalHangingPlate",
00536         logicRear,
00537         false,
00538         0);
00539 
00540     //water pipe
00541     solidWaterPipe = new G4Tubs("solidWaterPipe",
00542         0,
00543         (*besEMCGeometry).waterPipeDr,
00544         (*besEMCGeometry).BSCDz,
00545         0.*deg,
00546         360.*deg);
00547 
00548     logicWaterPipe = new G4LogicalVolume(solidWaterPipe,waterPipe,"logicalWaterPipe");
00549 
00550     physiWaterPipe = new G4PVPlacement(0,
00551         G4ThreeVector((*besEMCGeometry).cablePosX[0]-2*(*besEMCGeometry).cableDr,
00552           (*besEMCGeometry).cablePosY[0]-(*besEMCGeometry).cableDr-(*besEMCGeometry).waterPipeDr,
00553           0),
00554         logicWaterPipe,
00555         "physicalWaterPipe",
00556         logicBSCPhi,
00557         false,
00558         0);
00559     //---------------------------------------------------------------------------------
00560 
00561 
00562     //
00563     //Theta Cell
00564     //
00565     G4String nameCrystalAndCasing="CrystalAndCasing";
00566 
00567     G4int id=0;     //ID of crystals after distinguishing left and right
00568     for(i=startID;i<=thetaNbCrystals;i++) //================
00569     {
00570       ostringstream strSolidCasing;
00571       strSolidCasing << "solidBSCCasing"  << i-1;
00572       ostringstream strVolumeCasing;
00573       strVolumeCasing << "logicalBSCCasing"  << i-1;
00574       ostringstream strPhysiCasing;
00575       strPhysiCasing << "physicalBSCCasing"  << i-1;
00576 
00577       if(i>(*besEMCGeometry).BSCNbTheta)
00578       {
00579         id=i-(*besEMCGeometry).BSCNbTheta-1;
00580         solidBSCTheta = new G4Trap(strSolidCasing.str(),
00581             zHalfLength[id],
00582             thetaAxis[id],
00583             -phiAxis[id],
00584             yHalfLength1[id],
00585             xHalfLength2[id],
00586             xHalfLength1[id],
00587             -tanAlpha1[id],
00588             yHalfLength2[id],
00589             xHalfLength4[id],
00590             xHalfLength3[id],
00591             -tanAlpha2[id]);
00592 
00593         //G4cout<<"in EmcConstr1: "<<strSolidCasing.str()<<" x1="<<xHalfLength1[id]<<" y1="<<yHalfLength1[id]<<" theta="<<thetaAxis[id]
00594         //<<" phi="<<-phiAxis[id]<<" a1="<<-tanAlpha1[id]<<G4endl;
00595 
00596         logicBSCTheta = new G4LogicalVolume(solidBSCTheta,
00597             fCasingMaterial,
00598             strVolumeCasing.str());
00599 
00600         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1] = new G4RotationMatrix();
00601         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1]->rotateZ(-90*deg);
00602         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1]
00603           ->rotateX(-thetaPosition[id]);
00604 
00605 
00606         physiBSCTheta = 
00607           new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
00608               G4ThreeVector(xPosition[id],
00609                 yPosition[id],
00610                 zPosition[id]),
00611               strPhysiCasing.str(),
00612               logicBSCTheta,
00613               physiBSCPhi,
00614               false,
00615               i-1);
00616 
00617         if(logicBSCTheta)
00618         {
00619           G4VisAttributes* rightVisAtt= new G4VisAttributes(G4Colour(1.0,0.,0.));
00620           rightVisAtt->SetVisibility(true);
00621           logicBSCTheta->SetVisAttributes(rightVisAtt);
00622           logicBSCTheta->SetVisAttributes(G4VisAttributes::Invisible);
00623         }
00624 
00625         ostringstream strRear;
00626         strRear << "physicalRearBox_1_"  << i-1;
00627 
00628         physiRear = new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
00629             G4ThreeVector((*besEMCGeometry).rearBoxPosX[id],
00630               (*besEMCGeometry).rearBoxPosY[id],
00631               (*besEMCGeometry).rearBoxPosZ[id]),
00632             strRear.str(),
00633             logicRear,
00634             physiBSCPhi,
00635             false,
00636             i-1);
00637 
00638         ostringstream strGirder;
00639         strGirder << "solidOpenningCutGirder_1_"  << i-1;
00640         solidOCGirder = new G4Cons(strGirder.str(),
00641             (*besEMCGeometry).OCGirderRmin1[id],
00642             (*besEMCGeometry).BSCPhiRmax,
00643             (*besEMCGeometry).OCGirderRmin2[id],
00644             (*besEMCGeometry).BSCPhiRmax,
00645             (*besEMCGeometry).OCGirderDz[id]/2,
00646             360.*deg-(*besEMCGeometry).OCGirderAngle/2,
00647             (*besEMCGeometry).OCGirderAngle/2-da);
00648 
00649         ostringstream strVGirder;
00650         strVGirder << "logicalOpenningCutGirder_1_"  << i-1;
00651         logicOCGirder = new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder.str());
00652         logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
00653 
00654         ostringstream strPGirder;
00655         strPGirder << "physicalOpenningCutGirder_1_"  << i-1;
00656         physiOCGirder = new G4PVPlacement(0,
00657             G4ThreeVector(0,0,(*besEMCGeometry).OCGirderPosZ[id]),
00658             logicOCGirder,
00659             strPGirder.str(),
00660             logicBSCPhi,
00661             false,
00662             0);
00663 
00664         if(id<(*besEMCGeometry).BSCNbTheta-1)
00665         {
00666           G4double zLength = (*besEMCGeometry).OCGirderPosZ[id+1]
00667             -(*besEMCGeometry).OCGirderPosZ[id]
00668             -(*besEMCGeometry).OCGirderDz[id+1]/2-(*besEMCGeometry).OCGirderDz[id]/2;
00669           G4double zPosition = (*besEMCGeometry).OCGirderPosZ[id+1]
00670             -(*besEMCGeometry).OCGirderDz[id+1]/2-zLength/2;
00671 
00672           ostringstream strGirder2;
00673           strGirder2 << "solidOpenningCutGirder_2_"  << i-1;
00674           solidOCGirder = new G4Cons(strGirder2.str(),
00675               (*besEMCGeometry).OCGirderRmin2[id],
00676               (*besEMCGeometry).BSCPhiRmax,
00677               (*besEMCGeometry).OCGirderRmin1[id+1],
00678               (*besEMCGeometry).BSCPhiRmax,
00679               zLength/2,
00680               360.*deg-(*besEMCGeometry).OCGirderAngle/2,
00681               (*besEMCGeometry).OCGirderAngle/2-da);
00682 
00683           ostringstream strVGirder2;
00684           strVGirder2 << "logicalOpenningCutGirder_2_"  << i-1;
00685           logicOCGirder = new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder2.str());
00686           logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
00687 
00688           ostringstream strPGirder2;
00689           strPGirder2 << "physicalOpenningCutGirder_2_"  << i-1;
00690           physiOCGirder = new G4PVPlacement(0,
00691               G4ThreeVector(0,0,zPosition),
00692               logicOCGirder,
00693               strPGirder2.str(),
00694               logicBSCPhi,
00695               false,
00696               0);
00697         } 
00698 
00699         ostringstream strBSCCable;
00700         strBSCCable << "solidBSCCable_1_"  << i-1;
00701         solidCable = new G4Tubs(strBSCCable.str(),
00702             0,
00703             (*besEMCGeometry).cableDr,
00704             (*besEMCGeometry).cableLength[id]/2,
00705             0.*deg,
00706             360.*deg);
00707 
00708         ostringstream strVBSCCable;
00709         strVBSCCable << "logicalBSCCable_1_"  << i-1;
00710         logicCable = new G4LogicalVolume(solidCable,cable,strVBSCCable.str());
00711 
00712         ostringstream strPBSCCable;
00713         strPBSCCable << "physicalBSCCable_1_"  << i-1;
00714         physiCable = new G4PVPlacement(0,
00715             G4ThreeVector((*besEMCGeometry).cablePosX[id],
00716               (*besEMCGeometry).cablePosY[id],
00717               (*besEMCGeometry).cablePosZ[id]),
00718             logicCable,
00719             strPBSCCable.str(),
00720             logicBSCPhi,
00721             false,
00722             0);
00723         logicCable->SetVisAttributes(G4VisAttributes::Invisible);
00724       }
00725       else
00726       {
00727         id=(*besEMCGeometry).BSCNbTheta-i;
00728         solidBSCTheta = new G4Trap(strSolidCasing.str(),
00729             zHalfLength[id],
00730             thetaAxis[id],
00731             phiAxis[id],
00732             yHalfLength1[id],
00733             xHalfLength1[id],
00734             xHalfLength2[id],
00735             tanAlpha1[id],
00736             yHalfLength2[id],
00737             xHalfLength3[id],
00738             xHalfLength4[id],
00739             tanAlpha2[id]);
00740 
00741         //          G4cout<<"in EmcConstr2: "<<strSolidCasing.str()<<" x1="<<xHalfLength1[id]<<" y1="<<yHalfLength1[id]<<" theta="<<thetaAxis[id]
00742         //                <<" phi="<<phiAxis[id]<<" a1="<<tanAlpha1[id]<<G4endl;
00743 
00744         logicBSCTheta = new G4LogicalVolume(solidBSCTheta,
00745             fCasingMaterial,
00746             strVolumeCasing.str());
00747 
00748         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1] = new G4RotationMatrix();
00749         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1]->rotateZ(-90*deg);
00750         rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1]
00751           ->rotateX(-180*deg+thetaPosition[id]);
00752         physiBSCTheta = 
00753           new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
00754               G4ThreeVector(xPosition[id],
00755                 yPosition[id],
00756                 -zPosition[id]),
00757               strPhysiCasing.str(),
00758               logicBSCTheta,
00759               physiBSCPhi,
00760               false,
00761               i-1);
00762         if(logicBSCTheta)
00763         {
00764           G4VisAttributes* rightVisAtt= new G4VisAttributes(G4Colour(1.0,0.,0.));
00765           rightVisAtt->SetVisibility(true);
00766           logicBSCTheta->SetVisAttributes(rightVisAtt);
00767           logicBSCTheta->SetVisAttributes(G4VisAttributes::Invisible);
00768         }
00769 
00770         ostringstream strRear;
00771         strRear << "physicalRearBox_2_"  << i-1;
00772 
00773         physiRear = new G4PVPlacement(rotateMatrix[(*besEMCGeometry).BSCNbPhi+i-1],
00774             G4ThreeVector((*besEMCGeometry).rearBoxPosX[id],
00775               (*besEMCGeometry).rearBoxPosY[id],
00776               -(*besEMCGeometry).rearBoxPosZ[id]),
00777             strRear.str(),
00778             logicRear,
00779             physiBSCPhi,
00780             false,
00781             i-1);
00782 
00783         ostringstream strGirder;
00784         strGirder << "solidOpenningCutGirder_3_"  << i-1;
00785         solidOCGirder = new G4Cons(strGirder.str(),
00786             (*besEMCGeometry).OCGirderRmin2[id],
00787             (*besEMCGeometry).BSCPhiRmax,
00788             (*besEMCGeometry).OCGirderRmin1[id],
00789             (*besEMCGeometry).BSCPhiRmax,
00790             (*besEMCGeometry).OCGirderDz[id]/2,
00791             360.*deg-(*besEMCGeometry).OCGirderAngle/2,
00792             (*besEMCGeometry).OCGirderAngle/2-da);
00793 
00794         ostringstream strVGirder;
00795         strVGirder << "logicalOpenningCutGirder_3_"  << i-1;
00796         logicOCGirder = new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder.str());
00797         logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
00798 
00799         ostringstream strPGirder;
00800         strPGirder << "physicalOpenningCutGirder_3_"  << i-1;
00801         physiOCGirder = new G4PVPlacement(0,
00802             G4ThreeVector(0,0,-(*besEMCGeometry).OCGirderPosZ[id]),
00803             logicOCGirder,
00804             strPGirder.str(),
00805             logicBSCPhi,
00806             false,
00807             0);
00808 
00809         if(id<(*besEMCGeometry).BSCNbTheta-1)
00810         {
00811           G4double zLength = (*besEMCGeometry).OCGirderPosZ[id+1]-(*besEMCGeometry).OCGirderPosZ[id]
00812             -(*besEMCGeometry).OCGirderDz[id+1]/2-(*besEMCGeometry).OCGirderDz[id]/2;
00813           G4double zPosition = (*besEMCGeometry).OCGirderPosZ[id+1]-(*besEMCGeometry).OCGirderDz[id+1]/2-zLength/2;
00814 
00815           ostringstream strGirder2;
00816           strGirder2 << "solidOpenningCutGirder_4_"  << i-1;
00817           solidOCGirder = new G4Cons(strGirder2.str(),
00818               (*besEMCGeometry).OCGirderRmin1[id+1],
00819               (*besEMCGeometry).BSCPhiRmax,
00820               (*besEMCGeometry).OCGirderRmin2[id],
00821               (*besEMCGeometry).BSCPhiRmax,
00822               zLength/2,
00823               360.*deg-(*besEMCGeometry).OCGirderAngle/2,
00824               (*besEMCGeometry).OCGirderAngle/2-da);
00825 
00826           ostringstream strVGirder2;
00827           strVGirder2 << "logicalOpenningCutGirder_4_"  << i-1;
00828           logicOCGirder 
00829             = new G4LogicalVolume(solidOCGirder,stainlessSteel,strVGirder2.str());
00830           logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
00831 
00832           ostringstream strPGirder2;
00833           strPGirder2 << "physicalOpenningCutGirder_4_"  << i-1;
00834           physiOCGirder = new G4PVPlacement(0,
00835               G4ThreeVector(0,0,-zPosition),
00836               logicOCGirder,
00837               strPGirder2.str(),
00838               logicBSCPhi,
00839               false,
00840               0);
00841         }
00842 
00843         ostringstream strBSCCable;
00844         strBSCCable << "solidBSCCable_2_"  << i-1;
00845         solidCable = new G4Tubs(strBSCCable.str(),
00846             0,
00847             (*besEMCGeometry).cableDr,
00848             (*besEMCGeometry).cableLength[id]/2,
00849             0.*deg,
00850             360.*deg);
00851 
00852         ostringstream strVBSCCable;
00853         strVBSCCable << "logicalBSCCable_2_"  << i-1;
00854         logicCable = new G4LogicalVolume(solidCable,cable,strVBSCCable.str());
00855 
00856         ostringstream strPBSCCable;
00857         strPBSCCable << "physicalBSCCable_2_"  << i-1;
00858         physiCable = new G4PVPlacement(0,
00859             G4ThreeVector((*besEMCGeometry).cablePosX[id],
00860               (*besEMCGeometry).cablePosY[id],
00861               -(*besEMCGeometry).cablePosZ[id]),
00862             logicCable,
00863             strPBSCCable.str(),
00864             logicBSCPhi,
00865             false,
00866             0);
00867         logicCable->SetVisAttributes(G4VisAttributes::Invisible);
00868 
00869       }
00870 
00871       ostringstream strCrystal;
00872       strCrystal << "physicalCrystal"  << i-1;
00873       physiBSCCrystal = new G4PVParameterised(
00874           strCrystal.str(),
00875           logicBSCCrystal,
00876           physiBSCTheta,
00877           kZAxis,
00878           1,//for this method,it must be 1.
00879           crystalParam);
00880       (*besEMCGeometry).physiBSCCrystal[i]=physiBSCCrystal;
00881       //G4cout << (*besEMCGeometry).physiBSCCrystal[i] << G4endl;
00882       physiBSCCrystal->SetCopyNo(i);
00883 
00884 
00885       if(verboseLevel>4)
00886         G4cout << "BesEmcConstruction*****************************"<< G4endl
00887           << "point of crystal =" <<physiBSCCrystal << G4endl
00888           //           << "point of mother  =" <<physiBSCCrystal->GetMotherPhysical() << G4endl
00889           << "point of excepted=" <<physiBSCTheta << G4endl;
00890       //G4Exception("BesEMCConstruction::Construct() starting............");
00891     }
00892 
00893     //
00894     //always return the physical World
00895     //
00896     if(verboseLevel>0)PrintEMCParameters();
00897     //  return physiBSC;
00898 
00899     ConstructSPFrame(logicBSCWorld,besEMCGeometry);
00900     ConstructEndGeometry(logicEMC);
00901   }
00902 
00903   //Set vis attributes and sensitive detector 
00904   SetVisAndSD();
00905 
00906   //list geo tree
00907   if(logicEMC&&physiEMC&&verboseLevel>4){
00908     G4cout<<"logicEmc "<<logicEMC<<"  physiEmc "<<physiEMC<<G4endl;
00909     G4cout<<"list geo tree"<<G4endl;
00910 
00911     int NdaughterofEMC = logicEMC->GetNoDaughters();
00912 
00913     for(int i = 0; i < NdaughterofEMC; i++)
00914     {
00915       G4LogicalVolume *daughterofEmc = logicEMC->GetDaughter(i)->GetLogicalVolume();
00916       G4cout<<i<<"/"<<NdaughterofEMC<<" name: "<<daughterofEmc->GetName()<<" "<<daughterofEmc<<" shape: "<<daughterofEmc->GetSolid()->GetName()<<G4endl;
00917       int NdaughterofEmc_2 = daughterofEmc->GetNoDaughters();
00918       for(int j = 0; j < NdaughterofEmc_2; j++)
00919       {
00920         G4LogicalVolume *daughterofEmc_2 = daughterofEmc->GetDaughter(j)->GetLogicalVolume();
00921         G4cout<<"     --> "<<j<<"/"<<NdaughterofEmc_2<<" name: "<<daughterofEmc_2->GetName()<<" "<<daughterofEmc_2<<" shape: "<<daughterofEmc_2->GetSolid()->GetName()<<G4endl;
00922         int NdaughterofEmc_3 = daughterofEmc_2->GetNoDaughters();
00923         for(int k = 0; k < NdaughterofEmc_3; k++)
00924         {
00925           G4LogicalVolume *daughterofEmc_3 = daughterofEmc_2->GetDaughter(k)->GetLogicalVolume();
00926           G4cout<<"          --> "<<k<<"/"<<NdaughterofEmc_3<<" name: "<<daughterofEmc_3->GetName()<<" "<<daughterofEmc_3<<" shape: "<<daughterofEmc_3->GetSolid()->GetName()<<G4endl;
00927           int NdaughterofEmc_4 = daughterofEmc_3->GetNoDaughters();
00928           for(int m = 0; m < NdaughterofEmc_4; m++)
00929           {
00930             G4LogicalVolume *daughterofEmc_4 = daughterofEmc_3->GetDaughter(m)->GetLogicalVolume();
00931             G4cout<<"               --> "<<m<<"/"<<NdaughterofEmc_4<<" name: "<<daughterofEmc_4->GetName()<<" "<<daughterofEmc_4<<" shape: "<<daughterofEmc_4->GetSolid()->GetName()<<G4endl; 
00932             if(daughterofEmc_3->GetSolid()->GetName().contains("solidBSCCasing"))
00933             {
00934               G4Trap *Crystal = (G4Trap *)daughterofEmc_3->GetSolid();
00935               double hz = Crystal->GetZHalfLength();
00936               double hx1 = Crystal->GetXHalfLength1();
00937               double hx2 = Crystal->GetXHalfLength2();
00938               double hx3 = Crystal->GetXHalfLength3();
00939               double hx4 = Crystal->GetXHalfLength4();
00940               double hy1 = Crystal->GetYHalfLength1();
00941               double hy2 = Crystal->GetYHalfLength2();
00942               double tanalpha1 = Crystal->GetTanAlpha1();
00943               double tanalpha2 = Crystal->GetTanAlpha2();
00944               G4cout<<"                   --> "<<hx1<<" "<<hx2<<" "<<hx3<<" "<<hx4<<" "<<hy1<<" "<<hy2<<" "<<hz<<" "<<tanalpha1<<" "<<tanalpha2<<G4endl;
00945 
00946             }//if(SolidCrystal)
00947           }//4
00948         }//3
00949       }//2
00950     }//1
00951   }
00952 
00953 
00954 }
00955 
00956 void BesEmcConstruction::ConstructEndGeometry(G4LogicalVolume* logicEMC)
00957 {
00958   G4Material* fCrystalMaterial = G4Material::GetMaterial("Cesiumiodide");
00959   G4VisAttributes* crystalVisAtt= new G4VisAttributes(G4Colour(0.5,0,1.0));
00960   crystalVisAtt->SetVisibility(false);
00961   G4VisAttributes* endPhiVisAtt= new G4VisAttributes(G4Colour(0,1.0,0));
00962   endPhiVisAtt->SetVisibility(false);
00963   const G4double zoomConst = 0.995;
00964   const G4double da=0.001*deg;
00965 
00966   //world volume of endcap
00967   //east end
00968   solidEnd = new G4Cons("solidEndWorld",(*emcEnd).WorldRmin1,(*emcEnd).WorldRmax1,
00969       (*emcEnd).WorldRmin2,(*emcEnd).WorldRmax2,
00970       (*emcEnd).WorldDz/2,0.*deg,360.*deg);
00971   logicEnd = new G4LogicalVolume(solidEnd, G4Material::GetMaterial("Aluminium"), "logicalEndWorld", 0, 0, 0);
00972   physiEnd = new G4PVPlacement(0,               // no rotation
00973       G4ThreeVector(0,0,(*emcEnd).WorldZPosition),
00974       logicEnd,               // its logical volume
00975       "physicalEndWorld0",             // its name
00976       logicEMC,               // its mother  volume
00977       false,                  // no boolean operations
00978       0);                     // no field specific to volume
00979   if(logicEnd)
00980     logicEnd->SetVisAttributes(G4VisAttributes::Invisible);
00981 
00982 
00983   //west end
00984   G4RotationMatrix *rotateEnd = new G4RotationMatrix();
00985   rotateEnd->rotateY(180.*deg);
00986   physiEnd = new G4PVPlacement(rotateEnd,
00987       G4ThreeVector(0,0,-(*emcEnd).WorldZPosition),
00988       logicEnd,
00989       "physicalEndWorld2",
00990       logicEMC,
00991       false,
00992       2);
00993 
00995   // emc endcap sectors (east)                                            //
00997   //                                    20mm gap                          // 
00998   //                                      ||                              //
00999   //                                \   7 || 6   /                        //
01000   //                           -   8 \    ||    / 5   -                   //
01001   //                             -    \   ||   /    -                     //
01002   //                        _  9   -   \  ||  /   -   4  _                //
01003   //                          - _    -  \ || /  -    _ -                  //
01004   //                              - _  - \||/ -  _ -                      //
01005   //                         10        - -||- -         3                 //
01006   //                      ----------------||----------------              //
01007   //                         11        - -||- -         2                 //
01008   //                              _ -  - /||\ -  - _                      //
01009   //                          _ -    -  / || \  -    - _                  //
01010   //                        -  12  -   /  ||  \   -   1  -                //
01011   //                             -    /   ||   \    -                     //
01012   //                           -  13 /    ||    \  0  -                   //
01013   //                                /  14 || 15  \                        //
01014   //                                      ||                              //
01016 
01017   // 1/16 of endcap world,which has some symmetry
01018   // sector 0-6,8-14
01019   solidEndPhi = new G4Cons("solidEndPhi0",
01020       (*emcEnd).SectorRmin1,(*emcEnd).SectorRmax1,(*emcEnd).SectorRmin2,(*emcEnd).SectorRmax2,
01021       (*emcEnd).SectorDz/2,0.*deg,22.5*deg-da);
01022   logicEndPhi = new G4LogicalVolume(solidEndPhi, G4Material::GetMaterial("Air"), "logicalEndPhi0", 0, 0, 0);
01023   for(G4int i=0;i<14;i++)
01024   {
01025     if((i!=6)&&(i!=7))
01026     {
01027       G4RotationMatrix *rotatePhi = new G4RotationMatrix();
01028       rotatePhi->rotateZ(-i*22.5*deg+67.5*deg);
01029       ostringstream strEndPhi;
01030       strEndPhi << "physicalEndPhi"  << i;
01031       physiEndPhi = new G4PVPlacement(rotatePhi,//0,logicEndPhi,strEndPhi.str(),logicEnd,false,i);
01032           G4ThreeVector(0,0,(*emcEnd).SectorZPosition),logicEndPhi,strEndPhi.str(),logicEnd,false,i);
01033     }
01034   }
01035   if(logicEndPhi)
01036     logicEndPhi->SetVisAttributes(endPhiVisAtt);
01037 
01038   for(G4int i=0;i<35;i++)
01039   {
01040     ostringstream strEndCasing;
01041     strEndCasing << "solidEndCasing_0_"  << i;
01042 
01043     //-************tranform to new coodinate!    liangyt 2007.5.7  *******
01044     G4ThreeVector newfPnt[8];
01045     G4ThreeVector center(0.0, 0.0, 0.0);
01046     G4ThreeVector rotAngle(0.0, 0.0, 0.0);
01047 
01048     TransformToArb8( (*emcEnd).fPnt[i], newfPnt, center, rotAngle );
01049 
01050     emcEnd->Zoom(newfPnt,zoomConst);    //change emcEnd.fPnt[i] to newfPnt
01051 
01052     G4RotationMatrix *rotatePhiIrregBox = new G4RotationMatrix();
01053     rotatePhiIrregBox->rotateX(rotAngle.x());
01054     rotatePhiIrregBox->rotateY(rotAngle.y());
01055     rotatePhiIrregBox->rotateZ(rotAngle.z());
01056     //-*******************************************************************
01057 
01058     solidEndCasing = new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);   //liangyt
01059 
01060     ostringstream strVEndCasing;
01061     strVEndCasing << "logicalEndCasing_0_"  << i;
01062     logicEndCasing = new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
01063 
01064     ostringstream strPEndCasing;
01065     strPEndCasing << "physicalEndCasing_0_"  << i;
01066     physiEndCasing = new G4PVPlacement(rotatePhiIrregBox,center,
01067         logicEndCasing,strPEndCasing.str(),logicEndPhi,false,i);  //change with rot and pos now!
01068 
01069     ostringstream strEndCrystal;
01070     strEndCrystal << "solidEndCrystal_0_"  << i;
01071 
01072     emcEnd->ModifyForCasing((*emcEnd).zoomPoint,i);
01073     solidEndCrystal = new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
01074 
01075     ostringstream strVEndCrystal;
01076     strVEndCrystal << "logicalEndCrystal_0_"  << i;
01077     logicEndCrystal = new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
01078 
01079     ostringstream strPEndCrystal;
01080     strPEndCrystal << "physicalEndCrystal_0_"  << i;
01081     physiEndCrystal = new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,false,i);
01082 
01083     logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
01084     logicEndCrystal->SetVisAttributes(crystalVisAtt);
01085     logicEndCrystal->SetSensitiveDetector(besEMCSD);
01086   }
01087 
01088   // the top area which has 20 mm gap
01089   // sector 6,14
01090   solidEndPhi = new G4Cons("solidEndPhi1",
01091       (*emcEnd).WorldRmin1,(*emcEnd).WorldRmax1,(*emcEnd).WorldRmin2,(*emcEnd).WorldRmax2,  
01092       (*emcEnd).WorldDz/2,67.5*deg,22.5*deg-da);
01093   logicEndPhi = new G4LogicalVolume(solidEndPhi, G4Material::GetMaterial("Air"), "logicalEndPhi1", 0, 0, 0);
01094   for(G4int i=0;i<2;i++)
01095   {
01096     G4RotationMatrix *rotatePhi = new G4RotationMatrix();
01097     rotatePhi->rotateZ(-i*180.*deg);
01098     ostringstream strEndPhi;
01099     strEndPhi << "physicalEndPhi"  << i*8+6;
01100     physiEndPhi = new G4PVPlacement(rotatePhi,G4ThreeVector(0,0,(*emcEnd).SectorZPosition),
01101         logicEndPhi,strEndPhi.str(),logicEnd,false,i*8+6);
01102   }
01103   if(logicEndPhi)
01104     logicEndPhi->SetVisAttributes(endPhiVisAtt);
01105 
01106   for(G4int i=0;i<35;i++)
01107   {
01108     ostringstream strEndCasing;
01109     strEndCasing << "solidEndCasing_1_"  << i;
01110 
01111     //-************tranform to new coodinate!    liangyt 2007.5.7  *******
01112     G4ThreeVector newfPnt[8];
01113     G4ThreeVector center(0.0, 0.0, 0.0);
01114     G4ThreeVector rotAngle(0.0, 0.0, 0.0);
01115 
01116     TransformToArb8( (*emcEnd).fPnt1[i], newfPnt, center, rotAngle );
01117 
01118     emcEnd->Zoom(newfPnt,zoomConst);    //change emcEnd.fPnt[i] to newfPnt
01119 
01120     G4RotationMatrix *rotatePhiIrregBox = new G4RotationMatrix();
01121     rotatePhiIrregBox->rotateX(rotAngle.x());
01122     rotatePhiIrregBox->rotateY(rotAngle.y());
01123     rotatePhiIrregBox->rotateZ(rotAngle.z());
01124     //-*******************************************************************
01125 
01126     solidEndCasing = new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);
01127 
01128     ostringstream strVEndCasing;
01129     strVEndCasing << "logicalEndCasing_1_"  << i;
01130     logicEndCasing = new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
01131 
01132     ostringstream strPEndCasing;
01133     strPEndCasing << "physicalEndCasing_1_"  << i;
01134     physiEndCasing = new G4PVPlacement(rotatePhiIrregBox,center,
01135         logicEndCasing,strPEndCasing.str(),logicEndPhi,false,i);  //change with rot and pos now!
01136 
01137     ostringstream strEndCrystal;
01138     strEndCrystal << "solidEndCrystal_1_"  << i;
01139 
01140     emcEnd->ModifyForCasing((*emcEnd).zoomPoint,i);
01141     solidEndCrystal = new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
01142 
01143     ostringstream strVEndCrystal;
01144     strVEndCrystal << "logicalEndCrystal_1_"  << i;
01145     logicEndCrystal = new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
01146 
01147     ostringstream strPEndCrystal;
01148     strPEndCrystal << "physicalEndCrystal_1_"  << i;
01149     physiEndCrystal = new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,false,i);
01150 
01151     logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
01152     logicEndCrystal->SetVisAttributes(crystalVisAtt);
01153     logicEndCrystal->SetSensitiveDetector(besEMCSD);
01154   }
01155 
01156   (*emcEnd).ReflectX();
01157 
01158   // sector 7,15
01159   for(G4int i=0;i<35;i++)
01160     for (G4int j=0;j<8;j++)
01161       (*emcEnd).fPnt1[i][j].rotateZ(-90.*deg);
01162 
01163   solidEndPhi = new G4Cons("solidEndPhi2",
01164       (*emcEnd).WorldRmin1,(*emcEnd).WorldRmax1,(*emcEnd).WorldRmin2,(*emcEnd).WorldRmax2,
01165       (*emcEnd).WorldDz/2,0*deg,22.5*deg-da);
01166   logicEndPhi = new G4LogicalVolume(solidEndPhi, G4Material::GetMaterial("Air"), "logicalEndPhi2", 0, 0, 0);
01167   for(G4int i=0;i<2;i++)
01168   {
01169     G4RotationMatrix *rotatePhi = new G4RotationMatrix();
01170     rotatePhi->rotateZ(-i*180.*deg-90.*deg);
01171     ostringstream strEndPhi;
01172     strEndPhi << "physicalEndPhi"  << i*8+7;
01173     physiEndPhi = new G4PVPlacement(rotatePhi,G4ThreeVector(0,0,(*emcEnd).SectorZPosition),
01174         logicEndPhi,strEndPhi.str(),logicEnd,false,i*8+7);
01175   }
01176   if(logicEndPhi)
01177     logicEndPhi->SetVisAttributes(endPhiVisAtt);
01178 
01179   for(G4int i=0;i<35;i++)
01180   {
01181     ostringstream strEndCasing;
01182     strEndCasing << "solidEndCasing_2_"  << i;
01183 
01184     //-************tranform to new coodinate!    liangyt 2007.5.7  *******
01185     G4ThreeVector newfPnt[8];
01186     G4ThreeVector center(0.0, 0.0, 0.0);
01187     G4ThreeVector rotAngle(0.0, 0.0, 0.0);
01188 
01189     TransformToArb8( (*emcEnd).fPnt1[i], newfPnt, center, rotAngle );
01190 
01191     emcEnd->Zoom(newfPnt,zoomConst);    //change emcEnd.fPnt[i] to newfPnt
01192 
01193     G4RotationMatrix *rotatePhiIrregBox = new G4RotationMatrix();
01194     rotatePhiIrregBox->rotateX(rotAngle.x());
01195     rotatePhiIrregBox->rotateY(rotAngle.y());
01196     rotatePhiIrregBox->rotateZ(rotAngle.z());
01197     //-*******************************************************************
01198 
01199     solidEndCasing = new G4IrregBox(strEndCasing.str(),(*emcEnd).zoomPoint);
01200 
01201     ostringstream strVEndCasing;
01202     strVEndCasing << "logicalEndCasing_2_"  << i;
01203     logicEndCasing = new G4LogicalVolume(solidEndCasing,fCasingMaterial,strVEndCasing.str());
01204 
01205     ostringstream strPEndCasing;
01206     strPEndCasing << "physicalEndCasing_2_"  << i;
01207     physiEndCasing = new G4PVPlacement(rotatePhiIrregBox,center,
01208         logicEndCasing,strPEndCasing.str(),logicEndPhi,false,i);  //change with rot and pos now!
01209 
01210     ostringstream strEndCrystal;
01211     strEndCrystal << "solidEndCrystal_2_"  << i;
01212 
01213     emcEnd->ModifyForCasing((*emcEnd).zoomPoint,i);
01214     solidEndCrystal = new G4IrregBox(strEndCrystal.str(),(*emcEnd).cryPoint);
01215 
01216     ostringstream strVEndCrystal;
01217     strVEndCrystal << "logicalEndCrystal_2_"  << i;
01218     logicEndCrystal = new G4LogicalVolume(solidEndCrystal,fCrystalMaterial,strVEndCrystal.str());
01219 
01220     ostringstream strPEndCrystal;
01221     strPEndCrystal << "physicalEndCrystal_2_"  << i;
01222     physiEndCrystal = new G4PVPlacement(0,0,logicEndCrystal,strPEndCrystal.str(),logicEndCasing,false,i);
01223 
01224     logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
01225     logicEndCrystal->SetVisAttributes(crystalVisAtt);
01226     logicEndCrystal->SetSensitiveDetector(besEMCSD);
01227   }
01228 }
01229 
01230 G4int BesEmcConstruction::ComputeEndCopyNb(G4int num)
01231 {
01232   G4int copyNb;
01233   switch(num){
01234     case 30:
01235       copyNb = 5;
01236       break;
01237     case 31:
01238       copyNb = 6;
01239       break;
01240     case 32:
01241       copyNb = 14;
01242       break;
01243     case 33:
01244       copyNb = 15;
01245       break;
01246     case 34:
01247       copyNb = 16;
01248       break;
01249     default:
01250       copyNb = num;
01251       break;
01252   }
01253   return copyNb;
01254 }
01255 
01256 void BesEmcConstruction::ConstructSPFrame(G4LogicalVolume* logicEMC, BesEmcGeometry *besEMCGeometry)
01257 {
01258   G4double rmax=(*besEMCGeometry).BSCRmax+2.*mm;  //radius from 942mm to 940mm
01259   solidSupportBar = new G4Tubs("solidSupportBar0",
01260       rmax+(*besEMCGeometry).SPBarThickness1,
01261       rmax+(*besEMCGeometry).SPBarThickness+(*besEMCGeometry).SPBarThickness1,
01262       (*besEMCGeometry).BSCDz
01263       +(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz,
01264       0.*deg,
01265       360.*deg);
01266 
01267   logicSupportBar = new G4LogicalVolume(solidSupportBar,stainlessSteel,"logicalSupportBar0");
01268 
01269   physiSupportBar = new G4PVPlacement(0,0,logicSupportBar,"physicalSupportBar0",logicEMC,false,0);
01270 
01271   solidSupportBar1 = new G4Tubs("solidSupportBar1",
01272       rmax,
01273       rmax+(*besEMCGeometry).SPBarThickness1,
01274       (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3,
01275       (*besEMCGeometry).BSCPhiDphi-(*besEMCGeometry).SPBarDphi/2,
01276       (*besEMCGeometry).SPBarDphi);
01277 
01278   logicSupportBar1 = new G4LogicalVolume(solidSupportBar1,stainlessSteel,"logicalSupportBar1");
01279 
01280   for(G4int i=0;i<(*besEMCGeometry).BSCNbPhi/2;i++)
01281   {
01282     G4RotationMatrix *rotateSPBar = new G4RotationMatrix();
01283     rotateSPBar->rotateZ((*besEMCGeometry).BSCPhiDphi-i*2*(*besEMCGeometry).BSCPhiDphi);
01284     ostringstream strSupportBar1;
01285     strSupportBar1 << "physicalSupportBar1_"  << i;
01286     physiSupportBar1 = new G4PVPlacement(rotateSPBar,0,
01287         logicSupportBar1,strSupportBar1.str(),logicEMC,false,0);
01288   }
01289 
01290   //end ring
01291   solidEndRing = new G4Tubs("solidEndRing",
01292       (*besEMCGeometry).EndRingRmin,
01293       (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr/2,
01294       (*besEMCGeometry).EndRingDz/2,
01295       0.*deg,
01296       360.*deg);
01297 
01298   solidGear = new G4Tubs("solidGear",
01299       (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr/2,
01300       (*besEMCGeometry).EndRingRmin+(*besEMCGeometry).EndRingDr,
01301       (*besEMCGeometry).EndRingDz/2,
01302       0.*deg,
01303       (*besEMCGeometry).BSCPhiDphi);
01304 
01305   //taper ring
01306   solidTaperRing1 = new G4Tubs("solidTaperRing1",
01307       (*besEMCGeometry).TaperRingRmin1,
01308       (*besEMCGeometry).TaperRingRmin1+(*besEMCGeometry).TaperRingThickness1,
01309       (*besEMCGeometry).TaperRingInnerLength/2,
01310       0.*deg,
01311       360.*deg);
01312 
01313   solidTaperRing2 = new G4Cons("solidTaperRing2",
01314       (*besEMCGeometry).TaperRingRmin1,
01315       (*besEMCGeometry).TaperRingRmin1+(*besEMCGeometry).TaperRingDr,
01316       (*besEMCGeometry).TaperRingRmin2,
01317       (*besEMCGeometry).TaperRingRmin2+(*besEMCGeometry).TaperRingDr,
01318       (*besEMCGeometry).TaperRingDz/2,
01319       0.*deg,
01320       360.*deg);
01321 
01322   solidTaperRing3 = new G4Cons("solidTaperRing3",
01323       (*besEMCGeometry).BSCRmax2,
01324       (*besEMCGeometry).BSCRmax2+(*besEMCGeometry).TaperRingOuterLength1,
01325       (*besEMCGeometry).TaperRingRmin2+(*besEMCGeometry).TaperRingDr,
01326       (*besEMCGeometry).TaperRingRmin2+(*besEMCGeometry).TaperRingDr+(*besEMCGeometry).TaperRingOuterLength,
01327       (*besEMCGeometry).TaperRingThickness3/2,
01328       0.*deg,
01329       360.*deg);
01330 
01331   logicEndRing = new G4LogicalVolume(solidEndRing,stainlessSteel,"logicalEndRing");
01332   logicGear = new G4LogicalVolume(solidGear,stainlessSteel,"logicalGear");
01333   logicTaperRing1 = new G4LogicalVolume(solidTaperRing1,stainlessSteel,"logicalTaperRing1");
01334   logicTaperRing2 = new G4LogicalVolume(solidTaperRing2,stainlessSteel,"logicalTaperRing2");
01335   logicTaperRing3 = new G4LogicalVolume(solidTaperRing3,stainlessSteel,"logicalTaperRing3");
01336 
01337   for(G4int i=0;i<2;i++)
01338   {
01339     G4RotationMatrix *rotateSPRing = new G4RotationMatrix();
01340     G4double zEndRing,z1,z2,z3;
01341     if(i==0)
01342     {
01343       zEndRing = (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz/2;
01344       z1 = (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3
01345         -(*besEMCGeometry).TaperRingDz-(*besEMCGeometry).TaperRingInnerLength/2;
01346       z2 = (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3-(*besEMCGeometry).TaperRingDz/2;
01347       z3 = (*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3/2;
01348     }
01349     else
01350     {
01351       rotateSPRing->rotateY(180.*deg);
01352       zEndRing = -((*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3+(*besEMCGeometry).EndRingDz/2);
01353       z1 = -((*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3
01354           -(*besEMCGeometry).TaperRingDz-(*besEMCGeometry).TaperRingInnerLength/2);
01355       z2 = -((*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3-(*besEMCGeometry).TaperRingDz/2);
01356       z3 = -((*besEMCGeometry).BSCDz+(*besEMCGeometry).TaperRingThickness3/2);
01357     }
01358 
01359     ostringstream strEndRing;
01360     strEndRing << "physicalEndRing_"  << i;
01361     physiEndRing = new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,zEndRing),
01362         logicEndRing,strEndRing.str(),logicEMC,false,0);
01363 
01364     for(G4int j=0;j<(*besEMCGeometry).BSCNbPhi/2;j++)
01365     {
01366       G4RotationMatrix *rotateGear = new G4RotationMatrix();
01367       rotateGear->rotateZ((*besEMCGeometry).BSCPhiDphi/2-j*2*(*besEMCGeometry).BSCPhiDphi);
01368 
01369       ostringstream strGear;
01370       strGear << "physicalGear_"  << i << "_" <<j;
01371       physiGear = new G4PVPlacement(rotateGear,G4ThreeVector(0,0,zEndRing),
01372           logicGear,strGear.str(),logicEMC,false,0);
01373     }
01374 
01375     ostringstream strTaperRing1;
01376     strTaperRing1 << "physicalTaperRing1_"  << i;
01377     physiTaperRing1 = new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z1),
01378         logicTaperRing1,strTaperRing1.str(),logicEMC,false,0);
01379 
01380     ostringstream strTaperRing2;
01381     strTaperRing2 << "physicalTaperRing2_"  << i;
01382     physiTaperRing2 = new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z2),
01383         logicTaperRing2,strTaperRing2.str(),logicEMC,false,0);
01384 
01385     ostringstream strTaperRing3;
01386     strTaperRing3 << "physicalTaperRing3_"  << i;
01387     physiTaperRing3 = new G4PVPlacement(rotateSPRing,G4ThreeVector(0,0,z3),
01388         logicTaperRing3,strTaperRing3.str(),logicEMC,false,0);
01389   }
01390 }
01391 
01392 //Set vis attributes and sensitive detector
01393 void BesEmcConstruction::SetVisAndSD()
01394 {
01395   //-------------------------------------------------------------
01396   //Barrel
01397   G4VisAttributes* bscVisAtt= new G4VisAttributes(G4Colour(0.5,0.5,0.5));
01398   bscVisAtt->SetVisibility(false);
01399   logicEMC->SetVisAttributes(bscVisAtt);
01400   if(logicBSCWorld)
01401     logicBSCWorld->SetVisAttributes(G4VisAttributes::Invisible);
01402 
01403   if (logicBSCCrystal) {
01404     //G4cout<<"find BSCCrystal "<<G4endl;
01405     G4VisAttributes* crystalVisAtt= new G4VisAttributes(G4Colour(0,0,1.0));
01406     crystalVisAtt->SetVisibility(true);
01407     logicBSCCrystal->SetVisAttributes(crystalVisAtt);
01408     logicBSCCrystal->SetSensitiveDetector(besEMCSD);
01409   }
01410 
01411   if(logicBSCPhi) {
01412     G4VisAttributes* rightVisAtt= new G4VisAttributes(G4Colour(1.0,0.,1.0));
01413     rightVisAtt->SetVisibility(false);
01414     logicBSCPhi->SetVisAttributes(rightVisAtt);
01415   }
01416 
01417   if(logicRear)
01418     logicRear->SetVisAttributes(G4VisAttributes::Invisible);
01419   if(logicOrgGlass)
01420     logicOrgGlass->SetVisAttributes(G4VisAttributes::Invisible);
01421   if(logicRearCasing)
01422     logicRearCasing->SetVisAttributes(G4VisAttributes::Invisible);
01423   if(logicAlPlate) 
01424     logicAlPlate->SetVisAttributes(G4VisAttributes::Invisible);
01425   if(logicPD) {
01426     logicPD->SetVisAttributes(G4VisAttributes::Invisible);
01427     logicPD->SetSensitiveDetector(besEMCSD);
01428   }
01429   if(logicPreAmpBox)
01430     logicPreAmpBox->SetVisAttributes(G4VisAttributes::Invisible);
01431   if(logicAirInPABox)
01432     logicAirInPABox->SetVisAttributes(G4VisAttributes::Invisible);
01433   if(logicHangingPlate)
01434     logicHangingPlate->SetVisAttributes(G4VisAttributes::Invisible);
01435   if(logicWaterPipe)
01436     logicWaterPipe->SetVisAttributes(G4VisAttributes::Invisible);
01437 
01438   //-------------------------------------------------------------
01439   //Support system
01440   G4VisAttributes* ringVisAtt= new G4VisAttributes(G4Colour(0.5,0.25,0.));
01441   ringVisAtt->SetVisibility(false);
01442   if(logicSupportBar)
01443     logicSupportBar->SetVisAttributes(ringVisAtt);
01444   if(logicSupportBar1)
01445     logicSupportBar1->SetVisAttributes(ringVisAtt);
01446   if(logicEndRing)
01447     logicEndRing->SetVisAttributes(ringVisAtt);
01448   if(logicGear)
01449     logicGear->SetVisAttributes(ringVisAtt);
01450   if(logicTaperRing1)
01451     logicTaperRing1->SetVisAttributes(ringVisAtt);
01452   if(logicTaperRing2)
01453     logicTaperRing2->SetVisAttributes(ringVisAtt);
01454   if(logicTaperRing3)
01455     logicTaperRing3->SetVisAttributes(ringVisAtt);
01456 
01457   //-------------------------------------------------------------
01458   //Endcap
01459   G4VisAttributes* endPhiVisAtt= new G4VisAttributes(G4Colour(0,1.0,0));
01460   endPhiVisAtt->SetVisibility(false);
01461   if(logicEnd)
01462     logicEnd->SetVisAttributes(endPhiVisAtt);
01463 }
01464 
01465 //Get logical volumes from gdml
01466 void BesEmcConstruction::GetLogicalVolume()
01467 {
01468   //-------------------------------------------------------------
01469   //Barrel
01470   logicBSCWorld = FindLogicalVolume("logicalBSCWorld"); 
01471   logicBSCCrystal = FindLogicalVolume("logicalCrystal");
01472   logicBSCPhi = FindLogicalVolume("logicalBSCPhi");
01473   logicRear = FindLogicalVolume("logicalRearBox");
01474   logicOrgGlass = FindLogicalVolume("logicalOrganicGlass");
01475   logicRearCasing = FindLogicalVolume("logicalRearCasing");
01476   logicAlPlate = FindLogicalVolume("logicalAlPlate");
01477   logicPD = FindLogicalVolume("logicalPD");
01478   logicPreAmpBox = FindLogicalVolume("logicalPreAmpBox");
01479   logicAirInPABox = FindLogicalVolume("logicalAirInPABox");
01480   logicHangingPlate = FindLogicalVolume("logicalHangingPlate");
01481   logicWaterPipe = FindLogicalVolume("logicalWaterPipe");
01482 
01483   for(int i = 0; i < 44; i++){
01484     std::ostringstream osnameBSCCasing;
01485     osnameBSCCasing << "logicalBSCCasing"<<i;
01486     logicBSCTheta = FindLogicalVolume( osnameBSCCasing.str() );
01487     if(logicBSCTheta)
01488     {
01489       G4VisAttributes* rightVisAtt= new G4VisAttributes(G4Colour(1.0, 0.0,0.0));
01490       rightVisAtt->SetVisibility(false);
01491       logicBSCTheta->SetVisAttributes(rightVisAtt);
01492     }
01493 
01494     std::ostringstream osnameBSCCable1;
01495     osnameBSCCable1 << "logicalBSCCable_1_"<<i;
01496     logicCable = FindLogicalVolume( osnameBSCCable1.str() );
01497     if(logicCable)
01498       logicCable->SetVisAttributes(G4VisAttributes::Invisible);
01499 
01500     std::ostringstream osnameBSCCable2;
01501     osnameBSCCable2 << "logicalBSCCable_2_"<<i;
01502     logicCable = FindLogicalVolume( osnameBSCCable2.str() );
01503     if(logicCable)
01504       logicCable->SetVisAttributes(G4VisAttributes::Invisible);
01505 
01506     std::ostringstream osnameOCGirder1;
01507     osnameOCGirder1 <<"logicalOpenningCutGirder_1_"<<i;
01508     logicOCGirder = FindLogicalVolume( osnameOCGirder1.str() );
01509     if(logicOCGirder)
01510       logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
01511 
01512     std::ostringstream osnameOCGirder2;
01513     osnameOCGirder2 <<"logicalOpenningCutGirder_2_"<<i;
01514     logicOCGirder = FindLogicalVolume( osnameOCGirder2.str() );
01515     if(logicOCGirder)
01516       logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
01517 
01518     std::ostringstream osnameOCGirder3;
01519     osnameOCGirder3 <<"logicalOpenningCutGirder_3_"<<i;
01520     logicOCGirder = FindLogicalVolume( osnameOCGirder3.str() );
01521     if(logicOCGirder)
01522       logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
01523 
01524     std::ostringstream osnameOCGirder4;
01525     osnameOCGirder4 <<"logicalOpenningCutGirder_4_"<<i;
01526     logicOCGirder = FindLogicalVolume( osnameOCGirder4.str() );
01527     if(logicOCGirder)
01528       logicOCGirder->SetVisAttributes(G4VisAttributes::Invisible);
01529   }
01530 
01531   //-------------------------------------------------------------
01532   //Support system
01533   logicSupportBar = FindLogicalVolume("logicalSupportBar0");
01534   logicSupportBar1 = FindLogicalVolume("logicalSupportBar1");
01535   logicEndRing = FindLogicalVolume("logicalEndRing");
01536   logicGear = FindLogicalVolume("logicalGear");
01537   logicTaperRing1 = FindLogicalVolume("logicalTaperRing1");
01538   logicTaperRing2 = FindLogicalVolume("logicalTaperRing2");
01539   logicTaperRing3 = FindLogicalVolume("logicalTaperRing3");
01540 
01541   //-------------------------------------------------------------
01542   //Endcap
01543   logicEnd = FindLogicalVolume("logicalEndWorld");
01544 
01545   for(G4int sector=0;sector<3;sector++) {
01546     std::ostringstream osnameEndPhi;
01547     osnameEndPhi<<"logicalEndPhi"<<sector;
01548     logicEndPhi = FindLogicalVolume(osnameEndPhi.str());
01549     if(logicEndPhi) {
01550       logicEndPhi->SetVisAttributes(G4VisAttributes::Invisible);
01551     } else {
01552       G4cout<<"Can't find logicEndPhi!"<<G4endl;
01553     }
01554 
01555     for(G4int cryNb=0;cryNb<35;cryNb++) {
01556 
01557       std::ostringstream osnameEndCrystal;
01558       osnameEndCrystal<<"logicalEndCrystal_"<<sector<<"_"<<cryNb;
01559       logicEndCrystal = FindLogicalVolume( osnameEndCrystal.str() );
01560       if(logicEndCrystal) {
01561         logicEndCrystal->SetSensitiveDetector(besEMCSD);
01562         G4VisAttributes* crystalVisAtt
01563           = new G4VisAttributes(G4Colour(0.5,0,1.0));
01564         crystalVisAtt->SetVisibility(false);
01565         logicEndCrystal->SetVisAttributes(crystalVisAtt);
01566       } else {
01567         G4cout<<"Can't find: "<<osnameEndCrystal.str()<<G4endl;
01568       }
01569 
01570       std::ostringstream osnameEndCasing;
01571       osnameEndCasing<<"logicalEndCasing_"<<sector<<"_"<<cryNb;
01572       logicEndCasing = FindLogicalVolume( osnameEndCasing.str() );
01573       if(logicEndCasing) {
01574         logicEndCasing->SetVisAttributes(G4VisAttributes::Invisible);
01575       } else {
01576         G4cout<<"Can't find: "<<osnameEndCasing.str()<<G4endl;
01577       }
01578     }
01579   }
01580 }
01581 
01582 void BesEmcConstruction::DefineMaterials()
01583 {
01584   G4String name, symbol;             //a=mass of a mole;
01585   G4double a, z, density;            //z=mean number of protons;  
01586   //  G4int iz, n;                       //iz=number of protons  in an isotope;
01587   // n=number of nucleons in an isotope;
01588 
01589   G4int ncomponents, natoms;
01590   G4double fractionmass;
01591   //G4double abundance, fractionmass;
01592   //  G4double temperature, pressure;
01593 
01594   //for debug
01595   //  G4Exception("BesEmcConstruction::DefineMaterials() starting...........");
01596   //
01597   // define Elements
01598   //
01599   G4Element* H=G4Element::GetElement("Hydrogen");
01600   if(!H)
01601   {
01602     a = 1.01*g/mole;
01603     H = new G4Element(name="Hydrogen",symbol="H" , z= 1., a);
01604   }
01605   G4Element* C=G4Element::GetElement("Carbon");
01606   if(!C)
01607   {
01608     a = 12.01*g/mole;
01609     C = new G4Element(name="Carbon"  ,symbol="C" , z= 6., a);
01610   }
01611   G4Element* O=G4Element::GetElement("Oxygen");
01612   if(!O)
01613   {
01614     a = 16.00*g/mole;
01615     O = new G4Element(name="Oxygen"  ,symbol="O" , z= 8., a);  
01616   }
01617 
01618   density = 0.344*g/cm3;
01619   G4Material* Tyvek = new G4Material(name="M_Polyethylene", density, ncomponents=2);
01620   Tyvek->AddElement(C, natoms=1);
01621   Tyvek->AddElement(H, natoms=2);
01622 
01623   density = 1.39*g/cm3;
01624   G4Material* Mylar = new G4Material(name="M_PolyethyleneTerephthlate", density, ncomponents=3);
01625   Mylar->AddElement(C, natoms=5);
01626   Mylar->AddElement(H, natoms=4);
01627   Mylar->AddElement(O, natoms=2);
01628 
01629   density = 1.18*g/cm3;
01630   organicGlass = new G4Material(name="M_OrganicGlass", density, ncomponents=3);
01631   organicGlass->AddElement(C, natoms=5);
01632   organicGlass->AddElement(H, natoms=7);
01633   organicGlass->AddElement(O, natoms=2);
01634 
01635   G4Material *Fe = new G4Material(name="M_Iron", z=26., a=55.85*g/mole, density=7.87*g/cm3);
01636   G4Material *Cr = new G4Material(name="M_Chromium", z=24., a=52.00*g/mole, density=8.72*g/cm3);
01637   G4Material *Ni = new G4Material(name="M_Nickel", z=28., a=58.69*g/mole, density=8.72*g/cm3);
01638 
01639   stainlessSteel = new G4Material(name="M_Cr18Ni9", density=7.85*g/cm3, ncomponents=3);
01640   stainlessSteel->AddMaterial(Fe, fractionmass=73.*perCent);
01641   stainlessSteel->AddMaterial(Cr, fractionmass=18.*perCent);
01642   stainlessSteel->AddMaterial(Ni, fractionmass=9.*perCent);
01643 
01644   G4Material *H2O = G4Material::GetMaterial("Water");
01645   G4Material *Cu = G4Material::GetMaterial("Copper");
01646   G4double dWater = 1.*g/cm3; //density
01647   G4double dCopper = 8.96*g/cm3;
01648   G4double aWater = ((*besEMCGeometry).waterPipeDr-(*besEMCGeometry).waterPipeThickness)
01649     *((*besEMCGeometry).waterPipeDr-(*besEMCGeometry).waterPipeThickness);    //area
01650   G4double aCopper = (*besEMCGeometry).waterPipeDr*(*besEMCGeometry).waterPipeDr-aWater;
01651   density = (dWater*aWater+dCopper*aCopper)/(aWater+aCopper);
01652 
01653   waterPipe = new G4Material(name="M_WaterPipe", density, ncomponents=2);
01654   fractionmass = dWater*aWater/(dWater*aWater+dCopper*aCopper);
01655   waterPipe->AddMaterial(H2O, fractionmass);
01656   fractionmass = dCopper*aCopper/(dWater*aWater+dCopper*aCopper);
01657   waterPipe->AddMaterial(Cu, fractionmass);
01658 
01659   cable  = new G4Material(name="M_Cable", density=4.*g/cm3, ncomponents=1);
01660   cable->AddMaterial(Cu,1);
01661 
01662   //for debug
01663   //G4Exception("BesEmcConstruction::DefineMaterials() running one.........");
01664   //
01665   // predigest the casing of crystals to a mixture
01666   //
01667 
01668   G4Material* Al=G4Material::GetMaterial("Aluminium");
01669   if(Al==NULL)
01670   {
01671     Al = new G4Material(name="Aluminium", z=13., a=26.98*g/mole, density=2.700*g/cm3);
01672   }
01673 
01674   G4Material *Si=G4Material::GetMaterial("M_Silicon");
01675   if(Si==NULL)
01676   {
01677     Si = new G4Material(name="M_Silicon", z=14., a=28.0855*g/mole, density=2.33*g/cm3);
01678   }
01679 
01680 
01681   //for debug
01682   G4double totalThickness=(*besEMCGeometry).fTyvekThickness
01683     +(*besEMCGeometry).fAlThickness+(*besEMCGeometry).fMylarThickness;
01684   density = (Tyvek->GetDensity()*(*besEMCGeometry).fTyvekThickness+
01685       Al->GetDensity()*(*besEMCGeometry).fAlThickness+
01686       Mylar->GetDensity()*(*besEMCGeometry).fMylarThickness)
01687     /totalThickness;
01688   G4Material* Casing = new G4Material(name="M_Casing", density, ncomponents=3);
01689   Casing->AddMaterial(
01690       Tyvek,
01691       fractionmass=Tyvek->GetDensity()/density
01692       *(*besEMCGeometry).fTyvekThickness
01693       /totalThickness);
01694   Casing->AddMaterial(
01695       Al,
01696       fractionmass=Al->GetDensity()/density
01697       *(*besEMCGeometry).fAlThickness
01698       /totalThickness);
01699   Casing->AddMaterial(
01700       Mylar,
01701       fractionmass=Mylar->GetDensity()/density
01702       *(*besEMCGeometry).fMylarThickness
01703       /totalThickness);
01704   fCasingMaterial = Casing;
01705   rearCasingMaterial = Tyvek;
01706   //for debug
01707   //  G4Exception("BesEmcConstruction::DefineMaterials() running two........");
01708   fCrystalMaterial = G4Material::GetMaterial("Cesiumiodide");
01709 
01710 }
01711 
01712 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01713 void BesEmcConstruction::PrintEMCParameters()
01714 {
01715   G4cout << "-------------------------------------------------------"<< G4endl
01716     << "---> There are "
01717     << phiNbCrystals << "(max=" << (*besEMCGeometry).BSCNbPhi
01718     << ") crystals along phi direction and "
01719     << thetaNbCrystals << "(max=" << (*besEMCGeometry).BSCNbTheta
01720     << ") crystals along theta direction."<< G4endl
01721     << "The crystals have sizes of "
01722     << (*besEMCGeometry).BSCCryLength/cm << "cm(L) and "
01723     << (*besEMCGeometry).BSCYFront/cm << "cm(Y) with " 
01724     << fCrystalMaterial->GetName() <<"."<< G4endl
01725     << "The casing is layer of "
01726     << (*besEMCGeometry).fTyvekThickness/mm << "mm tyvek,"
01727     << (*besEMCGeometry).fAlThickness/mm << "mm aluminum and"
01728     << (*besEMCGeometry).fMylarThickness/mm << "mm mylar."<< G4endl
01729     << "-------------------------------------------------------"<< G4endl;
01730   G4cout << G4Material::GetMaterial("PolyethyleneTerephthlate") << G4endl
01731     << G4Material::GetMaterial("Casing") << G4endl
01732     << G4Material::GetMaterial("Polyethylene") << G4endl
01733     << "-------------------------------------------------------"<< G4endl;
01734 }
01735 
01736 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01737 
01738 void BesEmcConstruction::SetCrystalMaterial(G4String materialChoice)
01739 {
01740   // search the material by its name   
01741   G4Material* pttoMaterial = G4Material::GetMaterial(materialChoice);     
01742   if (pttoMaterial)
01743   {fCrystalMaterial = pttoMaterial;
01744     logicBSCCrystal->SetMaterial(pttoMaterial); 
01745     PrintEMCParameters();
01746   }             
01747 }
01748 
01749 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01750 
01751 void BesEmcConstruction::SetCasingMaterial(G4String materialChoice)
01752 {
01753   // search the material by its name 
01754   G4Material* pttoMaterial = G4Material::GetMaterial(materialChoice);  
01755   if (pttoMaterial)
01756   {fCasingMaterial = pttoMaterial;
01757     logicBSCTheta->SetMaterial(pttoMaterial);
01758     PrintEMCParameters();
01759   }             
01760 }
01761 
01762 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01763 
01764 void BesEmcConstruction::SetCasingThickness(G4ThreeVector val)
01765 {
01766   // change Gap thickness and recompute the calorimeter parameters
01767   (*besEMCGeometry).fTyvekThickness = val('X');
01768   (*besEMCGeometry).fAlThickness    = val('Y');
01769   (*besEMCGeometry).fMylarThickness = val('Z');
01770 }  
01771 
01772 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01773 
01774 void BesEmcConstruction::SetBSCRmin(G4double val)
01775 {
01776   (*besEMCGeometry).BSCRmin = val;
01777 }
01778 
01779 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01780 
01781 void BesEmcConstruction::SetBSCNbPhi(G4int val)
01782 {
01783   (*besEMCGeometry).BSCNbPhi = val;
01784 }
01785 
01786 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01787 
01788 void BesEmcConstruction::SetBSCNbTheta(G4int val)
01789 {
01790   (*besEMCGeometry).BSCNbTheta = val;
01791 }
01792 
01793 void BesEmcConstruction::SetStartIDTheta(G4int val)
01794 {
01795   startID = val;
01796 }
01797 
01798 
01799 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01800 
01801 void BesEmcConstruction::SetBSCCrystalLength(G4double val)
01802 {
01803   (*besEMCGeometry).BSCCryLength = val;
01804 }
01805 
01806 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01807 
01808 void BesEmcConstruction::SetBSCYFront0(G4double val)
01809 {
01810   (*besEMCGeometry).BSCYFront0 = val;
01811 }
01812 
01813 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01814 
01815 void BesEmcConstruction::SetBSCYFront(G4double val)
01816 {
01817   (*besEMCGeometry).BSCYFront = val;
01818 }
01819 
01820 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01821 
01822 void BesEmcConstruction::SetBSCPosition0(G4double val)
01823 {
01824   (*besEMCGeometry).BSCPosition0 = val;
01825 }
01826 
01827 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01828 
01829 void BesEmcConstruction::SetBSCPosition1(G4double val)
01830 {
01831   (*besEMCGeometry).BSCPosition1 = val;
01832 }
01833 
01834 
01835 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01836 
01837 void BesEmcConstruction::SetMagField(G4double fieldValue)
01838 {
01839   //apply a global uniform magnetic field along Z axis
01840   G4FieldManager* fieldMgr 
01841     = G4TransportationManager::GetTransportationManager()->GetFieldManager();
01842 
01843   if(magField) delete magField;         //delete the existing magn field
01844 
01845   if(fieldValue!=0.)                    // create a new one if non nul
01846   { magField = new G4UniformMagField(G4ThreeVector(0.,0.,fieldValue));        
01847     fieldMgr->SetDetectorField(magField);
01848     fieldMgr->CreateChordFinder(magField);
01849     fmagField=fieldValue;
01850   } else {
01851     magField = 0;
01852     fieldMgr->SetDetectorField(magField);
01853     fmagField=0.;
01854   }
01855 }
01856 
01857 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01858 //??????????????????????????????????????????????? 
01859 void BesEmcConstruction::UpdateGeometry()
01860 {
01861   ;//G4RunManager::GetRunManager()->DefineWorldVolume(BesDetectorConstruction::Construct());
01862 }
01863 
01864 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
01865 
01866   void
01867 BesEmcConstruction::TransformToArb8( const G4ThreeVector fPnt[8], G4ThreeVector newfPnt[8], G4ThreeVector &center, G4ThreeVector &rotAngle )
01868 {
01869   HepPoint3D point[8];
01870   center = G4ThreeVector(0.0, 0.0, 0.0);
01871   for (int i = 0; i < 8; i++) {
01872     point[i] = HepPoint3D( fPnt[i].x(), fPnt[i].y(), fPnt[i].z() );
01873     center += point[i];
01874   }
01875   center /= 8.0;
01876 
01877   HepPlane3D bottomPlane( point[4], point[5], point[6] );
01878   HepPoint3D centerProject = bottomPlane.point( center );
01879   Hep3Vector newZ = center - centerProject;
01880 
01881   rotAngle = RotAngleFromNewZ( newZ );
01882   G4RotationMatrix *g4Rot = new G4RotationMatrix();
01883   g4Rot->rotateX( rotAngle.x() );
01884   g4Rot->rotateY( rotAngle.y() );
01885 
01886   G4AffineTransform *transform = new G4AffineTransform( g4Rot, center );
01887   transform->Invert();
01888   for (int i = 0; i < 8; i++) {
01889     newfPnt[i] = transform->TransformPoint(fPnt[i]);
01890   }
01891   delete g4Rot;
01892   delete transform;
01893 }
01894 
01895 
01896   void 
01897 BesEmcConstruction::ThreeVectorTrans(G4ThreeVector fPnt[8], double x[8], double y[8], double z[8])
01898 {
01899   for (int i = 0; i < 8; i++) {
01900     x[i] = fPnt[i].x();
01901     y[i] = fPnt[i].y();
01902     z[i] = fPnt[i].z();
01903   }
01904 }
01905 
01906 
01907   Hep3Vector
01908 BesEmcConstruction::RotAngleFromNewZ( Hep3Vector newZ )
01909 {
01910   newZ.setMag(1.0);
01911   Hep3Vector x0(1, 0, 0), y0(0, 1, 0), z0(0, 0, 1);
01912   double dx, dy, dz = 0.0; 
01913 
01914   Hep3Vector a(0.0, newZ.y(), newZ.z());
01915   // three rotated angles, rotate dx angles(rad) by X axis,
01916   // then rotate dy by new Y axis, then dz by new Z axis;
01917   // all rotate in left hand system;
01918 
01919   // a, project of final newZ on original YZ plane, 0 <= theta < pi;
01920   if(a.mag() != 0.0) a.setMag(1.0);
01921   else cout << "newZ on X axis, a=(0,0,0)" << endl;
01922   dx = acos(a.dot(z0));
01923   if(a.dot(z0) == -1.0) dx = 0.0;
01924 
01925   // b, rotate dx angle of z0 around original X axis(x0), 
01926   Hep3Vector b(0, sin(dx), cos(dx));
01927   dy = acos(b.dot(newZ));
01928   if(newZ.x() > 0.0) dy = -dy;
01929 
01930   Hep3Vector v(dx, dy, dz);
01931   return v;
01932 }

Generated on Tue Nov 29 23:14:25 2016 for BOSS_7.0.2 by  doxygen 1.4.7