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 std::auto_ptr<CaloSubdetectorGeometry> CaloTowerHardcodeGeometryLoader::load() {
00009 CaloTowerGeometry* geom=new CaloTowerGeometry();
00010
00011 if( 0 == geom->cornersMgr() ) geom->allocateCorners (
00012 CaloTowerGeometry::k_NumberOfCellsForCorners ) ;
00013 if( 0 == geom->parMgr() ) geom->allocatePar (
00014 CaloTowerGeometry::k_NumberOfParametersPerShape*CaloTowerGeometry::k_NumberOfShapes,
00015 CaloTowerGeometry::k_NumberOfParametersPerShape ) ;
00016
00017 int nnn=0;
00018
00019 for (int ieta=-limits.lastHFRing(); ieta<=limits.lastHFRing(); ieta++) {
00020 if (ieta==0) continue;
00021 for (int iphi=1; iphi<=72; iphi++) {
00022 if (abs(ieta)>=limits.firstHFQuadPhiRing() && ((iphi-1)%4)==0) continue;
00023 if (abs(ieta)>=limits.firstHEDoublePhiRing() && ((iphi-1)%2)!=0) continue;
00024 ++nnn;
00025 }
00026 }
00027 if( geom->cornersMgr() == 0 ) geom->allocateCorners( nnn ) ;
00028 if( geom->parMgr() == 0 ) geom->allocatePar( 41, 3 ) ;
00029
00030 int n=0;
00031
00032 for (int ieta=-limits.lastHFRing(); ieta<=limits.lastHFRing(); ieta++) {
00033 if (ieta==0) continue;
00034 for (int iphi=1; iphi<=72; iphi++) {
00035 if (abs(ieta)>=limits.firstHFQuadPhiRing() && ((iphi-1)%4)==0) continue;
00036 if (abs(ieta)>=limits.firstHEDoublePhiRing() && ((iphi-1)%2)!=0) continue;
00037 geom->addCell(CaloTowerDetId(ieta,iphi),makeCell(ieta,iphi, geom));
00038 n++;
00039 }
00040 }
00041 edm::LogInfo("Geometry") << "CaloTowersHardcodeGeometry made " << n << " towers.";
00042 return std::auto_ptr<CaloSubdetectorGeometry>(geom);
00043 }
00044
00045 CaloCellGeometry*
00046 CaloTowerHardcodeGeometryLoader::makeCell( int ieta,
00047 int iphi,
00048 CaloSubdetectorGeometry* geom ) const {
00049 const double EBradius = 143.0;
00050 const double HOradius = 406.0+1.0;
00051 const double EEz = 320.0;
00052 const double HEz = 568.0;
00053 const double HFz = 1100.0;
00054 const double HFthick = 165;
00055
00056
00057 int etaRing=abs(ieta);
00058 int sign=(ieta>0)?(1):(-1);
00059 double eta1, eta2;
00060 if (etaRing>limits.lastHERing()) {
00061 eta1 = theHFEtaBounds[etaRing-limits.firstHFRing()];
00062 eta2 = theHFEtaBounds[etaRing-limits.firstHFRing()+1];
00063 } else {
00064 eta1 = theHBHEEtaBounds[etaRing-1];
00065 eta2 = theHBHEEtaBounds[etaRing];
00066 }
00067 double eta = 0.5*(eta1+eta2);
00068 double deta = (eta2-eta1);
00069
00070
00071 double dphi_nominal = 2.0*M_PI / limits.nPhiBins(1);
00072 double dphi_half = M_PI / limits.nPhiBins(etaRing);
00073
00074 double phi_low = dphi_nominal*(iphi-1);
00075 double phi = phi_low+dphi_half;
00076
00077 double x,y,z,thickness;
00078 bool alongZ=true;
00079 if (etaRing>limits.lastHERing()) {
00080 z=HFz;
00081 double r=z/sinh(eta);
00082 x=r * cos(phi);
00083 y=r * sin(phi);
00084 thickness=HFthick/tanh(eta);
00085 } else if (etaRing>17) {
00086 z=EEz;
00087 double r=z/sinh(eta);
00088 x=r * cos(phi);
00089 y=r * sin(phi);
00090 thickness=(HEz-EEz)/tanh(eta);
00091 } else {
00092 x=EBradius * cos(phi);
00093 y=EBradius * sin(phi);
00094 alongZ=false;
00095 z=EBradius * sinh(eta);
00096 thickness=(HOradius-EBradius) * cosh(eta);
00097 }
00098
00099 z*=sign;
00100 GlobalPoint point(x,y,z);
00101
00102 const double mysign ( !alongZ ? 1 : -1 ) ;
00103 std::vector<double> hh ;
00104 hh.reserve(5) ;
00105 hh.push_back( deta/2 ) ;
00106 hh.push_back( dphi_half ) ;
00107 hh.push_back( mysign*thickness/2. ) ;
00108
00109 hh.push_back( fabs( eta ) ) ;
00110 hh.push_back( fabs( z ) ) ;
00111
00112 return new calogeom::IdealObliquePrism(
00113 point,
00114 geom->cornersMgr(),
00115 CaloCellGeometry::getParmPtr( hh,
00116 geom->parMgr(),
00117 geom->parVecVec() ) );
00118 }