CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_13_patch3/src/Geometry/HcalTowerAlgo/src/CaloTowerHardcodeGeometryLoader.cc

Go to the documentation of this file.
00001 #include "Geometry/HcalTowerAlgo/interface/CaloTowerHardcodeGeometryLoader.h"
00002 #include "Geometry/HcalTowerAlgo/src/HcalHardcodeGeometryData.h"
00003 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
00004 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
00005 #include "Geometry/CaloGeometry/interface/IdealObliquePrism.h"
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 
00008 typedef CaloCellGeometry::CCGFloat CCGFloat ;
00009 
00010 std::auto_ptr<CaloSubdetectorGeometry> CaloTowerHardcodeGeometryLoader::load() {
00011   CaloTowerGeometry* geom=new CaloTowerGeometry();
00012 
00013   if( 0 == geom->cornersMgr() ) geom->allocateCorners ( 
00014      CaloTowerGeometry::k_NumberOfCellsForCorners ) ;
00015   if( 0 == geom->parMgr() ) geom->allocatePar (
00016      CaloTowerGeometry::k_NumberOfParametersPerShape*CaloTowerGeometry::k_NumberOfShapes,
00017      CaloTowerGeometry::k_NumberOfParametersPerShape ) ;
00018 
00019   int nnn=0;
00020   // simple loop
00021   for (int ieta=-limits.lastHFRing(); ieta<=limits.lastHFRing(); ieta++) {
00022     if (ieta==0) continue; // skip not existing eta=0 ring
00023     for (int iphi=1; iphi<=72; iphi++) {
00024       if (abs(ieta)>=limits.firstHFQuadPhiRing() && ((iphi-1)%4)==0) continue;
00025       if (abs(ieta)>=limits.firstHEDoublePhiRing() && ((iphi-1)%2)!=0) continue;
00026       ++nnn;
00027     }
00028   }
00029   if( geom->cornersMgr() == 0 ) geom->allocateCorners( nnn ) ; 
00030   if( geom->parMgr()     == 0 ) geom->allocatePar( 41, 3 ) ;
00031 
00032   int n=0;
00033   // simple loop
00034   for (int ieta=-limits.lastHFRing(); ieta<=limits.lastHFRing(); ieta++) {
00035     if (ieta==0) continue; // skip not existing eta=0 ring
00036     for (int iphi=1; iphi<=72; iphi++) {
00037       if (abs(ieta)>=limits.firstHFQuadPhiRing() && ((iphi-1)%4)==0) continue;
00038       if (abs(ieta)>=limits.firstHEDoublePhiRing() && ((iphi-1)%2)!=0) continue;
00039       makeCell(ieta,iphi, geom);
00040       n++;
00041     }
00042   }
00043   edm::LogInfo("Geometry") << "CaloTowersHardcodeGeometry made " << n << " towers.";
00044   return std::auto_ptr<CaloSubdetectorGeometry>(geom); 
00045 }
00046 
00047 void
00048 CaloTowerHardcodeGeometryLoader::makeCell( int ieta,
00049                                            int iphi,
00050                                            CaloSubdetectorGeometry* geom ) const {
00051   const double EBradius = 143.0; // cm
00052   const double HOradius = 406.0+1.0;
00053   const double EEz = 320.0; // rough (cm)
00054   const double HEz = 568.0; // back (cm)
00055   const double HFz = 1100.0;
00056   const double HFthick = 165;
00057   // Tower 17 is the last EB tower
00058 
00059   int etaRing=abs(ieta);
00060   int sign=(ieta>0)?(1):(-1);
00061   double eta1, eta2;
00062   if (etaRing>limits.lastHERing()) {
00063     eta1 = theHFEtaBounds[etaRing-limits.firstHFRing()];
00064     eta2 = theHFEtaBounds[etaRing-limits.firstHFRing()+1];
00065   } else {
00066     eta1 = theHBHEEtaBounds[etaRing-1];
00067     eta2 = theHBHEEtaBounds[etaRing];
00068   }
00069   double eta = 0.5*(eta1+eta2);
00070   double deta = (eta2-eta1);  
00071 
00072   // in radians
00073   double dphi_nominal = 2.0*M_PI / limits.nPhiBins(1); // always the same
00074   double dphi_half = M_PI / limits.nPhiBins(etaRing); // half-width
00075   
00076   double phi_low = dphi_nominal*(iphi-1); // low-edge boundaries are constant...
00077   double phi = phi_low+dphi_half;
00078 
00079   double x,y,z,thickness;
00080   bool alongZ=true;
00081   if (etaRing>limits.lastHERing()) { // forward
00082     z=HFz;
00083     double r=z/sinh(eta);
00084     x=r * cos(phi);
00085     y=r * sin(phi);
00086     thickness=HFthick/tanh(eta);
00087   } else if (etaRing>17) { // EE-containing
00088     z=EEz;
00089     double r=z/sinh(eta);
00090     x=r * cos(phi);
00091     y=r * sin(phi);
00092     thickness=(HEz-EEz)/tanh(eta);
00093   } else { // EB-containing
00094     x=EBradius * cos(phi);
00095     y=EBradius * sin(phi);
00096     alongZ=false;
00097     z=EBradius * sinh(eta);
00098     thickness=(HOradius-EBradius) * cosh(eta);
00099   }
00100 
00101   z*=sign;
00102   GlobalPoint point(x,y,z);
00103 
00104   const double mysign ( !alongZ ? 1 : -1 ) ;
00105   std::vector<CCGFloat> hh ;
00106   hh.reserve(5) ;
00107   hh.push_back( deta/2 ) ;
00108   hh.push_back( dphi_half ) ;
00109   hh.push_back( mysign*thickness/2. ) ;
00110 
00111   hh.push_back( fabs( eta ) ) ;
00112   hh.push_back( fabs( z ) ) ;
00113 
00114   geom->newCell( point, point, point,
00115                  CaloCellGeometry::getParmPtr( hh, 
00116                                                geom->parMgr(), 
00117                                                geom->parVecVec() ),
00118                  CaloTowerDetId( ieta, iphi ) ) ;
00119 }