Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:13:27

0001 #ifndef Point5MaterialMap_cc
0002 #define Point5MaterialMap_cc
0003 
0004 #include <iostream>
0005 #include <cmath>
0006 
0007 #include "GeneratorInterface/CosmicMuonGenerator/interface/CosmicMuonParameters.h"
0008 
0009 inline int inPlug(double vx, double vy, double vz, double PlugVx = PlugOnShaftVx, double PlugVz = PlugOnShaftVz) {
0010   if (vy > SurfaceOfEarth && vy < SurfaceOfEarth + PlugWidth) {
0011     if (vx > PlugVx - PlugXlength / 2. && vx < PlugVx + PlugXlength / 2. && vz > PlugVz - PlugZlength / 2. &&
0012         vz < PlugVz + PlugZlength / 2.)
0013       return Plug;
0014     if (vz >= PlugVz - PlugZlength / 2. - PlugNoseZlength && vz < PlugVz - PlugZlength / 2. &&
0015         vx > PlugVx - PlugNoseXlength / 2. && vx < PlugVx + PlugNoseXlength / 2.)
0016       return Plug;
0017   }
0018   return Unknown;
0019 }
0020 
0021 inline int inAirAfterPlug(double vx, double vy, double vz) {
0022   // particles above surface of earth
0023   if (vy >= SurfaceOfEarth)
0024     return Air;
0025 
0026   // CMS cavern (UXC 55)
0027   if (std::fabs(vz) < 26548. && sqrt((vx * 1.1576) * (vx * 1.1576) + vy * vy) < 15460. && vy > -8762)
0028     return Air;
0029 
0030   // access shaft (PX 56)
0031   if (vy > 0. && vy < (SurfaceOfEarth - 2250.) && sqrt(vx * vx + (vz - Z_PX56) * (vz - Z_PX56)) < 10250.)
0032     return Air;
0033 
0034   //surface hall ground floor
0035   if (vy >= SurfaceOfEarth - 2250. && vy < SurfaceOfEarth) {
0036     if (sqrt(vx * vx + (vz - Z_PX56) * (vz - Z_PX56)) < 10250. && vz - Z_PX56 > -7000. && vz - Z_PX56 < 7000.)
0037       return Air;
0038     if (vx > -2400. && vx < 2400. && vz - Z_PX56 >= -9800. && vz - Z_PX56 < -7000.)
0039       return Air;
0040   }
0041 
0042   // Shaft (PM 54)
0043   if (vy > 3233. && vy < (SurfaceOfEarth) &&
0044       sqrt((vx - 26600.) * (vx - 26600.) + (vz - 30100. - Z_PX56) * (vz - 30100. - Z_PX56)) < 6050.)
0045     //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56)) < 10050.)//@@@@@
0046     return Air;
0047 
0048   // Shaft (PM 56)
0049   if (vy > -3700. && vy < (SurfaceOfEarth) &&
0050       sqrt((vx - 18220.) * (vx - 18220.) + (vz + 24227. - Z_PX56) * (vz + 24227. - Z_PX56)) < 3550.)
0051     return Air;
0052 
0053   //Service cavern USC 55
0054   if (vz > -22050. && vz < 62050. && sqrt((vx - 29550.) * (vx - 29550.) + (vy - 3233.) * (vy - 3233.)) < 9050. &&
0055       vy > -3650.)
0056     return Air;
0057 
0058   //UXC55 cavern endcap beam openings
0059   if (((vz >= -29278. && vz < -26548.) || (vz >= 26548 && vz <= 29278.)) && sqrt(vy * vy + vx * vx) < 800.)
0060     return Air;
0061 
0062   //Pillar access between CMS collision cavern and service cavern
0063   if (vx > 10000. && vx < 20500. + 3050. &&  //TX54 galerie
0064       vz > 14460. && vz < 16260. &&          //estimated wall thickness 1m
0065       vy > -8680. && vy < -1507.)
0066     return Air;
0067 
0068   if (vx > 10000. && vx < 16500. &&  //TX54 galerie 30.9 degree opening to UXC55
0069       vy > -8680. && vy < -1507. &&  //@vx=13300 delta(vz)=1865., 2vx=1650. delta(vz)=0
0070       vz > 14460. - 1865. * (16500. - vx) / (3200.) && vz <= 14460.)
0071     return Air;
0072 
0073   if (vx > 13300. && vx < 20500. + 3050. &&  //TX54 galerie
0074       vz > 14460. && vz < 16260. &&          //estimated wall thickness 1m
0075       vy > -8680. && vy < -1507.)
0076     return Air;
0077 
0078   if (vx > 26600 - 6050. && vx < 20500. + 3050. &&  //TX54 going up in sewrvice cavern
0079       vz > 14460. && vz < 16260. &&                 //estimated wall thickness 1m
0080       vy >= -1507. && vy < 1000.)
0081     return Air;
0082 
0083   //R56, LHC beam tunnel East
0084   if (vz > -85000. && vz <= -29278. &&  //UJ57 junction to UXC55 cavern
0085       sqrt(vy * vy + (vx - 350.) * (vx - 350.)) < 1900. && vy > -1000.)
0086     return Air;
0087 
0088   //R54, LHC beam tunnel West
0089   if (vz >= 29278. && vz < 63000. &&  //UJ57 junction to UXC55 cavern
0090       sqrt(vy * vy + (vx - 350.) * (vx - 350.)) < 1900. && vy > -1000.)
0091     return Air;
0092 
0093   //UJ56 cavern
0094   if (vz > -58875. && vz < -33927. && sqrt(vy * vy + (vx - 4450.) * (vx - 4450.)) < 6750. &&
0095       vy > -1000. &&  //and beam shielding
0096       !(vx > 2250. && vx < 4250. && (vz > -(33927. + 18000.) || vz < -(33927. + 19500.))) &&
0097       !(vx >= 4250. && vx < 6650. && vz > -(33927. + 18000.) && vz < -(33927. + 16000.)))
0098     return Air;
0099 
0100   //connection between PM56 shaft and UJ56 cavern
0101   if (vx > 9000. && vx < 18220. &&
0102       sqrt((vy - 50.) * (vy - 50.) + (vz + 24227. - Z_PX56) * (vz + 24227. - Z_PX56)) < 3550. && vy > -1000.)
0103     return Air;
0104 
0105   return Unknown;
0106 }
0107 
0108 inline int inWallAfterAir(double vx, double vy, double vz) {
0109   // phase II surface building
0110   if (vy < SurfaceOfEarth && vy >= (SurfaceOfEarth - 2250.)) {
0111     if (std::fabs(vz - Z_PX56) < 30000. && std::fabs(vx) < 10950)
0112       return Wall;
0113     // foundation of crane
0114     if (std::fabs(vz - Z_PX56) < 9000. && std::fabs(vx) >= 10950 && std::fabs(vx) < 16950)
0115       return Wall;
0116   }
0117 
0118   // CMS cavern (UXC 55)
0119   if (std::fabs(vz) < 29278. && sqrt((vx * 1.1576) * (vx * 1.1576) + vy * vy) < 16830. && vy > -11762.)
0120     return Wall;
0121 
0122   // access shaft (PX 56)
0123   if (vy > 0. && vy < (SurfaceOfEarth - 2250.) &&  //t(shaft wall)=2150.
0124       sqrt(vx * vx + (vz - Z_PX56) * (vz - Z_PX56)) < 12400.)
0125     return Wall;
0126 
0127   // Shaft (PM 54)
0128   if (vy > 3233. && vy < (SurfaceOfEarth - 1000.) &&  //t~=t(PX56)/R(PX56)*R(PM54)
0129       sqrt((vx - 26600.) * (vx - 26600.) + (vz - 30100. - Z_PX56) * (vz - 30100. - Z_PX56)) <
0130           6050. + 2150. / 10250. * 6050.)
0131     return Wall;
0132   //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56))//@@@@@
0133   //< 10050.+2150./10250.*6050.) return Wall;//@@@@@
0134   else if (vy >= SurfaceOfEarth - 1000. && vy < SurfaceOfEarth &&  //t~=t(PX56)/R(PX56)*R(PM54)
0135            sqrt((vx - 26600.) * (vx - 26600.) + (vz - 30100. - Z_PX56) * (vz - 30100. - Z_PX56)) <
0136                6050. + 2150. / 10250. * 6050. + 1800.)
0137     return Wall;
0138   //sqrt((vx-5000.)*(vx-5000.) + (vz-30100.-Z_PX56)*(vz-30100.-Z_PX56))//@@@@@
0139   //< 10050.+2150./10250.*10050. +1800.) return Wall; //@@@@@@
0140 
0141   // Shaft (PM 56)
0142   if (vy > -5450. && vy < (SurfaceOfEarth - 1000.) &&  //t~=t(PX56)/R(PX56)*R(PM56)
0143       sqrt((vx - 18220.) * (vx - 18220.) + (vz + 24227. - Z_PX56) * (vz + 24227. - Z_PX56)) <
0144           3550. + 2150. / 10250. * 3550.)
0145     return Wall;
0146   else if (vy > SurfaceOfEarth - 1000. && vy < SurfaceOfEarth &&  //t~=t(PX56)/R(PX56)*R(PM56)
0147            sqrt((vx - 18220.) * (vx - 18220.) + (vz + 24227. - Z_PX56) * (vz + 24227. - Z_PX56)) <
0148                3550. + 2150. / 10250. * 3550. + 1800.)
0149     return Wall;
0150 
0151   //Service cavern USC 55
0152   if (vz > -(22050. + 1150.) && vz < (62050. + 1150.) &&
0153       sqrt((vx - 29550.) * (vx - 29550.) + (vy - 3233.) * (vy - 3233.)) < 9050. + 950. && vy > -3650. - 2000.)
0154     return Wall;  //8762=estimate, to be checked
0155 
0156   //Pillar between CMS collision cavern and service cavern
0157   if (vz > -29278. + 1000. && vz < 29278. + 1000.) {
0158     if (vy > -17985. && vy < 10410. && vx > 13300. && vx < 20500.)
0159       return Wall;
0160 
0161     if (vy > 0. && vy < 10410. && vx > 10000. && vx <= 13300.)
0162       return Wall;
0163 
0164     if (vy > -3650. - 2000. && vy < -3233. &&  // bottom edge between pillar and service cavern
0165         vx > 20000. && vx < 24000.)
0166       return Wall;
0167     if (vy > -11762. && vy < -5000. &&  // bottom edge between pillar and UXC55 cavern
0168         vx > 10500. && vx < 14000.)
0169       return Wall;
0170   }
0171 
0172   if (vy > -14000. && vy < -1450. &&  //TX54 galerie surrounding
0173       vz > 13460. && vz < 17260. && vx >= 20500. && vx < 24550.)
0174     return Wall;
0175 
0176   //R56, LHC beam tunnel East
0177   if (vz > -85000. && vz < -28510.) {  //UJ57 junction to UXC55 cavern
0178     if (sqrt(vy * vy + (vx - 350.) * (vx - 350.)) < 2250.)
0179       return Wall;
0180   }
0181 
0182   //R54, LHC beam tunnel West
0183   if (vz > 26550. && vz < 63000. &&  //UJ57 junction to UXC55 cavern
0184       sqrt(vy * vy + (vx - 350.) * (vx - 350.)) < 2250.)
0185     return Wall;
0186 
0187   //UJ56 cavern
0188   if (vz > -(58875. + 500.) && vz < -(33927. - 500.) && sqrt(vy * vy + (vx - 4450.) * (vx - 4450.)) < (6750. + 500.) &&
0189       vy > -3650.)
0190     return Wall;
0191 
0192   //connection between PM56 shaft and UJ56 cavern
0193   if (vx > 9000. && vx < 18220. &&
0194       sqrt((vy - 50.) * (vy - 50.) + (vz + 24227. - Z_PX56) * (vz + 24227. - Z_PX56)) < 3550. + 500. && vy > -3650.)
0195     return Wall;
0196 
0197   return Unknown;
0198 }
0199 
0200 inline int inClayOrRockAfterWall(double vx, double vy, double vz, double ClayWidth) {
0201   //So, it is not plug, air and wall, Check for clay
0202   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
0203     return Clay;
0204 
0205   //So, it is not plug, air, wall and clay, Check for rock
0206   if (vy < SurfaceOfEarth - ClayWidth)
0207     return Rock;
0208 
0209   return Unknown;
0210 }
0211 
0212 inline int inClayAfterWall(double vx, double vy, double vz, double ClayWidth) {
0213   //So, it is not plug, air and wall, Check for clay
0214   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
0215     return Clay;
0216 
0217   return Unknown;
0218 }
0219 
0220 inline int inRockAfterClay(double vx, double vy, double vz, double ClayWidth) {
0221   //So, it is not plug, air, wall and clay, Check for rock
0222   if (vy < SurfaceOfEarth - ClayWidth)
0223     return Rock;
0224 
0225   return Unknown;
0226 }
0227 
0228 inline int inMat(double vx,
0229                  double vy,
0230                  double vz,
0231                  double PlugVx = PlugOnShaftVx,
0232                  double PlugVz = PlugOnShaftVz,
0233                  double ClayWidth = DefaultClayWidth) {
0234   //check for Plug
0235   if (inPlug(vx, vy, vz, PlugVx, PlugVz))
0236     return Plug;
0237 
0238   //So, it is not plug, Check for air
0239   if (inAirAfterPlug(vx, vy, vz))
0240     return Air;
0241 
0242   //So, it is not plug and air, Check for wall
0243   if (inWallAfterAir(vx, vy, vz))
0244     return Wall;
0245 
0246   //So, it is not plug, air and wall, Check for clay
0247   if (vy >= SurfaceOfEarth - ClayWidth && vy < SurfaceOfEarth)
0248     return Clay;
0249 
0250   //So, it is not plug, air, wall and clay, Check for rock
0251   if (vy < SurfaceOfEarth - ClayWidth)
0252     return Rock;
0253 
0254   std::cout << "Point5MaterialMap.h: Warning! No Material recognised for point: vx=" << vx << " vy=" << vy
0255             << " vz=" << vz << std::endl;
0256   //Something went wrong
0257   return Unknown;
0258 }
0259 
0260 #endif