00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include "G4ThreeVector.hh"
00011
00012 #include "TrkExtAlg/ExtBesEmcGeometry.h"
00013 #include "TrkExtAlg/ExtBesEmcParameter.h"
00014 #include "G4ios.hh"
00015
00016 ExtBesEmcGeometry::ExtBesEmcGeometry()
00017 {
00018 verboseLevel = 0;
00019 }
00020
00021 ExtBesEmcGeometry::~ExtBesEmcGeometry()
00022 {;}
00023
00024 void ExtBesEmcGeometry::ReadEMCParameters()
00025 {
00026
00027 ExtBesEmcParameter& emcPara=ExtBesEmcParameter::GetInstance();
00028
00029
00030
00031 TaperRingDz = emcPara.GetTaperRingDz();
00032 TaperRingThickness1 = emcPara.GetTaperRingThickness1();
00033 TaperRingThickness2 = emcPara.GetTaperRingThickness2();
00034 TaperRingThickness3 = emcPara.GetTaperRingThickness3();
00035 TaperRingTheta = emcPara.GetTaperRingTheta()*deg;
00036 TaperRingInnerLength = emcPara.GetTaperRingInnerLength();
00037 TaperRingOuterLength = emcPara.GetTaperRingOuterLength();
00038
00039 BSCRmin = emcPara.GetBSCRmin();
00040 BSCDz = emcPara.GetBSCDz()-TaperRingThickness3;
00041
00042 BSCRmin1 = emcPara.GetBSCRmin1();
00043 BSCRmax1 = emcPara.GetBSCRmax1();
00044 BSCRmin2 = emcPara.GetBSCRmin2();
00045
00046 BSCDz1 = emcPara.GetBSCDz1()-TaperRingThickness3;
00047
00048 BSCAngleRotat = emcPara.GetBSCAngleRotat()*deg;
00049 BSCNbPhi = emcPara.GetBSCNbPhi();
00050 BSCNbTheta = emcPara.GetBSCNbTheta();
00051
00052 BSCCryLength = emcPara.GetCrystalLength();
00053 BSCYFront0 = emcPara.GetBSCYFront0();
00054 BSCYFront = emcPara.GetBSCYFront();
00055 BSCYFront1 = emcPara.GetBSCYFront1();
00056 BSCPosition0 = emcPara.GetBSCPosition0();
00057 BSCPosition1 = emcPara.GetBSCPosition1();
00058
00059 fTyvekThickness = emcPara.GetTyvekThickness();
00060 fAlThickness = emcPara.GetAlThickness();
00061 fMylarThickness = emcPara.GetMylarThickness();
00062
00063 rearBoxLength = emcPara.GetRearBoxLength();
00064 rearBoxDz = emcPara.GetRearBoxDz();
00065
00066 HangingPlateDz = emcPara.GetHangingPlateDz();
00067 OCGirderAngle = emcPara.GetOCGirderAngle()*deg;
00068
00069 rearCasingThickness = emcPara.GetRearCasingThickness();
00070
00071 orgGlassLengthX = emcPara.GetOrgGlassLengthX();
00072 orgGlassLengthY = emcPara.GetOrgGlassLengthY();
00073 orgGlassLengthZ = emcPara.GetOrgGlassLengthZ();
00074
00075 PDLengthX = emcPara.GetPDLengthX();
00076 PDLengthY = emcPara.GetPDLengthY();
00077 PDLengthZ = emcPara.GetPDLengthZ();
00078
00079 AlPlateDz = emcPara.GetAlPlateDz();
00080 PABoxDz = emcPara.GetPABoxDz();
00081 PABoxThickness = emcPara.GetPABoxThickness();
00082
00083 cableDr = emcPara.GetCableDr();
00084 waterPipeDr = emcPara.GetWaterPipeDr();
00085 waterPipeThickness= emcPara.GetWaterPipeThickness();
00086
00087 SPBarThickness = emcPara.GetSPBarThickness();
00088 SPBarThickness1 = emcPara.GetSPBarThickness1();
00089 SPBarwidth = emcPara.GetSPBarwidth();
00090
00091 EndRingDz = emcPara.GetEndRingDz();
00092 EndRingDr = emcPara.GetEndRingDr();
00093 EndRingRmin = emcPara.GetEndRingRmin();
00094
00095 for(G4int i=0;i<BSCNbTheta*6;i++)
00096 {
00097 zHalfLength[i] =0;
00098 thetaAxis[i] =0;
00099 phiAxis [i] =0;
00100 yHalfLength1[i] =0;
00101 xHalfLength1[i] =0;
00102 xHalfLength2[i] =0;
00103 tanAlpha1[i] =0;
00104 yHalfLength2[i] =0;
00105 xHalfLength3[i] =0;
00106 xHalfLength4[i] =0;
00107 tanAlpha2[i] =0;
00108 thetaPosition[i]=0;
00109 xPosition[i] =0;
00110 yPosition[i] =0;
00111 zPosition[i] =0;
00112 }
00113 }
00114
00115 void ExtBesEmcGeometry::PrintEMCParameters()
00116 {
00117
00118 ;
00119 }
00120
00121 void ExtBesEmcGeometry::ComputeEMCParameters()
00122 {
00123 ReadEMCParameters();
00124
00125
00126 G4double theta0;
00127 G4int i;
00128
00129
00130
00131
00132
00133 const G4double delta=1.*mm;
00134 TaperRingRmin1 = BSCRmax1-TaperRingThickness1;
00135 TaperRingRmin2 = TaperRingRmin1+TaperRingDz*tan(TaperRingTheta);
00136 TaperRingDr = TaperRingThickness2/cos(TaperRingTheta);
00137 TaperRingOuterLength1 = TaperRingOuterLength+TaperRingThickness3*tan(TaperRingTheta);
00138
00139 BSCRmax2 = TaperRingRmin2+TaperRingDr-TaperRingThickness3*tan(TaperRingTheta);
00140 BSCRmax = BSCRmin+33.8*cm;
00141 BSCPhiDphi= 360.*deg/BSCNbPhi;
00142 BSCPhiSphi= -BSCPhiDphi/2.;
00143 BSCPhiDz = BSCDz;
00144 BSCPhiRmax= BSCRmax-0.5*cm;
00145 BSCPhiRmin= BSCRmin/sin(90.*deg+BSCPhiDphi/2.)
00146 *sin(90.*deg+BSCPhiDphi/2.-BSCAngleRotat);
00147 SPBarDphi=2*atan(SPBarwidth/(2*BSCRmax));
00148
00149
00150 zHalfLength[0] = BSCCryLength/2.;
00151 yHalfLength1[0]= BSCYFront0/2.;
00152 xHalfLength1[0]= BSCPhiRmin*tan(BSCPhiDphi)/2.;
00153 xHalfLength2[0]= xHalfLength1[0];
00154 G4double rminprojected=BSCRmin*cos(BSCAngleRotat);
00155
00156 yHalfLength2[0]= yHalfLength1[0]
00157 +(BSCYFront0-BSCPosition0)/rminprojected*BSCCryLength/2.;
00158 xHalfLength3[0]= xHalfLength1[0]*(BSCPhiRmin+BSCCryLength)/BSCPhiRmin;
00159 xHalfLength4[0]= xHalfLength2[0]*(BSCPhiRmin+BSCCryLength)/BSCPhiRmin;
00160
00161 thetaPosition[0]= 90.*deg;
00162 xPosition[0] = BSCPhiRmin+BSCCryLength/2.;
00163 yPosition[0] =
00164 -(xHalfLength1[0]+xHalfLength3[0]+xHalfLength2[0]+xHalfLength4[0])/4.;
00165 zPosition[0] = (yHalfLength1[0]+yHalfLength2[0])/2.;
00166
00167 theta0= 90.*deg-atan((BSCYFront0-BSCPosition0)/rminprojected);
00168 if(verboseLevel>1)
00169 {
00170 G4cout << "------------------------------------>1"<< G4endl
00171 << "The direction of crystal is:"
00172 << theta0/deg <<"(deg)."<< G4endl;
00173 }
00174
00175 rearBoxPosX[0] = xPosition[0]+BSCCryLength/2.+rearBoxDz/2.;
00176 rearBoxPosY[0] = -(xHalfLength3[0]+xHalfLength4[0])/2.;
00177 rearBoxPosZ[0] = yHalfLength2[0];
00178
00179 OCGirderRmin1[0] = rearBoxPosX[0]+rearBoxDz/2.+delta;
00180 OCGirderRmin2[0] = rearBoxPosX[0]+rearBoxDz/2.+delta;
00181 OCGirderDz[0] = rearBoxLength;
00182 OCGirderPosZ[0] = rearBoxLength/2;
00183
00184 cableLength[0] = BSCDz-rearBoxPosZ[0];
00185 cablePosX[0] = BSCPhiRmax-cableDr-2.*mm;
00186 cablePosY[0] = rearBoxPosY[0]-3*cableDr;
00187 cablePosZ[0] = BSCDz-cableLength[0]/2;
00188
00189
00190 zHalfLength[1]= BSCCryLength/2.;
00191 yHalfLength1[1]= BSCYFront0/2.;
00192 G4double deltaR= BSCYFront0*cos(theta0);
00193 xHalfLength1[1]= xHalfLength1[0];
00194 xHalfLength2[1]= xHalfLength1[1]*(BSCPhiRmin+deltaR)/BSCPhiRmin;
00195 yHalfLength2[1]= yHalfLength1[1]
00196 +tan(theta0-atan(rminprojected/(BSCYFront0*(1.+1./sin(theta0))-
00197 BSCPosition1)))*BSCCryLength/2.;
00198 deltaR=deltaR+BSCCryLength*sin(theta0);
00199 xHalfLength4[1]= xHalfLength1[1]*(BSCPhiRmin+deltaR)/BSCPhiRmin;
00200 deltaR=deltaR-yHalfLength2[1]*2.*cos(theta0);
00201 xHalfLength3[1]= xHalfLength1[1]*(BSCPhiRmin+deltaR)/BSCPhiRmin;
00202
00203 thetaPosition[1]= theta0;
00204 xPosition[1] = BSCPhiRmin+BSCCryLength/2.*sin(theta0)
00205 + (3.*yHalfLength1[1]-yHalfLength2[1])/2.*cos(theta0);
00206 yPosition[1] =
00207 -(xHalfLength1[1]+xHalfLength3[1]+xHalfLength2[1]+xHalfLength4[1])/4.;
00208 zPosition[1] = yHalfLength1[0]*2.
00209 + (yHalfLength1[1]*2./tan(theta0)+BSCCryLength/2.)*cos(theta0)
00210 + (yHalfLength1[1]+yHalfLength2[1])/2.*sin(theta0);
00211
00212 rearBoxPosX[1] = xPosition[1]+(BSCCryLength/2.+rearBoxDz/2.)*sin(theta0);
00213 rearBoxPosY[1] = -(xHalfLength3[1]+xHalfLength4[1])/2.;
00214 rearBoxPosZ[1] = zPosition[1]+(BSCCryLength/2.+rearBoxDz/2.)*cos(theta0);
00215
00216 OCGirderRmin1[1] = rearBoxPosX[1]+rearBoxDz*sin(theta0)/2.+rearBoxLength*cos(theta0)/2+delta;
00217 OCGirderRmin2[1] = rearBoxPosX[1]+rearBoxDz*sin(theta0)/2.-rearBoxLength*cos(theta0)/2+delta;
00218 OCGirderDz[1] = rearBoxLength*sin(theta0);
00219 OCGirderPosZ[1] = rearBoxPosZ[1]+rearBoxDz*cos(theta0)/2.;
00220
00221 cableLength[1] = BSCDz-rearBoxPosZ[1];
00222 cablePosX[1] = cablePosX[0];
00223 cablePosY[1] = cablePosY[0]+2*cableDr;
00224 cablePosZ[1] = BSCDz-cableLength[1]/2;
00225
00226 theta0= theta0-atan((yHalfLength2[1]-yHalfLength1[1])*2./BSCCryLength);
00227
00228 for(i=2;i<BSCNbTheta;i++)
00229 {
00230 if(verboseLevel>1)
00231 {
00232 G4cout << "------------------------------------>"<<i<< G4endl
00233 << "The direction of crystal is:"
00234 << theta0/deg <<"(deg)." << G4endl;
00235 }
00236
00237 zHalfLength[i]=zHalfLength[0];
00238 yHalfLength1[i]=BSCYFront/2.;
00239 if(i==BSCNbTheta-1)
00240 {yHalfLength1[i]=BSCYFront1/2.;}
00241 deltaR=yHalfLength1[i]*2.*cos(theta0);
00242 xHalfLength1[i]=xHalfLength1[0];
00243 xHalfLength2[i]=xHalfLength1[i]/BSCPhiRmin*(BSCPhiRmin+deltaR);
00244 yHalfLength2[i]=yHalfLength1[i]
00245 *(1.+BSCCryLength/(rminprojected/sin(theta0)
00246 +yHalfLength1[i]*2./tan(theta0)));
00247 deltaR=deltaR+BSCCryLength*sin(theta0);
00248 xHalfLength4[i]=xHalfLength1[i]/BSCPhiRmin*(BSCPhiRmin+deltaR);
00249 deltaR=deltaR-yHalfLength2[i]*2.*cos(theta0);
00250 xHalfLength3[i]=xHalfLength1[i]/BSCPhiRmin*(BSCPhiRmin+deltaR);
00251
00252 thetaPosition[i]=theta0;
00253 xPosition[i]=BSCPhiRmin+zHalfLength[i]*sin(theta0)
00254 + (3.*yHalfLength1[i]-yHalfLength2[i])/2.*cos(theta0);
00255 yPosition[i]=
00256 -(xHalfLength1[i]+xHalfLength3[i]+xHalfLength2[i]+xHalfLength4[i])/4.;
00257 zPosition[i]=BSCPosition1+rminprojected/tan(theta0)
00258 +(2.*yHalfLength1[i]/tan(theta0)+zHalfLength[i])*cos(theta0)
00259 +(yHalfLength1[i]+yHalfLength2[i])/2.*sin(theta0);
00260
00261 rearBoxPosX[i] = xPosition[i]+(BSCCryLength/2.+rearBoxDz/2.)*sin(theta0);
00262 rearBoxPosY[i] = -(xHalfLength3[i]+xHalfLength4[i])/2.;
00263 rearBoxPosZ[i] = zPosition[i]+(BSCCryLength/2.+rearBoxDz/2.)*cos(theta0);
00264
00265 OCGirderRmin1[i] = rearBoxPosX[i]+rearBoxDz*sin(theta0)/2.+rearBoxLength*cos(theta0)/2+delta;
00266 OCGirderRmin2[i] = rearBoxPosX[i]+rearBoxDz*sin(theta0)/2.-rearBoxLength*cos(theta0)/2+delta;
00267 OCGirderDz[i] = rearBoxLength*sin(theta0);
00268 OCGirderPosZ[i] = rearBoxPosZ[i]+rearBoxDz*cos(theta0)/2.;
00269
00270 G4int xCable,yCable;
00271 xCable = i/4;
00272 yCable = i-4*(G4int)(i/4);
00273 cableLength[i] = BSCDz-(rearBoxPosZ[i]+rearBoxDz/2.*cos(theta0));
00274 cablePosX[i] = cablePosX[0]-xCable*2*cableDr;
00275 cablePosY[i] = cablePosY[0]+yCable*2*cableDr;
00276 cablePosZ[i] = BSCDz-cableLength[i]/2;
00277
00278 theta0=theta0-atan(2.*yHalfLength1[i]/(rminprojected/sin(theta0)
00279 +2.*yHalfLength1[i]/tan(theta0)));
00280
00281 }
00282 thetaPosition[BSCNbTheta]=theta0;
00283 if(verboseLevel>1)
00284 {
00285 G4cout << "------------------------------------>"<<i<< G4endl
00286 << "The direction of crystal is:"
00287 << theta0/deg <<"(deg)." << G4endl;
00288 }
00289
00290 G4double oop;
00291 for(i=0;i<BSCNbTheta;i++)
00292 {
00293 oop=sqrt((yHalfLength2[i]-yHalfLength1[i])*(yHalfLength2[i]
00294 -yHalfLength1[i])
00295 +(xHalfLength3[i]+xHalfLength4[i]-xHalfLength1[i]
00296 -xHalfLength2[i])*(xHalfLength3[i]+xHalfLength4[i]
00297 -xHalfLength1[i]-xHalfLength2[i])/4);
00298 thetaAxis[i]=atan(oop/BSCCryLength);
00299 phiAxis[i] =180.*deg+atan((yHalfLength2[i]-yHalfLength1[i])
00300 /(xHalfLength3[i]+xHalfLength4[i]
00301 -xHalfLength1[i]-xHalfLength2[i])*2.);
00302 tanAlpha2[i]=-(xHalfLength4[i]-xHalfLength3[i])/yHalfLength2[i]/2.;
00303 tanAlpha1[i]=-(xHalfLength2[i]-xHalfLength1[i])/yHalfLength1[i]/2.;
00304
00305
00306
00307
00308 }
00309
00310 }
00311
00312 void ExtBesEmcGeometry::ModifyForCasing()
00313 {
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330
00331
00332 G4double totalThickness=fTyvekThickness+fAlThickness+fMylarThickness;
00333 G4double delta=0.,angle1=0.*deg,angle2=0.*deg;
00334 G4double oop;
00335
00336 G4double rminprojected=BSCRmin*cos(BSCAngleRotat);
00337 rminprojected=BSCRmin;
00338
00339 G4int i;
00340 for(i=0;i<BSCNbTheta;i++)
00341 {
00342 zHalfLength[BSCNbTheta+i]=totalThickness/2.;
00343 yHalfLength1[BSCNbTheta+i]=yHalfLength1[i];
00344 yHalfLength2[BSCNbTheta+i]=yHalfLength1[i]
00345 +(yHalfLength2[i]-yHalfLength1[i])*totalThickness/BSCCryLength;
00346 xHalfLength1[BSCNbTheta+i]=xHalfLength1[i];
00347 xHalfLength2[BSCNbTheta+i]=xHalfLength2[i];
00348 xHalfLength1[BSCNbTheta*2+i]=xHalfLength3[BSCNbTheta+i]=
00349 xHalfLength1[i]*(BSCCryLength-totalThickness)/BSCCryLength
00350 +xHalfLength3[i]*totalThickness/BSCCryLength;
00351 xHalfLength2[BSCNbTheta*4+i]=xHalfLength4[BSCNbTheta+i]=
00352 xHalfLength2[i]*(BSCCryLength-totalThickness)/BSCCryLength
00353 +xHalfLength4[i]*totalThickness/BSCCryLength;
00354
00355 zHalfLength[BSCNbTheta*5+i]=zHalfLength[BSCNbTheta*4+i]=
00356 zHalfLength[BSCNbTheta*3+i]=zHalfLength[BSCNbTheta*2+i]=
00357 zHalfLength[i]-totalThickness/2.;
00358
00359 yHalfLength2[BSCNbTheta*2+i]=yHalfLength1[BSCNbTheta*2+i]=
00360 totalThickness/cos(thetaPosition[i+1]-thetaPosition[i])/2.;
00361 xHalfLength3[BSCNbTheta*2+i]=xHalfLength3[i];
00362 xHalfLength4[BSCNbTheta*2+i]=xHalfLength3[i]
00363 +(xHalfLength4[i]-xHalfLength3[i])*yHalfLength2[BSCNbTheta*2+i]
00364 /yHalfLength2[i];
00365 xHalfLength2[BSCNbTheta*2+i]=xHalfLength3[BSCNbTheta+i]
00366 +(xHalfLength4[BSCNbTheta+i]-xHalfLength3[BSCNbTheta+i])
00367 *yHalfLength1[BSCNbTheta*2+i]/yHalfLength2[BSCNbTheta*1+i];
00368
00369 yHalfLength2[BSCNbTheta*4+i]=yHalfLength1[BSCNbTheta*4+i]=
00370 totalThickness/2.;
00371 xHalfLength4[BSCNbTheta*4+i]=xHalfLength4[i];
00372 xHalfLength3[BSCNbTheta*4+i]=xHalfLength4[i]
00373 -(xHalfLength4[i]-xHalfLength3[i])*yHalfLength2[BSCNbTheta*4+i]
00374 /yHalfLength2[i];
00375 xHalfLength1[BSCNbTheta*4+i]=xHalfLength4[BSCNbTheta+i]
00376 -(xHalfLength4[BSCNbTheta+i]-xHalfLength3[BSCNbTheta+i])
00377 *yHalfLength1[BSCNbTheta*4+i]/yHalfLength2[BSCNbTheta*1+i];
00378
00379 delta=totalThickness/2.+yHalfLength1[BSCNbTheta*2+i];
00380 angle1=atan(yHalfLength2[i]/(xHalfLength4[i]-xHalfLength3[i]));
00381 angle2=atan(2.*(xHalfLength4[i]-xHalfLength2[i])*sin(angle1)
00382 /BSCCryLength);
00383
00384 yHalfLength1[BSCNbTheta*5+i]=yHalfLength1[BSCNbTheta*3+i]=
00385 yHalfLength1[i]-delta;
00386 yHalfLength2[BSCNbTheta*5+i]=yHalfLength2[BSCNbTheta*3+i]=
00387 yHalfLength2[i]-delta;
00388 xHalfLength4[BSCNbTheta*3+i]=xHalfLength3[BSCNbTheta*3+i]=
00389 xHalfLength2[BSCNbTheta*3+i]=xHalfLength1[BSCNbTheta*3+i]=
00390 totalThickness/cos(angle2)/sin(angle1)/2.;
00391 xHalfLength4[BSCNbTheta*5+i]=xHalfLength3[BSCNbTheta*5+i]=
00392 xHalfLength2[BSCNbTheta*5+i]=xHalfLength1[BSCNbTheta*5+i]=
00393 totalThickness/2.;
00394
00395 zHalfLength[i]=zHalfLength[i]-totalThickness/2.;
00396 yHalfLength1[i]=yHalfLength1[i]-delta;
00397 yHalfLength2[i]=yHalfLength2[i]-delta;
00398 delta=totalThickness*(1.+1./cos(angle2)/sin(angle1))/2.;
00399 xHalfLength1[i]=xHalfLength1[i]-delta;
00400 xHalfLength2[i]=xHalfLength2[i]-delta;
00401 xHalfLength3[i]=xHalfLength3[i]-delta;
00402 xHalfLength4[i]=xHalfLength4[i]-delta;
00403
00404 oop=sqrt((yHalfLength2[i]-yHalfLength1[i])*(yHalfLength2[i]
00405 -yHalfLength1[i])
00406 +(xHalfLength3[i]+xHalfLength4[i]-xHalfLength1[i]
00407 -xHalfLength2[i])*(xHalfLength3[i]+xHalfLength4[i]
00408 -xHalfLength1[i]-xHalfLength2[i])/4);
00409 thetaAxis[i]=atan(oop/BSCCryLength);
00410
00411 phiAxis[i] =180.*deg+atan((yHalfLength2[i]-yHalfLength1[i])
00412 /(xHalfLength3[i]+xHalfLength4[i]
00413 -xHalfLength1[i]-xHalfLength2[i])*2.);
00414
00415 oop=sqrt((yHalfLength2[BSCNbTheta+i]-yHalfLength1[BSCNbTheta+i])
00416 *(yHalfLength2[BSCNbTheta+i]-yHalfLength1[BSCNbTheta+i])
00417 +(xHalfLength3[BSCNbTheta+i]+xHalfLength4[BSCNbTheta+i]
00418 -xHalfLength1[BSCNbTheta+i]-xHalfLength2[BSCNbTheta+i])
00419 *(xHalfLength3[BSCNbTheta+i]+xHalfLength4[BSCNbTheta+i]
00420 -xHalfLength1[BSCNbTheta+i]-xHalfLength2[BSCNbTheta+i])/4);
00421 thetaAxis[BSCNbTheta+i]=atan(oop/totalThickness);
00422 phiAxis [BSCNbTheta+i]=
00423 -atan((yHalfLength2[BSCNbTheta+i]-yHalfLength1[BSCNbTheta+i])
00424 /(xHalfLength3[BSCNbTheta+i]+xHalfLength4[BSCNbTheta+i]
00425 -xHalfLength1[BSCNbTheta+i]-xHalfLength2[BSCNbTheta+i])*2.);
00426
00427 oop=sqrt((yHalfLength2[i]-yHalfLength1[i])*(yHalfLength2[i]
00428 -yHalfLength1[i])*4
00429 +(xHalfLength3[BSCNbTheta*2+i]+xHalfLength4[BSCNbTheta*2+i]
00430 -xHalfLength1[BSCNbTheta*2+i]-xHalfLength2[BSCNbTheta*2+i])
00431 *(xHalfLength3[BSCNbTheta*2+i]+xHalfLength4[BSCNbTheta*2+i]
00432 -xHalfLength1[BSCNbTheta*2+i]-xHalfLength2[BSCNbTheta*2+i])
00433 /4);
00434 thetaAxis[BSCNbTheta*2+i]=atan(oop/(zHalfLength[BSCNbTheta*2+i]*2));
00435 phiAxis [BSCNbTheta*2+i]=
00436 -atan((yHalfLength2[i]-yHalfLength1[i])
00437 /(xHalfLength3[BSCNbTheta*2+i]+xHalfLength4[BSCNbTheta*2+i]
00438 -xHalfLength1[BSCNbTheta*2+i]-xHalfLength2[BSCNbTheta*2+i])*4);
00439
00440 oop=sqrt((yHalfLength2[i]-yHalfLength1[i])*(yHalfLength2[i]
00441 -yHalfLength1[i])*4
00442 +(xHalfLength4[i]-xHalfLength2[i])
00443 *(xHalfLength4[i]-xHalfLength2[i])*4);
00444 thetaAxis[BSCNbTheta*3+i]=atan(oop/(zHalfLength[BSCNbTheta*3+i]*2));
00445 phiAxis [BSCNbTheta*3+i]=-atan((yHalfLength2[i]-yHalfLength1[i])
00446 /(xHalfLength4[i]-xHalfLength2[i]));
00447
00448 thetaAxis[BSCNbTheta*4+i]=
00449 atan((xHalfLength4[BSCNbTheta*4+i]+xHalfLength3[BSCNbTheta*4+i]
00450 -xHalfLength2[BSCNbTheta*4+i]-xHalfLength1[BSCNbTheta*4+i])/2.
00451 /(zHalfLength[BSCNbTheta*4+i]*2));
00452 phiAxis [BSCNbTheta*4+i]=0;
00453
00454 thetaAxis[BSCNbTheta*5+i]=atan((xHalfLength3[BSCNbTheta*5+i]
00455 -xHalfLength1[BSCNbTheta*5+i])
00456 /(zHalfLength[BSCNbTheta*5+i]*2));
00457 phiAxis [BSCNbTheta*5+i]=-90.*deg;
00458
00459 tanAlpha2[BSCNbTheta+i]=tanAlpha1[BSCNbTheta+i]=tanAlpha1[i]=
00460 -(xHalfLength2[i]-xHalfLength1[i])/yHalfLength1[i]/2.;
00461 tanAlpha1[BSCNbTheta*2+i]=(xHalfLength2[BSCNbTheta*2+i]
00462 -xHalfLength1[BSCNbTheta*2+i])
00463 /yHalfLength1[BSCNbTheta*2+i]/2.;
00464 tanAlpha1[BSCNbTheta*3+i]=tanAlpha1[i]*2.;
00465 tanAlpha1[BSCNbTheta*4+i]=(xHalfLength2[BSCNbTheta*4+i]
00466 -xHalfLength1[BSCNbTheta*4+i])
00467 /yHalfLength1[BSCNbTheta*4+i]/2.;
00468 tanAlpha1[BSCNbTheta*5+i]=(xHalfLength2[BSCNbTheta*5+i]
00469 -xHalfLength1[BSCNbTheta*5+i])
00470 /yHalfLength1[BSCNbTheta*5+i]/2.;
00471
00472 tanAlpha2[i]=-(xHalfLength4[i]-xHalfLength3[i])/yHalfLength2[i]/2.;
00473
00474 tanAlpha2[BSCNbTheta*2+i]=(xHalfLength4[BSCNbTheta*2+i]
00475 -xHalfLength3[BSCNbTheta*2+i])
00476 /yHalfLength2[BSCNbTheta*2+i]/2.;
00477 tanAlpha2[BSCNbTheta*3+i]=tanAlpha2[i]*2.;
00478 tanAlpha2[BSCNbTheta*4+i]=(xHalfLength4[BSCNbTheta*4+i]
00479 -xHalfLength3[BSCNbTheta*4+i])
00480 /yHalfLength2[BSCNbTheta*4+i]/2.;
00481 tanAlpha2[BSCNbTheta*5+i]=(xHalfLength4[BSCNbTheta*5+i]
00482 -xHalfLength3[BSCNbTheta*5+i])
00483 /yHalfLength2[BSCNbTheta*5+i]/2.;
00484
00485 zPosition[BSCNbTheta*5+i]=zPosition[BSCNbTheta*3+i]=zPosition[i]=
00486 zPosition[i]+totalThickness/2.*cos(thetaPosition[i])
00487 -yHalfLength1[BSCNbTheta*2+i]*sin(thetaPosition[i]);
00488 zPosition[i]=totalThickness/2.;
00489 xPosition[BSCNbTheta*5+i]=xPosition[BSCNbTheta*3+i]=xPosition[i]=
00490 xPosition[i]+totalThickness/2.*sin(thetaPosition[i])
00491 +totalThickness*(1./cos(thetaPosition[i+1]-thetaPosition[i])-1)/2.
00492 *cos(thetaPosition[i]);
00493 xPosition[i]=totalThickness*(1.-1./cos(angle2)/sin(angle1))/2.;
00494
00495 yPosition[i]=yPosition[i]
00496 +totalThickness*(1.-1./cos(angle2)/sin(angle1))/2.;
00497 yPosition[i]=yHalfLength1[BSCNbTheta*2+i]-totalThickness/2.;
00498 yPosition[BSCNbTheta*3+i]=yPosition[i]*2.+xHalfLength1[BSCNbTheta*3+i];
00499 yPosition[BSCNbTheta*5+i]=xHalfLength1[BSCNbTheta*5+i];
00500
00501 xPosition[BSCNbTheta+i]=BSCPhiRmin
00502 +zHalfLength[BSCNbTheta+i]*sin(thetaPosition[i])
00503 +(3.*yHalfLength1[BSCNbTheta+i]-yHalfLength2[BSCNbTheta+i])/2.
00504 *cos(thetaPosition[i]);
00505 yPosition[BSCNbTheta+i]=(xHalfLength1[BSCNbTheta+i]
00506 +xHalfLength3[BSCNbTheta+i]
00507 +xHalfLength2[BSCNbTheta+i]
00508 +xHalfLength4[BSCNbTheta+i])/4.;
00509 zPosition[BSCNbTheta+i]=BSCPosition1+rminprojected/tan(thetaPosition[i])
00510 +(2.*yHalfLength1[BSCNbTheta+i]/tan(thetaPosition[i])
00511 +zHalfLength[BSCNbTheta+i])*cos(thetaPosition[i])
00512 +(yHalfLength1[BSCNbTheta+i]+yHalfLength2[BSCNbTheta+i])/2.
00513 *sin(thetaPosition[i]);
00514
00515 xPosition[BSCNbTheta*2+i]=xPosition[i]
00516 +((yHalfLength1[i]+yHalfLength2[i])/2.+yHalfLength1[BSCNbTheta*2+i])
00517 *cos(thetaPosition[i]);
00518 zPosition[BSCNbTheta*2+i]=zPosition[i]
00519 -((yHalfLength1[i]+yHalfLength2[i])/2.+yHalfLength1[BSCNbTheta*2+i])
00520 *sin(thetaPosition[i]);
00521 yPosition[BSCNbTheta*2+i]=(xHalfLength1[BSCNbTheta*2+i]
00522 +xHalfLength3[BSCNbTheta*2+i]
00523 +xHalfLength2[BSCNbTheta*2+i]
00524 +xHalfLength4[BSCNbTheta*2+i])/4.;
00525
00526 xPosition[BSCNbTheta*4+i]=xPosition[i]
00527 -((yHalfLength1[i]+yHalfLength2[i])/2.+yHalfLength1[BSCNbTheta*4+i])
00528 *cos(thetaPosition[i]);
00529 zPosition[BSCNbTheta*4+i]=zPosition[i]
00530 -((yHalfLength1[i]+yHalfLength2[i])/2.+yHalfLength1[BSCNbTheta*4+i])
00531 *sin(thetaPosition[i]);
00532 yPosition[BSCNbTheta*4+i]=(xHalfLength1[BSCNbTheta*4+i]
00533 +xHalfLength3[BSCNbTheta*4+i]
00534 +xHalfLength2[BSCNbTheta*4+i]
00535 +xHalfLength4[BSCNbTheta*4+i])/4.;
00536
00537 }
00538
00539 if(verboseLevel>1)
00540 for(i=0;i<BSCNbTheta*6;i++)
00541 {
00542 G4cout << "The sizes of the " << i+1 << " crystal are:" << G4endl
00543 << "zHalfLength =" << zHalfLength[i]/cm << "(cm)," << G4endl
00544 << "thetaAxis =" << thetaAxis[i]/deg << "(deg),"<< G4endl
00545 << "phiAxis =" << phiAxis[i]/deg << "(deg),"<< G4endl
00546 << "yHalfLength1=" << yHalfLength1[i]/cm << "(cm)," << G4endl
00547 << "xHalfLength1=" << xHalfLength1[i]/cm << "(cm)," << G4endl
00548 << "xHalfLength2=" << xHalfLength2[i]/cm << "(cm)," << G4endl
00549 << "tanAlpha1 =" << tanAlpha1[i] << G4endl
00550 << "yHalfLength2=" << yHalfLength2[i]/cm << "(cm)," << G4endl
00551 << "xHalfLength3=" << xHalfLength3[i]/cm << "(cm)," << G4endl
00552 << "xHalfLength4=" << xHalfLength4[i]/cm << "(cm)," << G4endl
00553 << "tanAlpha2 =" << tanAlpha2[i] << "." << G4endl
00554 << "The position of the " << i+1 << " crystal is:" << G4endl
00555 << "(" << xPosition[i]/cm << ","
00556 << yPosition[i]/cm << ","
00557 << zPosition[i]/cm << ")cm" << G4endl;
00558 }
00559
00560 }
00561
00562 void ExtBesEmcGeometry::SetCasingThickness(G4ThreeVector val)
00563 {
00564
00565 fTyvekThickness = val('X');
00566 fAlThickness = val('Y');
00567 fMylarThickness = val('Z');
00568 }
00569
00570 G4double ExtBesEmcGeometry::GetXPosition(G4int NbCrystal)
00571 {
00572 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00573 {
00574 return xPosition[NbCrystal];
00575 }
00576 else
00577 {
00578 return 0.;
00579 }
00580 }
00581
00582 G4double ExtBesEmcGeometry::GetYPosition(G4int NbCrystal)
00583 {
00584 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00585 {
00586 return yPosition[NbCrystal];
00587 }
00588 else
00589 {
00590 return 0.;
00591 }
00592 }
00593
00594 G4double ExtBesEmcGeometry::GetZPosition(G4int NbCrystal)
00595 {
00596 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00597 {
00598 return zPosition[NbCrystal];
00599 }
00600 else
00601 {
00602 return 0.;
00603 }
00604 }
00605
00606 G4double ExtBesEmcGeometry::GetThetaPosition(G4int NbCrystal)
00607 {
00608 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00609 {
00610 return thetaPosition[NbCrystal];
00611 }
00612 else
00613 {
00614 return 0.;
00615 }
00616 }
00617
00618 G4double ExtBesEmcGeometry::GetZHalfLength(G4int NbCrystal)
00619 {
00620 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00621 {
00622 return zHalfLength[NbCrystal];
00623 }
00624 else
00625 {
00626 return 0.;
00627 }
00628 }
00629
00630 G4double ExtBesEmcGeometry::GetThetaAxis(G4int NbCrystal)
00631 {
00632 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00633 {
00634 return thetaAxis[NbCrystal];
00635 }
00636 else
00637 {
00638 return 0.;
00639 }
00640 }
00641
00642 G4double ExtBesEmcGeometry::GetPhiAxis(G4int NbCrystal)
00643 {
00644 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00645 {
00646 return phiAxis[NbCrystal];
00647 }
00648 else
00649 {
00650 return 0.;
00651 }
00652 }
00653
00654 G4double ExtBesEmcGeometry::GetYHalfLength1(G4int NbCrystal)
00655 {
00656 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00657 {
00658 return yHalfLength1[NbCrystal];
00659 }
00660 else
00661 {
00662 return 0.;
00663 }
00664 }
00665
00666 G4double ExtBesEmcGeometry::GetXHalfLength1(G4int NbCrystal)
00667 {
00668 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00669 {
00670 return xHalfLength1[NbCrystal];
00671 }
00672 else
00673 {
00674 return 0.;
00675 }
00676 }
00677
00678 G4double ExtBesEmcGeometry::GetXHalfLength2(G4int NbCrystal)
00679 {
00680 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00681 {
00682 return xHalfLength2[NbCrystal];
00683 }
00684 else
00685 {
00686 return 0.;
00687 }
00688 }
00689
00690 G4double ExtBesEmcGeometry::GetTanAlpha1(G4int NbCrystal)
00691 {
00692 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00693 {
00694 return tanAlpha1[NbCrystal];
00695 }
00696 else
00697 {
00698 return 0.;
00699 }
00700 }
00701
00702 G4double ExtBesEmcGeometry::GetYHalfLength2(G4int NbCrystal)
00703 {
00704 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00705 {
00706 return yHalfLength2[NbCrystal];
00707 }
00708 else
00709 {
00710 return 0.;
00711 }
00712 }
00713
00714 G4double ExtBesEmcGeometry::GetXHalfLength3(G4int NbCrystal)
00715 {
00716 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00717 {
00718 return xHalfLength3[NbCrystal];
00719 }
00720 else
00721 {
00722 return 0.;
00723 }
00724 }
00725
00726 G4double ExtBesEmcGeometry::GetXHalfLength4(G4int NbCrystal)
00727 {
00728 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00729 {
00730 return xHalfLength4[NbCrystal];
00731 }
00732 else
00733 {
00734 return 0.;
00735 }
00736 }
00737
00738 G4double ExtBesEmcGeometry::GetTanAlpha2(G4int NbCrystal)
00739 {
00740 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00741 {
00742 return tanAlpha2[NbCrystal];
00743 }
00744 else
00745 {
00746 return 0.;
00747 }
00748 }
00749
00750 G4VPhysicalVolume* ExtBesEmcGeometry::GetPhysiBSCCrystal(G4int NbCrystal)
00751 {
00752 if(NbCrystal>=0&&NbCrystal<BSCNbTheta*6)
00753 {
00754 return physiBSCCrystal[NbCrystal];
00755 }
00756 else
00757 {
00758 return NULL;
00759 }
00760 }
