DetLayers/src/MuonCSCDetLayerGeometryBuilder.cc

0001 #include <RecoMuon/DetLayers/src/MuonCSCDetLayerGeometryBuilder.h>
0002
0003 #include <RecoMuon/DetLayers/interface/MuRingForwardDoubleLayer.h>
0004 #include <RecoMuon/DetLayers/interface/MuDetRing.h>
0005 #include <DataFormats/MuonDetId/interface/CSCDetId.h>
0006 #include <Geometry/CommonDetUnit/interface/GeomDet.h>
0007
0008 #include <Utilities/General/interface/precomputed_value_sort.h>
0009 #include <Geometry/CommonDetUnit/interface/DetSorting.h>
0010
0011 #include <FWCore/MessageLogger/interface/MessageLogger.h>
0012
0013 #include <iostream>
0014
0015 using namespace std;
0016
0017 pair<vector<DetLayer*>, vector<DetLayer*> > MuonCSCDetLayerGeometryBuilder::buildLayers(const CSCGeometry& geo) {
0018   vector<DetLayer*> result[2];  // one for each endcap
0019
0020   for (int i = 0; i < 2; i++) {
0021     int endcap = i + 1;
0022
0023     // ME/1/1a (= station 1, ring 4) and ME/1/1b (= station 1, ring 1)
0024     {
0025       vector<int> rings;
0026       rings.push_back(4);
0027       rings.push_back(1);
0028
0029       MuRingForwardDoubleLayer* layer = buildLayer(endcap, 1, rings, geo);
0030       if (layer)
0031         result[i].push_back(layer);
0032     }
0033
0034     // ME/1/2 and 1/3 (= station 1, ring 2 and 3)
0035     {
0036       vector<int> rings;
0037       rings.push_back(2);
0038       rings.push_back(3);
0039
0040       MuRingForwardDoubleLayer* layer = buildLayer(endcap, 1, rings, geo);
0041       if (layer)
0042         result[i].push_back(layer);
0043     }
0044
0045     // Stations 2,3,4
0046     for (int station = 2; station <= CSCDetId::maxStationId(); station++) {
0047       vector<int> rings;
0048       for (int ring = CSCDetId::minRingId(); ring <= CSCDetId::maxRingId(); ring++) {
0049         rings.push_back(ring);
0050       }
0051       MuRingForwardDoubleLayer* layer = buildLayer(endcap, station, rings, geo);
0052       if (layer)
0053         result[i].push_back(layer);
0054     }
0055   }
0056   pair<vector<DetLayer*>, vector<DetLayer*> > res_pair(result[0], result[1]);
0057   return res_pair;
0058 }
0059
0060 MuRingForwardDoubleLayer* MuonCSCDetLayerGeometryBuilder::buildLayer(int endcap,
0061                                                                      int station,
0062                                                                      vector<int>& rings,
0063                                                                      const CSCGeometry& geo) {
0064   const std::string metname = "Muon|RecoMuon|RecoMuonDetLayers|MuonCSCDetLayerGeometryBuilder";
0065   MuRingForwardDoubleLayer* result = nullptr;
0066
0067   vector<const ForwardDetRing*> frontRings, backRings;
0068
0069   for (vector<int>::iterator ring = rings.begin(); ring != rings.end(); ring++) {
0070     vector<const GeomDet*> frontGeomDets, backGeomDets;
0071     for (int chamber = CSCDetId::minChamberId(); chamber <= CSCDetId::maxChamberId(); chamber++) {
0072       CSCDetId detId(endcap, station, (*ring), chamber, 0);
0073       const GeomDet* geomDet = geo.idToDet(detId);
0074       // we sometimes loop over more chambers than there are in ring
0075       bool isInFront = isFront(station, *ring, chamber);
0076       if (geomDet != nullptr) {
0077         if (isInFront) {
0078           frontGeomDets.push_back(geomDet);
0079         } else {
0080           backGeomDets.push_back(geomDet);
0081         }
0082         LogTrace(metname) << "get CSC chamber " << CSCDetId(endcap, station, (*ring), chamber, 0)
0083                           << " at R=" << geomDet->position().perp() << ", phi=" << geomDet->position().phi()
0084                           << ", z= " << geomDet->position().z() << " isFront? " << isInFront;
0085       }
0086     }
0087
0088     if (!backGeomDets.empty()) {
0089       backRings.push_back(makeDetRing(backGeomDets));
0090     }
0091
0092     if (!frontGeomDets.empty()) {
0093       frontRings.push_back(makeDetRing(frontGeomDets));
0094       assert(!backGeomDets.empty());
0095       float frontz = frontRings[0]->position().z();
0096       float backz = backRings[0]->position().z();
0097       assert(fabs(frontz) < fabs(backz));
0098     }
0099   }
0100
0101   // How should they be sorted?
0102   //    precomputed_value_sort(muDetRods.begin(), muDetRods.end(), geomsort::ExtractZ<GeometricSearchDet,float>());
0103   result = new MuRingForwardDoubleLayer(frontRings, backRings);
0104   LogTrace(metname) << "New MuRingForwardLayer with " << frontRings.size() << " and " << backRings.size()
0105                     << " rings, at Z " << result->position().z() << " R1: " << result->specificSurface().innerRadius()
0106                     << " R2: " << result->specificSurface().outerRadius();
0107   return result;
0108 }
0109
0110 bool MuonCSCDetLayerGeometryBuilder::isFront(int station, int ring, int chamber) {
0111   bool result = false;
0112
0113   bool isOverlapping = !(station == 1 && ring == 3);
0114   // not overlapping means back
0115   if (isOverlapping) {
0116     bool isEven = (chamber % 2 == 0);
0117     // odd chambers are bolted to the iron, which faces
0118     // forward in 1&2, backward in 3&4, so...
0119     result = (station < 3) ? isEven : !isEven;
0120   }
0121   return result;
0122 }
0123
0124 MuDetRing* MuonCSCDetLayerGeometryBuilder::makeDetRing(vector<const GeomDet*>& geomDets) {
0125   const std::string metname = "Muon|RecoMuon|RecoMuonDetLayers|MuonCSCDetLayerGeometryBuilder";
0126
0127   precomputed_value_sort(geomDets.begin(), geomDets.end(), geomsort::DetPhi());
0128   MuDetRing* result = new MuDetRing(geomDets);
0129   LogTrace(metname) << "New MuDetRing with " << geomDets.size() << " chambers at z=" << result->position().z()
0130                     << " R1: " << result->specificSurface().innerRadius()
0131                     << " R2: " << result->specificSurface().outerRadius();
0132   return result;
0133 }