CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
RPCtoDTTranslator.cc
Go to the documentation of this file.
1 //-------------------------------------------------
2 //
3 // Class: RPCtoDTTranslator
4 //
5 // RPCtoDTTranslator
6 //
7 //
8 // Author :
9 // G. Flouris U Ioannina Mar. 2015
10 // modifications: G Karathanasis U Athens
11 //--------------------------------------------------
12 
13 #include <iostream>
14 #include <iomanip>
15 #include <iterator>
16 #include <cmath>
17 #include <map>
18 
25 
26 using namespace std;
27 
28 RPCtoDTTranslator::RPCtoDTTranslator(RPCDigiCollection const& inrpcDigis) : m_rpcDigis{inrpcDigis} {}
29 
30 namespace {
31  constexpr int max_rpc_bx = 2;
32  constexpr int min_rpc_bx = -2;
33 
34  struct rpc_hit {
35  int bx;
36  int station;
37  int sector;
38  int wheel;
39  RPCDetId detid;
40  int strip;
41  int roll;
42  int layer;
43  rpc_hit(int pbx, int pstation, int psector, int pwheel, RPCDetId pdet, int pstrip, int proll, int player)
44  : bx(pbx),
45  station(pstation),
46  sector(psector),
47  wheel(pwheel),
48  detid(pdet),
49  strip(pstrip),
50  roll(proll),
51  layer(player) {}
52  };
53 
54  //Need to shift the index so that index 0
55  // corresponds to min_rpc_bx
56  class BxToHit {
57  public:
58  BxToHit() : m_hits{} {} //zero initializes
59 
60  static bool outOfRange(int iBX) { return (iBX > max_rpc_bx or iBX < min_rpc_bx); }
61 
62  int& operator[](int iBX) { return m_hits[iBX - min_rpc_bx]; }
63 
64  size_t size() const { return m_hits.size(); }
65 
66  private:
67  std::array<int, max_rpc_bx - min_rpc_bx + 1> m_hits;
68  };
69 } // namespace
70 
71 void RPCtoDTTranslator::run(const RPCGeometry& rpcGeometry) {
72  std::vector<L1MuDTChambPhDigi> l1ttma_out;
73  std::vector<L1MuDTChambPhDigi> l1ttma_hits_out;
74 
75  std::vector<rpc_hit> vrpc_hit_layer1, vrpc_hit_layer2, vrpc_hit_st3, vrpc_hit_st4;
76 
78  for (auto chamber = m_rpcDigis.begin(); chamber != m_rpcDigis.end(); ++chamber) {
79  RPCDetId detid = (*chamber).first;
80  for (auto digi = (*chamber).second.first; digi != (*chamber).second.second; ++digi) {
81  if (detid.region() != 0)
82  continue; //Region = 0 Barrel
83  if (BxToHit::outOfRange(digi->bx()))
84  continue;
85  if (detid.layer() == 1)
86  vrpc_hit_layer1.emplace_back(digi->bx(),
87  detid.station(),
88  detid.sector(),
89  detid.ring(),
90  detid,
91  digi->strip(),
92  detid.roll(),
93  detid.layer());
94  if (detid.station() == 3)
95  vrpc_hit_st3.emplace_back(digi->bx(),
96  detid.station(),
97  detid.sector(),
98  detid.ring(),
99  detid,
100  digi->strip(),
101  detid.roll(),
102  detid.layer());
103  if (detid.layer() == 2)
104  vrpc_hit_layer2.emplace_back(digi->bx(),
105  detid.station(),
106  detid.sector(),
107  detid.ring(),
108  detid,
109  digi->strip(),
110  detid.roll(),
111  detid.layer());
112  if (detid.station() == 4)
113  vrpc_hit_st4.emplace_back(digi->bx(),
114  detid.station(),
115  detid.sector(),
116  detid.ring(),
117  detid,
118  digi->strip(),
119  detid.roll(),
120  detid.layer());
121  }
122  }
123 
124  vector<int> vcluster_size;
125  int cluster_id = -1;
126  int itr = 0;
127  // int hits[5][4][12][2][5][3][100]= {{{{{{{0}}}}}}};
128  std::map<RPCHitCleaner::detId_Ext, int> hits;
129  int cluster_size = 0;
130  for (auto chamber = m_rpcDigis.begin(); chamber != m_rpcDigis.end(); ++chamber) {
131  RPCDetId detid = (*chamber).first;
132  int strip_n1 = -10000;
133  int bx_n1 = -10000;
134  if (detid.region() != 0)
135  continue; //Region = 0 Barrel
136  for (auto digi = (*chamber).second.first; digi != (*chamber).second.second; ++digi) {
137  if (fabs(digi->bx()) > 3)
138  continue;
139  //Create cluster ids and store their size
140  //if((digi->strip()+1!=strip_n1)|| digi->bx()!=bx_n1){
141  if (abs(digi->strip() - strip_n1) != 1 || digi->bx() != bx_n1) {
142  if (itr != 0)
143  vcluster_size.push_back(cluster_size);
144  cluster_size = 0;
145  cluster_id++;
146  }
147  itr++;
148  cluster_size++;
150  //hits[(detid.ring()+2)][(detid.station()-1)][(detid.sector()-1)][(detid.layer()-1)][(digi->bx()+2)][detid.roll()-1][digi->strip()]= cluster_id ;
151  RPCHitCleaner::detId_Ext tmp{detid, digi->bx(), digi->strip()};
152  hits[tmp] = cluster_id;
154  strip_n1 = digi->strip();
155  bx_n1 = digi->bx();
156  }
157  }
158  vcluster_size.push_back(cluster_size);
159 
160  for (int wh = -2; wh <= 2; wh++) {
161  for (int sec = 1; sec <= 12; sec++) {
162  for (int st = 1; st <= 4; st++) {
163  int rpcbx = 0;
164  std::vector<int> delta_phib;
165  bool found_hits = false;
166  std::vector<int> rpc2dt_phi, rpc2dt_phib;
168  int itr1 = 0;
169  for (unsigned int l1 = 0; l1 < vrpc_hit_layer1.size(); l1++) {
170  RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer1[l1].detid, vrpc_hit_layer1[l1].bx, vrpc_hit_layer1[l1].strip};
171  int id = hits[tmp];
172  int phi1 = radialAngle(vrpc_hit_layer1[l1].detid, rpcGeometry, vrpc_hit_layer1[l1].strip);
173  if (vcluster_size[id] == 2 && itr1 == 0) {
174  itr1++;
175  continue;
176  }
177  if (vcluster_size[id] == 2 && itr1 == 1) {
178  itr1 = 0;
179  phi1 = phi1 + (radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip));
180  phi1 /= 2;
181  }
182  int itr2 = 0;
183  for (unsigned int l2 = 0; l2 < vrpc_hit_layer2.size(); l2++) {
184  if (vrpc_hit_layer1[l1].station != st || vrpc_hit_layer2[l2].station != st)
185  continue;
186  if (vrpc_hit_layer1[l1].sector != sec || vrpc_hit_layer2[l2].sector != sec)
187  continue;
188  if (vrpc_hit_layer1[l1].wheel != wh || vrpc_hit_layer2[l2].wheel != wh)
189  continue;
190  if (vrpc_hit_layer1[l1].bx != vrpc_hit_layer2[l2].bx)
191  continue;
192  RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer2[l2].detid, vrpc_hit_layer2[l2].bx, vrpc_hit_layer2[l2].strip};
193  int id = hits[tmp];
194 
195  if (vcluster_size[id] == 2 && itr2 == 0) {
196  itr2++;
197  continue;
198  }
199 
200  //int phi1 = radialAngle(vrpc_hit_layer1[l1].detid, c, vrpc_hit_layer1[l1].strip) ;
201  int phi2 = radialAngle(vrpc_hit_layer2[l2].detid, rpcGeometry, vrpc_hit_layer2[l2].strip);
202  if (vcluster_size[id] == 2 && itr2 == 1) {
203  itr2 = 0;
204  phi2 = phi2 + (radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip));
205  phi2 /= 2;
206  }
207  int average = ((phi1 + phi2) / 2) << 2; //10-bit->12-bit
208  rpc2dt_phi.push_back(average); //Convert and store to 12-bit
209  //int xin = localX(vrpc_hit_layer1[l1].detid, c, vrpc_hit_layer1[l1].strip);
210  //int xout = localX(vrpc_hit_layer2[l2].detid, c, vrpc_hit_layer2[l2].strip);
211  //cout<<(phi1<<2)<<" "<<l1<<" "<<vrpc_hit_layer1[l1].station<<endl;
212  //cout<<(phi2<<2)<<" "<<l1<<" "<<vrpc_hit_layer1[l1].station<<endl;
213  int xin = localXX((phi1 << 2), 1, vrpc_hit_layer1[l1].station);
214  int xout = localXX((phi2 << 2), 2, vrpc_hit_layer2[l2].station);
215  if (vcluster_size[id] == 2 && itr2 == 1) {
216  int phi1_n1 = radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip);
217  int phi2_n1 = radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip);
218  xin += localXX((phi1_n1 << 2), 1, vrpc_hit_layer1[l1].station);
219  xout += localXX((phi2_n1 << 2), 2, vrpc_hit_layer2[l2].station);
220  xin /= 2;
221  xout /= 2;
222  }
223  //cout<<">>"<<xin<<" "<<xout<<endl;
224  int phi_b = bendingAngle(xin, xout, average);
225  //cout<<"phib "<<phi_b<<endl;
226  rpc2dt_phib.push_back(phi_b);
227  //delta_phib to find the highest pt primitve
228  delta_phib.push_back(abs(phi_b));
229  found_hits = true;
230  rpcbx = vrpc_hit_layer1[l1].bx;
231  }
232  }
233  if (found_hits) {
234  //cout<<"found_hits"<<endl;
235  int min_index = std::distance(delta_phib.begin(), std::min_element(delta_phib.begin(), delta_phib.end())) + 0;
236  //cout<<min_index<<endl;
237  l1ttma_out.emplace_back(rpcbx, wh, sec - 1, st, rpc2dt_phi[min_index], rpc2dt_phib[min_index], 3, 0, 0, 2);
238  }
240  BxToHit hit;
241  itr1 = 0;
242  for (unsigned int l1 = 0; l1 < vrpc_hit_layer1.size(); l1++) {
243  if (vrpc_hit_layer1[l1].station != st || st > 2 || vrpc_hit_layer1[l1].sector != sec ||
244  vrpc_hit_layer1[l1].wheel != wh)
245  continue;
246  //int id = hits[vrpc_hit_layer1[l1].wheel+2][(vrpc_hit_layer1[l1].station-1)][(vrpc_hit_layer1[l1].sector-1)][(vrpc_hit_layer1[l1].layer-1)][(vrpc_hit_layer1[l1].bx+2)][vrpc_hit_layer1[l1].roll-1][vrpc_hit_layer1[l1].strip];
247  RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer1[l1].detid, vrpc_hit_layer1[l1].bx, vrpc_hit_layer1[l1].strip};
248  int id = hits[tmp];
249  if (vcluster_size[id] == 2 && itr1 == 0) {
250  itr1++;
251  continue;
252  }
253  int phi2 = radialAngle(vrpc_hit_layer1[l1].detid, rpcGeometry, vrpc_hit_layer1[l1].strip);
254  phi2 = phi2 << 2;
255  if (vcluster_size[id] == 2 && itr1 == 1) {
256  itr1 = 0;
257  phi2 = phi2 + (radialAngle(vrpc_hit_layer1[l1 - 1].detid, rpcGeometry, vrpc_hit_layer1[l1 - 1].strip) << 2);
258  phi2 /= 2;
259  }
260 
261  l1ttma_hits_out.emplace_back(
262  vrpc_hit_layer1[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_layer1[l1].bx], 0, 2);
263  hit[vrpc_hit_layer1[l1].bx]++;
264  }
265  itr1 = 0;
266  for (unsigned int l2 = 0; l2 < vrpc_hit_layer2.size(); l2++) {
267  if (vrpc_hit_layer2[l2].station != st || st > 2 || vrpc_hit_layer2[l2].sector != sec ||
268  vrpc_hit_layer2[l2].wheel != wh)
269  continue;
270  RPCHitCleaner::detId_Ext tmp{vrpc_hit_layer2[l2].detid, vrpc_hit_layer2[l2].bx, vrpc_hit_layer2[l2].strip};
271  int id = hits[tmp];
272  // int id = hits[vrpc_hit_layer2[l2].wheel+2][(vrpc_hit_layer2[l2].station-1)][(vrpc_hit_layer2[l2].sector-1)][(vrpc_hit_layer2[l2].layer-1)][(vrpc_hit_layer2[l2].bx+2)][vrpc_hit_layer2[l2].roll-1][vrpc_hit_layer2[l2].strip];
273  if (vcluster_size[id] == 2 && itr1 == 0) {
274  itr1++;
275  continue;
276  }
277  int phi2 = radialAngle(vrpc_hit_layer2[l2].detid, rpcGeometry, vrpc_hit_layer2[l2].strip);
278  phi2 = phi2 << 2;
279  if (vcluster_size[id] == 2 && itr1 == 1) {
280  itr1 = 0;
281  phi2 = phi2 + (radialAngle(vrpc_hit_layer2[l2 - 1].detid, rpcGeometry, vrpc_hit_layer2[l2 - 1].strip) << 2);
282  phi2 /= 2;
283  }
284  l1ttma_hits_out.emplace_back(
285  vrpc_hit_layer2[l2].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_layer2[l2].bx], 0, 2);
286  hit[vrpc_hit_layer2[l2].bx]++;
287  }
288  itr1 = 0;
289 
290  for (unsigned int l1 = 0; l1 < vrpc_hit_st3.size(); l1++) {
291  if (st != 3 || vrpc_hit_st3[l1].station != 3 || vrpc_hit_st3[l1].wheel != wh ||
292  vrpc_hit_st3[l1].sector != sec)
293  continue;
294  RPCHitCleaner::detId_Ext tmp{vrpc_hit_st3[l1].detid, vrpc_hit_st3[l1].bx, vrpc_hit_st3[l1].strip};
295  int id = hits[tmp];
296  //int id = hits[vrpc_hit_st3[l1].wheel+2][(vrpc_hit_st3[l1].station-1)][(vrpc_hit_st3[l1].sector-1)][(vrpc_hit_st3[l1].layer-1)][(vrpc_hit_st3[l1].bx+2)][vrpc_hit_st3[l1].roll-1][vrpc_hit_st3[l1].strip];
297  if (vcluster_size[id] == 2 && itr1 == 0) {
298  itr1++;
299  continue;
300  }
301  int phi2 = radialAngle(vrpc_hit_st3[l1].detid, rpcGeometry, vrpc_hit_st3[l1].strip);
302  phi2 = phi2 << 2;
303  if (vcluster_size[id] == 2 && itr1 == 1) {
304  itr1 = 0;
305  phi2 = phi2 + (radialAngle(vrpc_hit_st3[l1 - 1].detid, rpcGeometry, vrpc_hit_st3[l1 - 1].strip) << 2);
306  phi2 /= 2;
307  }
308  l1ttma_hits_out.emplace_back(
309  vrpc_hit_st3[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_st3[l1].bx], 0, 2);
310  hit[vrpc_hit_st3[l1].bx]++;
311  }
312  itr1 = 0;
313 
314  for (unsigned int l1 = 0; l1 < vrpc_hit_st4.size(); l1++) {
315  if (st != 4 || vrpc_hit_st4[l1].station != 4 || vrpc_hit_st4[l1].wheel != wh ||
316  vrpc_hit_st4[l1].sector != sec)
317  continue;
318  //int id = hits[vrpc_hit_st4[l1].wheel+2][(vrpc_hit_st4[l1].station-1)][(vrpc_hit_st4[l1].sector-1)][(vrpc_hit_st4[l1].layer-1)][(vrpc_hit_st4[l1].bx+2)][vrpc_hit_st4[l1].roll-1][vrpc_hit_st4[l1].strip];
319  RPCHitCleaner::detId_Ext tmp{vrpc_hit_st4[l1].detid, vrpc_hit_st4[l1].bx, vrpc_hit_st4[l1].strip};
320  int id = hits[tmp];
321  if (vcluster_size[id] == 2 && itr1 == 0) {
322  itr1++;
323  continue;
324  }
325  int phi2 = radialAngle(vrpc_hit_st4[l1].detid, rpcGeometry, vrpc_hit_st4[l1].strip);
326  phi2 = phi2 << 2;
327  if (vcluster_size[id] == 2 && itr1 == 1) {
328  itr1 = 0;
329  phi2 = phi2 + (radialAngle(vrpc_hit_st4[l1 - 1].detid, rpcGeometry, vrpc_hit_st4[l1 - 1].strip) << 2);
330  phi2 /= 2;
331  }
332  l1ttma_hits_out.emplace_back(
333  vrpc_hit_st4[l1].bx, wh, sec - 1, st, phi2, 0, 3, hit[vrpc_hit_st4[l1].bx], 0, 2);
334  hit[vrpc_hit_st4[l1].bx]++;
335  //l1ttma_out.push_back(rpc2dt_out);
336 
337  //break;
338  }
339  }
340  }
341  }
343  m_rpcdt_translated.setContainer(l1ttma_out);
345  m_rpchitsdt_translated.setContainer(l1ttma_hits_out);
346 }
347 
349 int RPCtoDTTranslator::radialAngle(RPCDetId detid, const RPCGeometry& rpcGeometry, int strip) {
350  int radialAngle;
351  // from phiGlobal to radialAngle of the primitive in sector sec in [1..12] <- RPC scheme
352 
353  const RPCRoll* roll = rpcGeometry.roll(detid);
354  GlobalPoint stripPosition = roll->toGlobal(roll->centreOfStrip(strip));
355 
356  double globalphi = stripPosition.phi();
357  int sector = (roll->id()).sector();
358  if (sector == 1)
359  radialAngle = int(globalphi * 1024);
360  else {
361  if (globalphi >= 0)
362  radialAngle = int((globalphi - (sector - 1) * Geom::pi() / 6.) * 1024);
363  else
364  radialAngle = int((globalphi + (13 - sector) * Geom::pi() / 6.) * 1024);
365  }
366  return radialAngle;
367 }
368 
370 int RPCtoDTTranslator::localX(RPCDetId detid, const RPCGeometry& rpcGeometry, int strip) {
371  const RPCRoll* roll = rpcGeometry.roll(detid);
372 
374  GlobalPoint p1cmPRG = roll->toGlobal(LocalPoint(1, 0, 0));
375  GlobalPoint m1cmPRG = roll->toGlobal(LocalPoint(-1, 0, 0));
376  float phip1cm = p1cmPRG.phi();
377  float phim1cm = m1cmPRG.phi();
378  int direction = (phip1cm - phim1cm) / abs(phip1cm - phim1cm);
380 
381  return direction * roll->centreOfStrip(strip).x();
382 }
383 
384 int RPCtoDTTranslator::bendingAngle(int xin, int xout, int phi) {
385  // use chamber size and max angle in hw units 512
386  int atanv = (int)(atan((xout - xin) / 34.6) * 512);
387  if (atanv > 512)
388  return 512;
389  int rvalue = atanv - phi / 8;
390  return rvalue;
391 }
392 
394  //R[stat][layer] - radius of rpc station/layer from center of CMS
395  double R[2][2] = {{410.0, 444.8}, {492.7, 527.3}};
396  double rvalue = R[station - 1][layer - 1] * tan(phi / 4096.);
397  return rvalue;
398 }
LocalPoint centreOfStrip(int strip) const
Definition: RPCRoll.cc:26
tuple array
Definition: mps_check.py:216
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventIDconst &, edm::Timestampconst & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
static int bendingAngle(int, int, int)
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
RPCtoDTTranslator(const RPCDigiCollection &inrpcDigis)
static int localX(RPCDetId, const RPCGeometry &, int)
function - will be replaced by LUTs(?)
const RPCDigiCollection & m_rpcDigis
constexpr std::array< uint8_t, layerIndexSize > layer
RPCDetId id() const
Definition: RPCRoll.cc:16
int roll() const
Definition: RPCDetId.h:92
void setContainer(Phi_Container inputSegments)
L1MuDTChambPhContainer m_rpcdt_translated
Output PhContainer.
int ring() const
Definition: RPCDetId.h:59
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
T operator[](int i) const
void run(const RPCGeometry &)
int layer() const
Definition: RPCDetId.h:85
L1MuDTChambPhContainer m_rpchitsdt_translated
int sector() const
Sector id: the group of chambers at same phi (and increasing r)
Definition: RPCDetId.h:81
int average
Definition: PDRates.py:138
static int localXX(int, int, int)
constexpr double pi()
Definition: Pi.h:31
tmp
align.sh
Definition: createJobs.py:716
T x() const
Definition: PV3DBase.h:59
tuple size
Write out results.
static int radialAngle(RPCDetId, const RPCGeometry &, int)
function - will be replaced by LUTs(?)
const RPCRoll * roll(RPCDetId id) const
Return a roll given its id.
Definition: RPCGeometry.cc:50
int region() const
Region id: 0 for Barrel, +/-1 For +/- Endcap.
Definition: RPCDetId.h:53
int station() const
Definition: RPCDetId.h:78