CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MuonRPCDetLayerGeometryBuilder.cc
Go to the documentation of this file.
2 
9 
13 
15 
16 #include <iostream>
17 
18 using namespace std;
19 
20 namespace rpcdetlayergeomsort {
21  template <class T, class Scalar = typename T::Scalar>
24  Scalar operator()(const T* p) const {return fabs(p->specificSurface().innerRadius());}
25  Scalar operator()(const T& p) const {return fabs(p.specificSurface().innerRadius());}
26  };
27 }
28 
29 
31 }
32 
33 
34 // Builds the forward (first) and backward (second) layers
35 pair<vector<DetLayer*>, vector<DetLayer*> >
37 
38  vector<DetLayer*> result[2];
39 
40  for (int endcap = -1; endcap<=1; endcap+=2) {
41  int iendcap = (endcap==1) ? 0 : 1; // +1: forward, -1: backward
42 
43  // ME 1
44  int firstStation=1;
45 
46  // ME 1/1
47  for (int layer = RPCDetId::minLayerId; layer <= RPCDetId::maxLayerId; ++layer) {
48  vector<int> rolls;
49  std::vector<int> rings;
50  int FirstStationRing = 1;
51  rings.push_back(FirstStationRing);
52  for(int roll = RPCDetId::minRollId+1;
53  roll <= RPCDetId::maxRollId; ++roll) {
54  rolls.push_back(roll);
55  }
56 
57 
58 
59  MuRingForwardDoubleLayer* ringLayer = buildLayer(endcap, rings,
60  firstStation , layer,
61  rolls, geo);
62  if (ringLayer) result[iendcap].push_back(ringLayer);
63 
64  }
65 
66  // ME 1/2 and ME1/3
67  for(int layer = RPCDetId::minLayerId; layer <= RPCDetId::maxLayerId; ++layer) {
68  vector<int> rolls;
69  std::vector<int> rings;
70  for(int ring = 2; ring <= 3; ++ring) {
71  rings.push_back(ring);
72  }
73  for(int roll = RPCDetId::minRollId+1; roll <= RPCDetId::maxRollId;
74  ++roll) {
75  rolls.push_back(roll);
76  }
77 
78  MuRingForwardDoubleLayer* ringLayer = buildLayer(endcap, rings, firstStation , layer, rolls, geo);
79  if (ringLayer) result[iendcap].push_back(ringLayer);
80  }
81 
82 
83  // ME 2 and ME 3
84  for(int station = 2; station <= RPCDetId::maxStationId; ++station) {
85  for(int layer = RPCDetId::minLayerId; layer <= RPCDetId::maxLayerId; ++layer) {
86  vector<int> rolls;
87  std::vector<int> rings;
89  rings.push_back(ring);
90  }
91  for(int roll = RPCDetId::minRollId+1; roll <= RPCDetId::maxRollId; ++roll) {
92  rolls.push_back(roll);
93  }
94 
95  MuRingForwardDoubleLayer* ringLayer = buildLayer(endcap, rings, station, layer, rolls, geo);
96  if (ringLayer) result[iendcap].push_back(ringLayer);
97  }
98  }
99 
100  }
101  pair<vector<DetLayer*>, vector<DetLayer*> > res_pair(result[0], result[1]);
102  return res_pair;
103 
104 }
105 
106 
107 
110  int layer,
111  vector<int>& rolls,
112  const RPCGeometry& geo) {
113 
114  const std::string metname = "Muon|RPC|RecoMuon|RecoMuonDetLayers|MuonRPCDetLayerGeometryBuilder";
115 
116  vector<const ForwardDetRing*> frontRings, backRings;
117 
118 
119  for (std::vector<int>::iterator ring=rings.begin(); ring<rings.end();++ring){
120  for (vector<int>::iterator roll = rolls.begin(); roll!=rolls.end(); ++roll) {
121  vector<const GeomDet*> frontDets, backDets;
122  for(int sector = RPCDetId::minSectorForwardId; sector <= RPCDetId::maxSectorForwardId; ++sector) {
123  for(int subsector = RPCDetId::minSubSectorForwardId; subsector <= RPCDetId::maxSectorForwardId; ++subsector) {
124  RPCDetId rpcId(endcap,*ring, station,sector,layer,subsector, (*roll));
125  bool isInFront = isFront(rpcId);
126  const GeomDet* geomDet = geo.idToDet(rpcId);
127  if (geomDet) {
128 
129  if(isInFront)
130  {
131  frontDets.push_back(geomDet);
132  }
133  else
134  {
135  backDets.push_back(geomDet);
136  }
137  LogTrace(metname) << "get RPC Endcap roll "
138  << rpcId
139  << (isInFront ? "front" : "back ")
140  << " at R=" << geomDet->position().perp()
141  << ", phi=" << geomDet->position().phi()
142  << ", Z=" << geomDet->position().z();
143 
144  }
145  }
146  }
147  if (frontDets.size()!=0) {
148  precomputed_value_sort(frontDets.begin(), frontDets.end(), geomsort::DetPhi());
149  frontRings.push_back(new MuDetRing(frontDets));
150  LogTrace(metname) << "New front ring with " << frontDets.size()
151  << " chambers at z="<< frontRings.back()->position().z();
152  }
153  if (backDets.size()!=0) {
154  precomputed_value_sort(backDets.begin(), backDets.end(), geomsort::DetPhi());
155  backRings.push_back(new MuDetRing(backDets));
156  LogTrace(metname) << "New back ring with " << backDets.size()
157  << " chambers at z="<< backRings.back()->position().z();
158  }
159  }
160  }
161 
163 
164  if(!backRings.empty() || !frontRings.empty())
165  {
167  precomputed_value_sort(frontRings.begin(), frontRings.end(), SortRingByInnerR());
168  precomputed_value_sort(backRings.begin(), backRings.end(), SortRingByInnerR());
169 
170  result = new MuRingForwardDoubleLayer(frontRings, backRings);
171  LogTrace(metname) << "New layer with " << frontRings.size()
172  << " front rings and " << backRings.size()
173  << " back rings, at Z " << result->position().z();
174  }
175 
176  return result;
177 }
178 
179 
180 vector<DetLayer*>
182  const std::string metname = "Muon|RPC|RecoMuon|RecoMuonDetLayers|MuonRPCDetLayerGeometryBuilder";
183 
184  vector<DetLayer*> detlayers;
185  vector<MuRodBarrelLayer*> result;
186  int region =0;
187 
189 
190  vector<const GeomDet*> geomDets;
191 
192  for(int layer=RPCDetId::minLayerId; layer<= RPCDetId::maxLayerId;++layer){
193  for(int sector = RPCDetId::minSectorId; sector <= RPCDetId::maxSectorId; sector++) {
194  for(int subsector = RPCDetId::minSubSectorId; subsector <= RPCDetId::maxSubSectorId; subsector++) {
195  for(int wheel = RPCDetId::minRingBarrelId; wheel <= RPCDetId::maxRingBarrelId; wheel++) {
196  for(int roll=RPCDetId::minRollId+1; roll <= RPCDetId::maxRollId; roll++){
197  const GeomDet* geomDet = geo.idToDet(RPCDetId(region,wheel,station,sector,layer,subsector,roll));
198  if (geomDet) {
199  geomDets.push_back(geomDet);
200  LogTrace(metname) << "get RPC Barrel roll " << RPCDetId(region,wheel,station,sector,layer,subsector,roll)
201  << " at R=" << geomDet->position().perp()
202  << ", phi=" << geomDet->position().phi() ;
203  }
204  }
205  }
206  }
207  }
208  }
209  makeBarrelLayers(geomDets, result);
210  }
211 
212 
213  for(vector<MuRodBarrelLayer*>::const_iterator it = result.begin(); it != result.end(); it++)
214  detlayers.push_back((DetLayer*)(*it));
215 
216  return detlayers;
217 }
218 
219 
220 void MuonRPCDetLayerGeometryBuilder::makeBarrelLayers(vector<const GeomDet *> & geomDets,
221  vector<MuRodBarrelLayer*> & result)
222 {
223  const std::string metname = "Muon|RPC|RecoMuon|RecoMuonDetLayers|MuonRPCDetLayerGeometryBuilder";
224 
225  //Sort in R
226  precomputed_value_sort(geomDets.begin(), geomDets.end(),geomsort::DetR());
227 
228  // Clusterize in phi - phi0
229  float resolution(25); // cm
230  float r0 = float(geomDets.front()->position().perp());
231  float rMin = - float(resolution);
232  float rMax = float(geomDets.back()->position().perp()) - r0 + resolution;
233 
234  ClusterizingHistogram hisR( int((rMax-rMin)/resolution) + 1,
235  rMin, rMax);
236 
237  vector<const GeomDet*>::iterator first = geomDets.begin();
238  vector<const GeomDet*>::iterator last = geomDets.end();
239 
240  for (vector<const GeomDet*>::iterator i=first; i!=last; i++){
241  hisR.fill(float((*i)->position().perp())-r0);
242  LogTrace(metname) << "R " << float((*i)->position().perp())-r0;
243  }
244  vector<float> rClust = hisR.clusterize(resolution);
245 
246  // LogTrace(metname) << " Found " << phiClust.size() << " clusters in Phi, ";
247 
248  vector<const GeomDet*>::iterator layerStart = first;
249  vector<const GeomDet*>::iterator separ = first;
250 
251  for (unsigned int i=0; i<rClust.size(); i++) {
252  float rSepar;
253  if (i<rClust.size()-1) {
254  rSepar = (rClust[i] + rClust[i+1])/2.f;
255  } else {
256  rSepar = rMax;
257  }
258 
259  // LogTrace(metname) << " cluster " << i
260  // << " phisepar " << phiSepar <<endl;
261  while (separ < last && float((*separ)->position().perp())-r0 < rSepar ) {
262  // LogTrace(metname) << " roll at dphi: " << float((*separ)->position().phi())-phi0;
263  separ++;
264  }
265 
266  if (int(separ-layerStart) > 0) {
267  // we have a layer in R. Now separate it into rods
268  vector<const DetRod*> rods;
269  vector<const GeomDet*> layerDets(layerStart, separ);
270  makeBarrelRods(layerDets, rods);
271 
272  if (rods.size()!=0) {
273  result.push_back(new MuRodBarrelLayer(rods));
274  LogTrace(metname) << " New MuRodBarrelLayer with " << rods.size()
275  << " rods, at R " << result.back()->specificSurface().radius();
276  }
277  }
278  layerStart = separ;
279  }
280 }
281 
282 
283 void
284 MuonRPCDetLayerGeometryBuilder::makeBarrelRods(vector<const GeomDet *> & geomDets,
285  vector<const DetRod*> & result)
286 {
287  const std::string metname = "Muon|RPC|RecoMuon|RecoMuonDetLayers|MuonRPCDetLayerGeometryBuilder";
288 
289  //Sort in phi
290  precomputed_value_sort(geomDets.begin(), geomDets.end(),geomsort::DetPhi());
291 
292  // Clusterize in phi - phi0
293  float resolution(0.01); // rad
294  float phi0 = float(geomDets.front()->position().phi());
295  float phiMin = - float(resolution);
296  float phiMax = float(geomDets.back()->position().phi()) - phi0 + resolution;
297 
298  ClusterizingHistogram hisPhi( int((phiMax-phiMin)/resolution) + 1,
299  phiMin, phiMax);
300 
301  vector<const GeomDet*>::iterator first = geomDets.begin();
302  vector<const GeomDet*>::iterator last = geomDets.end();
303 
304  for (vector<const GeomDet*>::iterator i=first; i!=last; i++){
305  hisPhi.fill(float((*i)->position().phi())-phi0);
306  LogTrace(metname) << "C " << float((*i)->position().phi())-phi0;
307  }
308  vector<float> phiClust = hisPhi.clusterize(resolution);
309 
310  // LogTrace(metname) << " Found " << phiClust.size() << " clusters in Phi, ";
311 
312  vector<const GeomDet*>::iterator rodStart = first;
313  vector<const GeomDet*>::iterator separ = first;
314 
315  for (unsigned int i=0; i<phiClust.size(); i++) {
316  float phiSepar;
317  if (i<phiClust.size()-1) {
318  phiSepar = (phiClust[i] + phiClust[i+1])/2.f;
319  } else {
320  phiSepar = phiMax;
321  }
322 
323  // LogTrace(metname) << " cluster " << i
324  // << " phisepar " << phiSepar <<endl;
325  while (separ < last && float((*separ)->position().phi())-phi0 < phiSepar ) {
326  // LogTrace(metname) << " roll at dphi: " << float((*separ)->position().phi())-phi0;
327  separ++;
328  }
329 
330  if (int(separ-rodStart) > 0) {
331  result.push_back(new MuDetRod(rodStart, separ));
332  LogTrace(metname) << " New MuDetRod with " << int(separ-rodStart)
333  << " rolls at R=" << (*rodStart)->position().perp()
334  << ", phi=" << float((*rodStart)->position().phi());
335 
336  }
337  rodStart = separ;
338  }
339 }
340 
341 
343 {
344  // ME1/2 is always in back
345  // if(rpcId.station() == 1 && rpcId.ring() == 2) return false;
346 
347  bool result = false;
348  int ring = rpcId.ring();
349  int station = rpcId.station();
350  // 20 degree rings are a little weird! not anymore from 17x
351  if(ring == 1 && station > 1)
352  {
353  // RE2/1 RE3/1 Upscope Geometry
354  /* goes (sector) (subsector) 1/3
355  1 1 back // front
356  1 2 front // back
357  1 3 front // front
358  2 1 front // back
359  2 2 back // from
360  2 3 back // back
361 
362  */
363  result = (rpcId.subsector() != 2);
364  if(rpcId.sector()%2 == 0) result = !result;
365  return result;
366  }
367  else
368  {
369  // 10 degree rings have odd subsectors in front
370  result = (rpcId.subsector()%2 == 0);
371  }
372  return result;
373 }
374 
375 
static const int maxStationId
Definition: RPCDetId.h:149
int i
Definition: DBlmapReader.cc:9
ExtractPhi< GeomDet, float > DetPhi
Definition: DetSorting.h:43
T perp() const
Definition: PV3DBase.h:71
static const int maxLayerId
Definition: RPCDetId.h:159
double Scalar
Definition: Definitions.h:27
const std::string metname
static const int maxRingForwardId
Definition: RPCDetId.h:143
static const int minRingBarrelId
Definition: RPCDetId.h:144
Geom::Phi< T > phi() const
Definition: PV3DBase.h:68
static const int minSubSectorId
Definition: RPCDetId.h:161
static MuRingForwardDoubleLayer * buildLayer(int endcap, std::vector< int > rings, int station, int layer, std::vector< int > &rolls, const RPCGeometry &geo)
static const int minRollId
Definition: RPCDetId.h:168
static const int minSectorId
Definition: RPCDetId.h:151
const Surface::PositionType & position() const
The position (origin of the R.F.)
Definition: GeomDet.h:41
static const int minSubSectorForwardId
Definition: RPCDetId.h:165
ExtractR< GeomDet, float > DetR
Definition: DetSorting.h:23
int ring() const
Definition: RPCDetId.h:76
static const int maxSubSectorId
Definition: RPCDetId.h:162
T z() const
Definition: PV3DBase.h:63
tuple result
Definition: query.py:137
static const int minSectorForwardId
Definition: RPCDetId.h:155
static bool isFront(const RPCDetId &rpcId)
double f[11][100]
static const int maxRollId
Definition: RPCDetId.h:169
bool first
Definition: L1TdeRCT.cc:94
#define LogTrace(id)
static const int minStationId
Definition: RPCDetId.h:148
static void makeBarrelRods(std::vector< const GeomDet * > &geomDets, std::vector< const DetRod * > &result)
static void makeBarrelLayers(std::vector< const GeomDet * > &geomDets, std::vector< MuRodBarrelLayer * > &result)
static const int maxRingBarrelId
Definition: RPCDetId.h:145
static const int maxSectorId
Definition: RPCDetId.h:152
int sector() const
Sector id: the group of chambers at same phi (and increasing r)
Definition: RPCDetId.h:106
static std::pair< std::vector< DetLayer * >, std::vector< DetLayer * > > buildEndcapLayers(const RPCGeometry &geo)
int subsector() const
SubSector id : some sectors are divided along the phi direction in subsectors (from 1 to 4 in Barrel...
Definition: RPCDetId.h:118
static const int minRingForwardId
Definition: RPCDetId.h:142
void precomputed_value_sort(RandomAccessIterator begin, RandomAccessIterator end, const Extractor &extr)
static const int maxSectorForwardId
Definition: RPCDetId.h:156
long double T
static const int minLayerId
Definition: RPCDetId.h:158
virtual const GeomDet * idToDet(DetId) const
Definition: RPCGeometry.cc:52
int station() const
Definition: RPCDetId.h:100
static std::vector< DetLayer * > buildBarrelLayers(const RPCGeometry &geo)
Builds the barrel layers. Result vector is sorted inside-out.