CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/GeneratorInterface/CosmicMuonGenerator/src/Point5MaterialMap.cc

Go to the documentation of this file.
00001 #ifndef Point5MaterialMap_cc
00002 #define Point5MaterialMap_cc
00003 
00004 #include <iostream>
00005 #include <cmath>
00006 
00007 #include "GeneratorInterface/CosmicMuonGenerator/interface/CosmicMuonParameters.h"
00008 
00009 
00010 
00011 inline int inPlug(double vx, double vy, double vz,
00012                   double PlugVx = PlugOnShaftVx, double PlugVz = PlugOnShaftVz) {
00013   if (vy > SurfaceOfEarth && vy < SurfaceOfEarth + PlugWidth) {
00014     if (vx > PlugVx - PlugXlength/2. && vx < PlugVx + PlugXlength/2. &&
00015         vz > PlugVz - PlugZlength/2. && vz < PlugVz + PlugZlength/2.) return Plug;
00016     if (vz >= PlugVz - PlugZlength/2. - PlugNoseZlength && vz < PlugVz - PlugZlength/2. &&
00017         vx > PlugVx - PlugNoseXlength/2. && vx < PlugVx + PlugNoseXlength/2.) return Plug;
00018   }
00019   return Unknown;
00020 }
00021 
00022 
00023 inline int inAirAfterPlug(double vx, double vy, double vz) {
00024   // particles above surface of earth
00025   if (vy >= SurfaceOfEarth) return Air;
00026   
00027   // CMS cavern (UXC 55)
00028   if (std::fabs(vz) < 26548. && sqrt((vx*1.1576)*(vx*1.1576) + vy*vy) < 15460. &&
00029       vy > -8762) return Air;
00030   
00031   // access shaft (PX 56)
00032   if (vy > 0. && vy < (SurfaceOfEarth-2250.) && 
00033       sqrt(vx*vx + (vz-Z_PX56)*(vz-Z_PX56)) < 10250.) return Air;
00034   
00035   //surface hall ground floor
00036   if (vy >= SurfaceOfEarth-2250. && vy < SurfaceOfEarth) {
00037     if (sqrt(vx*vx + (vz-Z_PX56)*(vz-Z_PX56)) < 10250.
00038         && vz-Z_PX56 > -7000. && vz-Z_PX56 < 7000.) return Air;
00039     if (vx > -2400. && vx < 2400. && vz-Z_PX56 >= -9800. && vz-Z_PX56 < -7000. )
00040       return Air;
00041   }
00042   
00043     // Shaft (PM 54)
00044   if (vy > 3233. && vy < (SurfaceOfEarth) && 
00045       sqrt((vx-26600.)*(vx-26600.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56)) < 6050.)
00046     //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56)) < 10050.)//@@@@@
00047     return Air;
00048   
00049   
00050   // Shaft (PM 56)
00051   if (vy > -3700. && vy < (SurfaceOfEarth) &&
00052       sqrt((vx-18220.)*(vx-18220.) + (vz+24227.-Z_PX56)*(vz+24227.-Z_PX56)) < 3550.)
00053       return Air;
00054   
00055   //Service cavern USC 55
00056   if (vz > -22050. && vz < 62050. &&
00057       sqrt((vx-29550.)*(vx-29550.) + (vy-3233.)*(vy-3233.)) < 9050. && 
00058       vy > -3650.) return Air;
00059   
00060   //UXC55 cavern endcap beam openings
00061   if (((vz >= -29278. && vz < -26548.) || (vz >= 26548 && vz <= 29278.)) &&
00062       sqrt(vy*vy + vx*vx) < 800.) return Air; 
00063   
00064   //Pillar access between CMS collision cavern and service cavern
00065   if (vx > 10000. && vx < 20500.+3050. && //TX54 galerie
00066       vz > 14460. && vz < 16260. &&//estimated wall thickness 1m
00067       vy > -8680. && vy < -1507.) return Air;
00068   
00069   if (vx > 10000. && vx < 16500. && //TX54 galerie 30.9 degree opening to UXC55
00070       vy > -8680. && vy < -1507. && //@vx=13300 delta(vz)=1865., 2vx=1650. delta(vz)=0
00071       vz > 14460.-1865.*(16500.-vx)/(3200.) && vz <= 14460.) return Air;
00072   
00073   if (vx > 13300. && vx < 20500.+3050. && //TX54 galerie
00074       vz > 14460. && vz < 16260. && //estimated wall thickness 1m
00075       vy > -8680. && vy < -1507.) return Air;
00076   
00077   if (vx > 26600-6050. && vx < 20500.+3050. && //TX54 going up in sewrvice cavern
00078       vz > 14460. && vz < 16260. && //estimated wall thickness 1m
00079       vy >= -1507. && vy < 1000.) return Air;
00080   
00081   //R56, LHC beam tunnel East
00082   if (vz > -85000. && vz <= -29278. && //UJ57 junction to UXC55 cavern
00083       sqrt(vy*vy + (vx-350.)*(vx-350.)) < 1900. &&
00084       vy > -1000.) return Air;
00085   
00086   //R54, LHC beam tunnel West
00087   if (vz >= 29278. && vz < 63000. && //UJ57 junction to UXC55 cavern
00088       sqrt(vy*vy + (vx-350.)*(vx-350.)) < 1900. &&
00089       vy > -1000.) return Air;
00090   
00091   //UJ56 cavern
00092   if (vz > -58875. && vz < -33927. &&
00093       sqrt(vy*vy + (vx-4450.)*(vx-4450.)) < 6750. && vy > -1000. && //and beam shielding
00094       !(vx > 2250. && vx < 4250. && (vz > -(33927.+18000.) || vz < -(33927.+19500.))) && 
00095       !(vx >= 4250. && vx < 6650. && vz > -(33927.+18000.) && vz < -(33927.+16000.)))
00096     return Air;
00097   
00098   //connection between PM56 shaft and UJ56 cavern
00099   if (vx > 9000. && vx < 18220. &&
00100       sqrt((vy-50.)*(vy-50.)+(vz+24227.-Z_PX56)*(vz+24227.-Z_PX56)) < 3550.
00101       && vy > -1000.) return Air;
00102   
00103   return Unknown;
00104   
00105 }
00106 
00107 
00108 inline int inWallAfterAir(double vx, double vy, double vz) {
00109   // phase II surface building
00110   if (vy < SurfaceOfEarth && vy >= (SurfaceOfEarth-2250.)) {
00111     if (std::fabs(vz-Z_PX56) < 30000. && std::fabs(vx) < 10950) return Wall;
00112     // foundation of crane
00113     if (std::fabs(vz-Z_PX56) < 9000. && std::fabs(vx) >= 10950 && std::fabs(vx) < 16950) 
00114       return Wall;
00115   }
00116   
00117   // CMS cavern (UXC 55)
00118   if (std::fabs(vz) < 29278. && sqrt((vx*1.1576)*(vx*1.1576) + vy*vy) < 16830. &&
00119       vy > -11762.) return Wall;
00120   
00121   // access shaft (PX 56)
00122   if (vy > 0. && vy < (SurfaceOfEarth-2250.) && //t(shaft wall)=2150.
00123       sqrt(vx*vx + (vz-Z_PX56)*(vz-Z_PX56)) < 12400.) return Wall;
00124   
00125   // Shaft (PM 54)
00126   if (vy > 3233. && vy < (SurfaceOfEarth-1000.) && //t~=t(PX56)/R(PX56)*R(PM54)
00127       sqrt((vx-26600.)*(vx-26600.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56)) 
00128       < 6050.+2150./10250.*6050.) return Wall;
00129   //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56))//@@@@@ 
00130   //< 10050.+2150./10250.*6050.) return Wall;//@@@@@
00131   else if (vy >= SurfaceOfEarth-1000. && vy < SurfaceOfEarth && //t~=t(PX56)/R(PX56)*R(PM54)
00132            sqrt((vx-26600.)*(vx-26600.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56)) 
00133            < 6050.+2150./10250.*6050. +1800.) return Wall;
00134   //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56))//@@@@@ 
00135   //< 10050.+2150./10250.*10050. +1800.) return Wall; //@@@@@@
00136   
00137   // Shaft (PM 56)
00138   if (vy > -5450. && vy < (SurfaceOfEarth-1000.) && //t~=t(PX56)/R(PX56)*R(PM56)
00139       sqrt((vx-18220.)*(vx-18220.) + (vz+24227.-Z_PX56)*(vz+24227.-Z_PX56)) 
00140       < 3550.+2150./10250.*3550.) return Wall;
00141   else if (vy > SurfaceOfEarth-1000. && vy < SurfaceOfEarth && //t~=t(PX56)/R(PX56)*R(PM56)
00142            sqrt((vx-18220.)*(vx-18220.) + (vz+24227.-Z_PX56)*(vz+24227.-Z_PX56)) 
00143            < 3550.+2150./10250.*3550. +1800.) return Wall;
00144   
00145   //Service cavern USC 55
00146   if (vz > -(22050.+1150.) && vz < (62050.+1150.) &&
00147       sqrt((vx-29550.)*(vx-29550.) + (vy-3233.)*(vy-3233.)) < 9050.+950. &&
00148       vy > -3650.-2000.) return Wall; //8762=estimate, to be checked 
00149   
00150   //Pillar between CMS collision cavern and service cavern
00151   if (vz > -29278.+1000. && vz < 29278.+1000.) {
00152     if (vy > -17985. && vy < 10410. && vx > 13300. && vx < 20500.) 
00153       return Wall;
00154     
00155     if (vy > 0. && vy < 10410. && vx > 10000. && vx <= 13300.)
00156       return Wall;
00157     
00158     if (vy > -3650.-2000. && vy < -3233. && // bottom edge between pillar and service cavern
00159         vx > 20000. && vx < 24000.) return Wall;
00160     if (vy > -11762. && vy < -5000. && // bottom edge between pillar and UXC55 cavern
00161         vx > 10500. && vx < 14000.) return Wall;
00162   }
00163   
00164   if (vy > -14000. && vy < -1450. && //TX54 galerie surrounding
00165       vz > 13460. && vz < 17260. && vx >= 20500. && vx < 24550.) 
00166     return Wall;
00167   
00168   //R56, LHC beam tunnel East
00169   if (vz > -85000. && vz < -28510.) { //UJ57 junction to UXC55 cavern
00170     if (sqrt(vy*vy + (vx-350.)*(vx-350.)) < 2250.) return Wall;
00171   }
00172   
00173   //R54, LHC beam tunnel West
00174   if (vz > 26550. && vz < 63000. && //UJ57 junction to UXC55 cavern
00175       sqrt(vy*vy + (vx-350.)*(vx-350.)) < 2250.) return Wall;
00176   
00177   //UJ56 cavern
00178   if (vz > -(58875.+500.) && vz < -(33927.-500.) &&
00179       sqrt(vy*vy + (vx-4450.)*(vx-4450.)) < (6750.+500.) && vy > -3650.) 
00180     return Wall;
00181   
00182   //connection between PM56 shaft and UJ56 cavern
00183   if (vx > 9000. && vx < 18220. &&
00184       sqrt((vy-50.)*(vy-50.)+(vz+24227.-Z_PX56)*(vz+24227.-Z_PX56)) < 3550.+500. 
00185       && vy > -3650.) return Wall;
00186   
00187   return Unknown;
00188 }
00189 
00190 
00191 
00192 inline int inClayOrRockAfterWall(double vx, double vy, double vz, double ClayWidth) {
00193   
00194   //So, it is not plug, air and wall, Check for clay
00195   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
00196     return Clay;
00197   
00198   //So, it is not plug, air, wall and clay, Check for rock
00199   if (vy < SurfaceOfEarth - ClayWidth)
00200     return Rock;
00201   
00202   return Unknown;
00203   
00204 }
00205 
00206 
00207 
00208 inline int inClayAfterWall(double vx, double vy, double vz, double ClayWidth) {
00209   
00210   //So, it is not plug, air and wall, Check for clay
00211   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
00212     return Clay;
00213   
00214   return Unknown;
00215   
00216 }
00217 
00218 
00219 
00220 inline int inRockAfterClay(double vx, double vy, double vz, double ClayWidth) {
00221   
00222   //So, it is not plug, air, wall and clay, Check for rock
00223   if (vy < SurfaceOfEarth - ClayWidth)
00224     return Rock;
00225   
00226   return Unknown;
00227   
00228 }
00229 
00230 
00231 
00232 
00233 
00234 inline int inMat(double vx, double vy, double vz,
00235                  double PlugVx = PlugOnShaftVx, double PlugVz = PlugOnShaftVz,
00236                  double ClayWidth = DefaultClayWidth) {
00237   
00238   //check for Plug
00239   if (inPlug(vx, vy, vz, PlugVx, PlugVz)) return Plug;
00240   
00241   //So, it is not plug, Check for air
00242   if (inAirAfterPlug(vx, vy, vz)) return Air;
00243   
00244   //So, it is not plug and air, Check for wall
00245   if (inWallAfterAir(vx, vy, vz)) return Wall;
00246   
00247   //So, it is not plug, air and wall, Check for clay
00248   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
00249     return Clay;
00250   
00251   //So, it is not plug, air, wall and clay, Check for rock
00252   if (vy < SurfaceOfEarth - ClayWidth)
00253     return Rock;
00254   
00255   
00256   std::cout << "Point5MaterialMap.h: Warning! No Material recognised for point: vx="
00257             << vx << " vy=" << vy << " vz=" << vz << std::endl;
00258   //Something went wrong
00259   return Unknown;
00260   
00261 }
00262 
00263 
00264 
00265 
00266 #endif