Reference Manual

Go to the documentation of this file.
00001 #include "Geometry/HcalTowerAlgo/interface/HcalFlexiHardcodeGeometryLoader.h"
00002 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
00003 #include "Geometry/CaloGeometry/interface/IdealObliquePrism.h"
00004 #include "Geometry/CaloGeometry/interface/IdealZPrism.h"
00005 #include "Geometry/HcalTowerAlgo/interface/HcalGeometry.h"
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 
00008 #include <vector>
00009 
00010 typedef CaloCellGeometry::CCGFloat CCGFloat ;
00011 //#define DebugLog
00012 // ==============> Loader Itself <==========================
00013 
00014 HcalFlexiHardcodeGeometryLoader::HcalFlexiHardcodeGeometryLoader(const edm::ParameterSet& ps) {
00015 
00016   MAX_HCAL_PHI = 72;
00017   DEGREE2RAD = M_PI / 180.;
00018 
00019   edm::ParameterSet ps0 = ps.getParameter<edm::ParameterSet>("HcalReLabel");
00020   bool relabel_= ps0.getUntrackedParameter<bool>("RelabelHits",false);
00021   if (relabel_) {
00022     edm::ParameterSet ps1 = ps0.getUntrackedParameter<edm::ParameterSet>("RelabelRules");
00023     m_segmentation.resize(29);
00024     for (int i=0; i<29; i++) {
00025       char name[10];
00026       snprintf(name,10,"Eta%d",i+1);
00027       if (i>0) {
00028         m_segmentation[i]=ps1.getUntrackedParameter<std::vector<int> >(name,m_segmentation[i-1]);
00029       } else {
00030         m_segmentation[i]=ps1.getUntrackedParameter<std::vector<int> >(name);
00031       }
00032 #ifdef DebugLog
00033       std::cout << name;
00034       for (unsigned int k=0; k<m_segmentation[i].size(); ++k) {
00035         std::cout << " [" << k << "] " << m_segmentation[i][k];
00036       }
00037       std::cout << std::endl;
00038 #endif
00039     }
00040   }
00041 
00042 }
00043 
00044 CaloSubdetectorGeometry* HcalFlexiHardcodeGeometryLoader::load(const HcalTopology& fTopology) {
00045   CaloSubdetectorGeometry* hcalGeometry = new HcalGeometry (fTopology);
00046   if( 0 == hcalGeometry->cornersMgr() ) hcalGeometry->allocateCorners ( fTopology.ncells() );
00047   if( 0 == hcalGeometry->parMgr() ) hcalGeometry->allocatePar (hcalGeometry->numberOfShapes(),
00048                                                                HcalGeometry::k_NumberOfParametersPerShape ) ;
00049   if (fTopology.mode() == HcalTopologyMode::H2) {  // TB geometry
00050     fillHBHO (hcalGeometry, makeHBCells(fTopology), true);
00051     fillHBHO (hcalGeometry, makeHOCells(), false);
00052     fillHE (hcalGeometry, makeHECells_H2());
00053  } else { // regular geometry
00054     fillHBHO (hcalGeometry, makeHBCells(fTopology), true);
00055     fillHBHO (hcalGeometry, makeHOCells(), false);
00056     fillHF (hcalGeometry, makeHFCells());
00057     fillHE (hcalGeometry, makeHECells(fTopology));
00058   }
00059   return hcalGeometry;
00060 }
00061 
00062 
00063 // ----------> HB <-----------
00064 std::vector <HcalFlexiHardcodeGeometryLoader::HBHOCellParameters> HcalFlexiHardcodeGeometryLoader::makeHBCells (const HcalTopology & topology) {
00065 
00066   const float HBRMIN = 181.1;
00067   const float HBRMAX = 288.8;
00068     
00069   float normalDepths[2] = {HBRMIN, HBRMAX};
00070   float ring15Depths[3] = {HBRMIN, 258.4, HBRMAX};
00071   float ring16Depths[3] = {HBRMIN, 190.4, 232.6};
00072   float layerDepths[18] = {HBRMIN, 188.7, 194.7, 200.7, 206.7, 212.7, 218.7,
00073                            224.7, 230.7, 236.7, 242.7, 249.3, 255.9, 262.5,
00074                            269.1, 275.7, 282.3, HBRMAX};
00075   float slhcDepths[4]   = {HBRMIN, 214., 239., HBRMAX};
00076 #ifdef DebugLog
00077   std::cout <<"FlexiGeometryLoader called for "<< topology.mode() << ":" << HcalTopologyMode::SLHC << std::endl;
00078 #endif
00079   std::vector <HcalFlexiHardcodeGeometryLoader::HBHOCellParameters> result;
00080   for(int iring = 1; iring <= 16; ++iring) {
00081     std::vector<float> depths;
00082     if (topology.mode() != HcalTopologyMode::SLHC) {
00083       if (iring == 15) {
00084         for (int i=0; i<3; ++i) depths.push_back(ring15Depths[i]);
00085       } else if (iring == 16) {
00086         for (int i=0; i<3; ++i) depths.push_back(ring16Depths[i]);
00087       } else {
00088         for (int i=0; i<2; ++i) depths.push_back(normalDepths[i]);
00089       }
00090     } else {
00091       if (m_segmentation.size() >= (unsigned int)(iring)) {
00092         int depth = m_segmentation[iring-1][0];
00093         depths.push_back(layerDepths[depth]);
00094         int layer = 1;
00095         for (unsigned int i=1; i<m_segmentation[iring-1].size(); ++i) {
00096           if (depth != m_segmentation[iring-1][i]) {
00097             depth = m_segmentation[iring-1][i];
00098             layer = i;
00099             if (iring != 16 || depth < 3)
00100               depths.push_back(layerDepths[depth]);
00101           }
00102           if (i >= 17) break;
00103         }
00104         if (layer <= 17) depths.push_back(HBRMAX);
00105       } else {
00106         for (int i=0; i<4; ++i) {
00107           if (iring != 16 || i < 3) {
00108             depths.push_back(slhcDepths[i]);
00109           }
00110         }
00111       }
00112     }
00113     unsigned int ndepth=depths.size()-1;
00114     unsigned int startingDepth=1;
00115     float etaMin=(iring-1)*0.087;
00116     float etaMax=iring*0.087;
00117     // topology.depthBinInformation(HcalBarrel, iring, ndepth, startingDepth);
00118 #ifdef DebugLog
00119     std::cout << "HBRing " << iring << " eta " << etaMin << ":" << etaMax << " depths " << ndepth << ":" << startingDepth;
00120     for (unsigned int i=0; i<depths.size(); ++i) std::cout << ":" << depths[i];
00121     std::cout << "\n";
00122 #endif
00123     for (unsigned int idepth = startingDepth; idepth <= ndepth; ++idepth) {
00124       float rmin = depths[idepth-1];
00125       float rmax = depths[idepth];
00126 #ifdef DebugLog
00127       std::cout << "HB " << idepth << " R " << rmin << ":" << rmax << "\n";
00128 #endif
00129       result.push_back(HcalFlexiHardcodeGeometryLoader::HBHOCellParameters(iring, (int)idepth, 1, 1, 5, rmin, rmax, etaMin, etaMax));
00130     }
00131   }
00132   return result;
00133 }
00134 
00135 
00136 
00137 // ----------> HO <-----------
00138 std::vector <HcalFlexiHardcodeGeometryLoader::HBHOCellParameters> HcalFlexiHardcodeGeometryLoader::makeHOCells () {
00139   const float HORMIN0 = 390.0;
00140   const float HORMIN1 = 412.6;
00141   const float HORMAX = 413.6;
00142   
00143   HcalFlexiHardcodeGeometryLoader::HBHOCellParameters cells [] = {
00144     // eta, depth, firstPhi, stepPhi, deltaPhi, rMin, rMax, etaMin, etaMax
00145     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 1, 4, 1, 1, 5, HORMIN0, HORMAX, 0.087*0, 0.087*1),
00146     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 2, 4, 1, 1, 5, HORMIN0, HORMAX, 0.087*1, 0.087*2),
00147     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 3, 4, 1, 1, 5, HORMIN0, HORMAX, 0.087*2, 0.087*3),
00148     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 4, 4, 1, 1, 5, HORMIN0, HORMAX, 0.087*3, 0.3075),
00149     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 5, 4, 1, 1, 5, HORMIN1, HORMAX, 0.3395,  0.087*5),
00150     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 6, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*5, 0.087*6),
00151     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 7, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*6, 0.087*7),
00152     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 8, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*7, 0.087*8),
00153     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters ( 9, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*8, 0.087*9),
00154     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (10, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*9,  0.8494),
00155     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (11, 4, 1, 1, 5, HORMIN1, HORMAX, 0.873, 0.087*11),
00156     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (12, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*11, 0.087*12),
00157     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (13, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*12, 0.087*13),
00158     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (14, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*13, 0.087*14),
00159     HcalFlexiHardcodeGeometryLoader::HBHOCellParameters (15, 4, 1, 1, 5, HORMIN1, HORMAX, 0.087*14, 0.087*15)
00160   };
00161   int nCells = sizeof(cells)/sizeof(HcalFlexiHardcodeGeometryLoader::HBHOCellParameters);
00162   std::vector <HcalFlexiHardcodeGeometryLoader::HBHOCellParameters> result;
00163   result.reserve (nCells);
00164   for (int i = 0; i < nCells; ++i) result.push_back (cells[i]);
00165   return result;
00166 }
00167 
00168 
00169 //
00170 // Convert constants to appropriate cells
00171 //
00172 void HcalFlexiHardcodeGeometryLoader::fillHBHO (CaloSubdetectorGeometry* fGeometry, const std::vector <HcalFlexiHardcodeGeometryLoader::HBHOCellParameters>& fCells, bool fHB) {
00173 
00174   for (size_t iCell = 0; iCell < fCells.size(); ++iCell) {
00175     const HcalFlexiHardcodeGeometryLoader::HBHOCellParameters& param = fCells[iCell];
00176     for (int iPhi = param.phiFirst; iPhi <= MAX_HCAL_PHI; iPhi += param.phiStep) {
00177       for (int iside = -1; iside <= 1; iside += 2) { // both detector sides are identical
00178         HcalDetId hid (fHB ? HcalBarrel : HcalOuter, param.eta*iside, iPhi, param.depth);
00179         float phiCenter = ((iPhi-1)*360./MAX_HCAL_PHI + 0.5*param.dphi) * DEGREE2RAD; // middle of the cell
00180         float etaCenter = 0.5*(param.etaMin + param.etaMax);
00181         float x = param.rMin* cos (phiCenter);
00182         float y = param.rMin* sin (phiCenter);
00183         float z = iside * param.rMin * sinh(etaCenter);
00184         // make cell geometry
00185         GlobalPoint refPoint (x,y,z); // center of the cell's face
00186         std::vector<CCGFloat> cellParams;
00187         cellParams.reserve (5);
00188         cellParams.push_back (0.5 * (param.etaMax - param.etaMin)); // deta_half
00189         cellParams.push_back (0.5 * param.dphi * DEGREE2RAD);  // dphi_half
00190         cellParams.push_back (0.5 * (param.rMax - param.rMin) * cosh (etaCenter)); // dr_half
00191         cellParams.push_back ( fabs( refPoint.eta() ) ) ;
00192         cellParams.push_back ( fabs( refPoint.z() ) ) ;
00193 #ifdef DebugLog
00194         std::cout << "HcalFlexiHardcodeGeometryLoader::fillHBHO-> " << hid << hid.ieta() << '/' << hid.iphi() << '/' << hid.depth() << refPoint << '/' << cellParams [0] << '/' << cellParams [1] << '/' << cellParams [2] << std::endl;
00195 #endif
00196         fGeometry->newCell(refPoint,  refPoint,  refPoint, 
00197                            CaloCellGeometry::getParmPtr(cellParams, 
00198                                                         fGeometry->parMgr(), 
00199                                                         fGeometry->parVecVec()),
00200                            hid ) ;
00201       }
00202     }
00203   }
00204 }
00205 
00206 
00207 // ----------> HE <-----------
00208 std::vector<HcalFlexiHardcodeGeometryLoader::HECellParameters> HcalFlexiHardcodeGeometryLoader::makeHECells (const HcalTopology & topology) {
00209 
00210   std::vector <HcalFlexiHardcodeGeometryLoader::HECellParameters> result;
00211   const float HEZMIN = 400.458;
00212   const float HEZMID = 436.168;
00213   const float HEZMAX = 549.268;
00214   float normalDepths[3] = {HEZMIN, HEZMID, HEZMAX};
00215   float tripleDepths[4] = {HEZMIN, 418.768, HEZMID, HEZMAX};
00216   float slhcDepths[5]   = {HEZMIN, 418.768, HEZMID, 493., HEZMAX};
00217   float ring16Depths[2] = {418.768,470.968};
00218   float ring16slhcDepths[3] = {418.768, 450., 470.968};
00219   float ring17Depths[2] = {409.698,514.468};
00220   float ring17slhcDepths[5] = {409.698, 435., 460., 495., 514.468};
00221   float ring18Depths[3] = {391.883,427.468,540.568};
00222   float ring18slhcDepths[5] = {391.883, 439.,  467., 504. , 540.568};
00223   float etaBounds[] = {0.087*15, 0.087*16, 0.087*17, 0.087*18,  0.087*19,
00224                        1.74, 1.83,  1.93, 2.043, 2.172, 2.322, 2.500,
00225                        2.650, 2.868, 3.000};
00226   float layerDepths[19] = {HEZMIN, 408.718, 416.978, 425.248, 433.508, 441.768,
00227                            450.038,458.298, 466.558, 474.828, 483.088, 491.348,
00228                            499.618,507.878, 516.138, 524.398, 532.668, 540.928,
00229                            HEZMAX};
00230 
00231   // count by ring - 16
00232   for(int iringm16=0; iringm16 <= 13; ++iringm16) {
00233     int iring = iringm16 + 16;
00234     std::vector<float> depths;
00235     unsigned int startingDepth = 1;
00236     if (topology.mode() != HcalTopologyMode::SLHC) {
00237       if (iring == 16)     
00238         {for (int i=0; i<2; ++i) depths.push_back(ring16Depths[i]); startingDepth = 3;}
00239       else if (iring == 17) 
00240         for (int i=0; i<2; ++i) depths.push_back(ring17Depths[i]);
00241       else if (iring == 18) 
00242         for (int i=0; i<3; ++i) depths.push_back(ring18Depths[i]);
00243       else if (iring == topology.lastHERing()) 
00244         for (int i=0; i<3; ++i) depths.push_back(tripleDepths[i]);
00245       else if (iring >= topology.firstHETripleDepthRing())
00246         for (int i=0; i<4; ++i) depths.push_back(tripleDepths[i]);
00247       else
00248         for (int i=0; i<3; ++i) depths.push_back(normalDepths[i]);
00249     } else {
00250       if (m_segmentation.size() >= (unsigned int)(iring)) {
00251         int depth = m_segmentation[iring-1][0];
00252         if (iring == 16)      depths.push_back(ring16Depths[0]);
00253         else if (iring == 17) depths.push_back(ring17Depths[0]);
00254         else if (iring == 18) depths.push_back(ring18Depths[0]);
00255         else                  depths.push_back(layerDepths[depth]);
00256         int layer = 1;
00257         float lastDepth = depths[0];
00258         for (unsigned int i=1; i<m_segmentation[iring-1].size(); ++i) {
00259           if (depth != m_segmentation[iring-1][i]) {
00260             depth = m_segmentation[iring-1][i];
00261             layer = i;
00262             if (layerDepths[depth] > lastDepth && (iring != 16 || depth > 3)) {
00263               depths.push_back(layerDepths[depth]);
00264               lastDepth = layerDepths[depth];
00265             }
00266           }
00267         }
00268         if (layer <= 17) depths.push_back(HEZMAX);
00269         if (iring == 16) startingDepth = 3;
00270       } else {
00271         if (iring == 16)     {for (int i=0; i<3; ++i) depths.push_back(ring16slhcDepths[i]); startingDepth = 3;}
00272         else if (iring == 17) for (int i=0; i<5; ++i) depths.push_back(ring17slhcDepths[i]);
00273         else if (iring == 18) for (int i=0; i<5; ++i) depths.push_back(ring18slhcDepths[i]);
00274         else                  for (int i=0; i<5; ++i) depths.push_back(slhcDepths[i]);
00275       }
00276     }
00277     float etamin = etaBounds[iringm16];
00278     float etamax = etaBounds[iringm16+1];
00279     unsigned int ndepth = depths.size()-1;
00280     //    topology.depthBinInformation(HcalEndcap, iring, ndepth, startingDepth);
00281 #ifdef DebugLog
00282     std::cout << "HERing " << iring << " eta " << etamin << ":" << etamax << " depths " << ndepth << ":" << startingDepth;
00283     for (unsigned int i=0; i<depths.size(); ++i) std::cout << ":" << depths[i];
00284     std::cout << "\n";
00285 #endif
00286     for (unsigned int idepth = 0; idepth < ndepth; ++idepth) {
00287       int depthIndex = (int)(idepth + startingDepth);
00288       float zmin = depths[idepth];
00289       float zmax = depths[idepth+1];
00290       if (depthIndex <= 7) {
00291 #ifdef DebugLog
00292         std::cout << "HE Depth " << idepth << ":" << depthIndex << " Z " << zmin << ":" << zmax << "\n";
00293 #endif
00294         int stepPhi = (iring >= topology.firstHEDoublePhiRing() ? 2 : 1);
00295         int deltaPhi =  (iring >= topology.firstHEDoublePhiRing() ? 10 : 5);
00296         if (topology.mode() != HcalTopologyMode::SLHC &&
00297             iring == topology.lastHERing()-1 && idepth == ndepth-1) {
00298 #ifdef DebugLog
00299           std::cout << "HE iEta " << iring << " Depth " << depthIndex << " Eta " << etamin << ":" << etaBounds[iringm16+2] << std::endl;
00300 #endif
00301           result.push_back(HcalFlexiHardcodeGeometryLoader::HECellParameters(iring, depthIndex, 1, stepPhi, deltaPhi, zmin, zmax, etamin, etaBounds[iringm16+2]));
00302         } else {
00303 #ifdef DebugLog
00304           std::cout << "HE iEta " << iring << " Depth " << depthIndex << " Eta " << etamin << ":" << etamax << std::endl;
00305 #endif
00306           result.push_back(HcalFlexiHardcodeGeometryLoader::HECellParameters(iring, depthIndex, 1, stepPhi, deltaPhi, zmin, zmax, etamin, etamax));
00307         }
00308       }
00309     }
00310   }
00311 
00312   return result;
00313 }
00314 
00315 
00316 // ----------> HE @ H2 <-----------
00317 std::vector <HcalFlexiHardcodeGeometryLoader::HECellParameters> HcalFlexiHardcodeGeometryLoader::makeHECells_H2 () {
00318 
00319   const float HEZMIN_H2 = 400.715;
00320   const float HEZMID_H2 = 436.285;
00321   const float HEZMAX_H2 = 541.885;
00322     
00323   HcalFlexiHardcodeGeometryLoader::HECellParameters cells [] = {
00324     // eta, depth, firstPhi, stepPhi, deltaPhi, zMin, zMax, etaMin, etaMax
00325     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 16, 3, 1, 1, 5, 409.885,   462.685,   1.305, 1.373),
00326     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 17, 1, 1, 1, 5, HEZMIN_H2, 427.485,   1.373, 1.444),
00327     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 17, 2, 1, 1, 5, 427.485,   506.685,   1.373, 1.444),
00328     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 18, 1, 1, 1, 5, HEZMIN_H2, HEZMID_H2, 1.444, 1.521),
00329     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 18, 2, 1, 1, 5, HEZMID_H2, 524.285,   1.444, 1.521),
00330     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 19, 1, 1, 1, 5, HEZMIN_H2, HEZMID_H2, 1.521, 1.603),
00331     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 19, 2, 1, 1, 5, HEZMID_H2, HEZMAX_H2, 1.521, 1.603),
00332     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 20, 1, 1, 1, 5, HEZMIN_H2, HEZMID_H2, 1.603, 1.693),
00333     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 20, 2, 1, 1, 5, HEZMID_H2, HEZMAX_H2, 1.603, 1.693),
00334     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 21, 1, 1, 2, 5, HEZMIN_H2, HEZMID_H2, 1.693, 1.79),
00335     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 21, 2, 1, 2, 5, HEZMID_H2, HEZMAX_H2, 1.693, 1.79),
00336     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 22, 1, 1, 2,10, HEZMIN_H2, HEZMID_H2, 1.79, 1.88),
00337     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 22, 2, 1, 2,10, HEZMID_H2, HEZMAX_H2, 1.79, 1.88),
00338     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 23, 1, 1, 2,10, HEZMIN_H2, HEZMID_H2, 1.88, 1.98),
00339     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 23, 2, 1, 2,10, HEZMID_H2, HEZMAX_H2, 1.88, 1.98),
00340     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 24, 1, 1, 2,10, HEZMIN_H2, 418.685,   1.98, 2.09),
00341     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 24, 2, 1, 2,10, 418.685,   HEZMID_H2, 1.98, 2.09),
00342     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 24, 3, 1, 2,10, HEZMID_H2, HEZMAX_H2, 1.98, 2.09),
00343     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 25, 1, 1, 2,10, HEZMIN_H2, 418.685,   2.09, 2.21),
00344     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 25, 2, 1, 2,10, 418.685,   HEZMID_H2, 2.09, 2.21),
00345     HcalFlexiHardcodeGeometryLoader::HECellParameters ( 25, 3, 1, 2,10, HEZMID_H2, HEZMAX_H2, 2.09, 2.21)
00346   };
00347   int nCells = sizeof(cells)/sizeof(HcalFlexiHardcodeGeometryLoader::HECellParameters);
00348   std::vector <HcalFlexiHardcodeGeometryLoader::HECellParameters> result;
00349   result.reserve (nCells);
00350   for (int i = 0; i < nCells; ++i) result.push_back (cells[i]);
00351   return result;
00352 }
00353 
00354 // ----------> HF <-----------
00355 std::vector <HcalFlexiHardcodeGeometryLoader::HFCellParameters> HcalFlexiHardcodeGeometryLoader::makeHFCells () {
00356 
00357   const float HFZMIN1 = 1115.;
00358   const float HFZMIN2 = 1137.;
00359   const float HFZMAX = 1280.1;
00360     
00361   HcalFlexiHardcodeGeometryLoader::HFCellParameters cells [] = {
00362     // eta, depth, firstPhi, stepPhi, deltaPhi, zMin, zMax, rMin, rMax
00363     HcalFlexiHardcodeGeometryLoader::HFCellParameters (29, 1, 1, 2, 10, HFZMIN1, HFZMAX,116.2,130.0),
00364     HcalFlexiHardcodeGeometryLoader::HFCellParameters (29, 2, 1, 2, 10, HFZMIN2, HFZMAX,116.2,130.0),
00365     HcalFlexiHardcodeGeometryLoader::HFCellParameters (30, 1, 1, 2, 10, HFZMIN1, HFZMAX, 97.5,116.2),
00366     HcalFlexiHardcodeGeometryLoader::HFCellParameters (30, 2, 1, 2, 10, HFZMIN2, HFZMAX, 97.5,116.2),
00367     HcalFlexiHardcodeGeometryLoader::HFCellParameters (31, 1, 1, 2, 10, HFZMIN1, HFZMAX, 81.8, 97.5),
00368     HcalFlexiHardcodeGeometryLoader::HFCellParameters (31, 2, 1, 2, 10, HFZMIN2, HFZMAX, 81.8, 97.5),
00369     HcalFlexiHardcodeGeometryLoader::HFCellParameters (32, 1, 1, 2, 10, HFZMIN1, HFZMAX, 68.6, 81.8),
00370     HcalFlexiHardcodeGeometryLoader::HFCellParameters (32, 2, 1, 2, 10, HFZMIN2, HFZMAX, 68.6, 81.8),
00371     HcalFlexiHardcodeGeometryLoader::HFCellParameters (33, 1, 1, 2, 10, HFZMIN1, HFZMAX, 57.6, 68.6),
00372     HcalFlexiHardcodeGeometryLoader::HFCellParameters (33, 2, 1, 2, 10, HFZMIN2, HFZMAX, 57.6, 68.6),
00373     HcalFlexiHardcodeGeometryLoader::HFCellParameters (34, 1, 1, 2, 10, HFZMIN1, HFZMAX, 48.3, 57.6),
00374     HcalFlexiHardcodeGeometryLoader::HFCellParameters (34, 2, 1, 2, 10, HFZMIN2, HFZMAX, 48.3, 57.6),
00375     HcalFlexiHardcodeGeometryLoader::HFCellParameters (35, 1, 1, 2, 10, HFZMIN1, HFZMAX, 40.6, 48.3),
00376     HcalFlexiHardcodeGeometryLoader::HFCellParameters (35, 2, 1, 2, 10, HFZMIN2, HFZMAX, 40.6, 48.3),
00377     HcalFlexiHardcodeGeometryLoader::HFCellParameters (36, 1, 1, 2, 10, HFZMIN1, HFZMAX, 34.0, 40.6),
00378     HcalFlexiHardcodeGeometryLoader::HFCellParameters (36, 2, 1, 2, 10, HFZMIN2, HFZMAX, 34.0, 40.6),
00379     HcalFlexiHardcodeGeometryLoader::HFCellParameters (37, 1, 1, 2, 10, HFZMIN1, HFZMAX, 28.6, 34.0),
00380     HcalFlexiHardcodeGeometryLoader::HFCellParameters (37, 2, 1, 2, 10, HFZMIN2, HFZMAX, 28.6, 34.0),
00381     HcalFlexiHardcodeGeometryLoader::HFCellParameters (38, 1, 1, 2, 10, HFZMIN1, HFZMAX, 24.0, 28.6),
00382     HcalFlexiHardcodeGeometryLoader::HFCellParameters (38, 2, 1, 2, 10, HFZMIN2, HFZMAX, 24.0, 28.6),
00383     HcalFlexiHardcodeGeometryLoader::HFCellParameters (39, 1, 1, 2, 10, HFZMIN1, HFZMAX, 20.1, 24.0),
00384     HcalFlexiHardcodeGeometryLoader::HFCellParameters (39, 2, 1, 2, 10, HFZMIN2, HFZMAX, 20.1, 24.0),
00385     HcalFlexiHardcodeGeometryLoader::HFCellParameters (40, 1, 3, 4, 20, HFZMIN1, HFZMAX, 16.9, 20.1),
00386     HcalFlexiHardcodeGeometryLoader::HFCellParameters (40, 2, 3, 4, 20, HFZMIN2, HFZMAX, 16.9, 20.1),
00387     HcalFlexiHardcodeGeometryLoader::HFCellParameters (41, 1, 3, 4, 20, HFZMIN1, HFZMAX, 12.5, 16.9),
00388     HcalFlexiHardcodeGeometryLoader::HFCellParameters (41, 2, 3, 4, 20, HFZMIN2, HFZMAX, 12.5, 16.9)
00389   };
00390   int nCells = sizeof(cells)/sizeof(HcalFlexiHardcodeGeometryLoader::HFCellParameters);
00391   std::vector <HcalFlexiHardcodeGeometryLoader::HFCellParameters> result;
00392   result.reserve (nCells);
00393   for (int i = 0; i < nCells; ++i) result.push_back (cells[i]);
00394   return result;
00395 }
00396   
00397 void HcalFlexiHardcodeGeometryLoader::fillHE (CaloSubdetectorGeometry* fGeometry, const std::vector <HcalFlexiHardcodeGeometryLoader::HECellParameters>& fCells) {
00398 
00399   for (size_t iCell = 0; iCell < fCells.size(); ++iCell) {
00400     const HcalFlexiHardcodeGeometryLoader::HECellParameters& param = fCells[iCell];
00401     for (int iPhi = param.phiFirst; iPhi <= MAX_HCAL_PHI; iPhi += param.phiStep) {
00402       for (int iside = -1; iside <= 1; iside += 2) { // both detector sides are identical
00403         HcalDetId hid (HcalEndcap, param.eta*iside, iPhi, param.depth);
00404         float phiCenter = ((iPhi-1)*360./MAX_HCAL_PHI + 0.5*param.dphi) * DEGREE2RAD; // middle of the cell
00405         float etaCenter = 0.5 * (param.etaMin + param.etaMax);
00406 
00407         float perp = param.zMin / sinh (etaCenter);
00408         float x = perp * cos (phiCenter);
00409         float y = perp * sin (phiCenter);
00410         float z = iside * param.zMin;
00411         // make cell geometry
00412         GlobalPoint refPoint (x,y,z); // center of the cell's face
00413         std::vector<CCGFloat> cellParams;
00414         cellParams.reserve (5);
00415         cellParams.push_back (0.5 * (param.etaMax - param.etaMin)); //deta_half
00416         cellParams.push_back (0.5 * param.dphi * DEGREE2RAD);  // dphi_half
00417         cellParams.push_back (-0.5 * (param.zMax - param.zMin) / tanh (etaCenter)); // dz_half, "-" means edges in Z
00418         cellParams.push_back ( fabs( refPoint.eta() ) ) ;
00419         cellParams.push_back ( fabs( refPoint.z() ) ) ;
00420 #ifdef DebugLog
00421         std::cout << "HcalFlexiHardcodeGeometryLoader::fillHE-> " << hid << refPoint << '/' << cellParams [0] << '/' << cellParams [1] << '/' << cellParams [2] << std::endl;
00422 #endif
00423         fGeometry->newCell(refPoint,  refPoint,  refPoint, 
00424                            CaloCellGeometry::getParmPtr(cellParams, 
00425                                                         fGeometry->parMgr(), 
00426                                                         fGeometry->parVecVec()),
00427                            hid ) ;
00428       }
00429     }
00430   }
00431 }
00432 
00433 void HcalFlexiHardcodeGeometryLoader::fillHF (CaloSubdetectorGeometry* fGeometry, const std::vector <HcalFlexiHardcodeGeometryLoader::HFCellParameters>& fCells) {
00434 
00435   for (size_t iCell = 0; iCell < fCells.size(); ++iCell) {
00436     const HcalFlexiHardcodeGeometryLoader::HFCellParameters& param = fCells[iCell];
00437     for (int iPhi = param.phiFirst; iPhi <= MAX_HCAL_PHI; iPhi += param.phiStep) {
00438       for (int iside = -1; iside <= 1; iside += 2) { // both detector sides are identical
00439         HcalDetId hid (HcalForward, param.eta*iside, iPhi, param.depth);
00440         float phiCenter = ((iPhi-1)*360./MAX_HCAL_PHI + 0.5*param.dphi) * DEGREE2RAD; // middle of the cell
00441         GlobalPoint inner (param.rMin, 0, param.zMin);
00442         GlobalPoint outer (param.rMax, 0, param.zMin);
00443         float iEta = inner.eta();
00444         float oEta = outer.eta();
00445         float etaCenter = 0.5 * ( iEta + oEta );
00446         
00447         float perp = param.zMin / sinh (etaCenter);
00448         float x = perp * cos (phiCenter);
00449         float y = perp * sin (phiCenter);
00450         float z = iside * param.zMin;
00451         // make cell geometry
00452         GlobalPoint refPoint (x,y,z); // center of the cell's face
00453         std::vector<CCGFloat> cellParams;
00454         cellParams.reserve (5);
00455         cellParams.push_back (0.5 * ( iEta - oEta )); // deta_half
00456         cellParams.push_back (0.5 * param.dphi * DEGREE2RAD);  // dphi_half
00457         cellParams.push_back (0.5 * (param.zMax - param.zMin)); // dz_half
00458         cellParams.push_back ( fabs( refPoint.eta()));
00459         cellParams.push_back ( fabs( refPoint.z() ) ) ;
00460 #ifdef DebugLog
00461         std::cout << "HcalFlexiHardcodeGeometryLoader::fillHF-> " << hid << refPoint << '/' << cellParams [0] << '/' << cellParams [1] << '/' << cellParams [2] << std::endl;
00462 #endif  
00463         fGeometry->newCell(refPoint,  refPoint,  refPoint, 
00464                            CaloCellGeometry::getParmPtr(cellParams, 
00465                                                         fGeometry->parMgr(), 
00466                                                         fGeometry->parVecVec()),
00467                            hid ) ;
00468       }
00469     }
00470   }
00471 }
/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_7/src/Geometry/HcalTowerAlgo/src/HcalFlexiHardcodeGeometryLoader.cc
