CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes
DTTFFEDSim Class Reference

#include <DTTFFEDSim.h>

Inheritance diagram for DTTFFEDSim:
edm::stream::EDProducer<>

Public Member Functions

 DTTFFEDSim (const edm::ParameterSet &pset)
 Constructor. More...
 
bool fillRawData (edm::Event &e, FEDRawDataCollection &data)
 Generate and fill FED raw data for a full event. More...
 
void produce (edm::Event &e, const edm::EventSetup &c) override
 Produce digis out of raw data. More...
 
 ~DTTFFEDSim () override
 Destructor. More...
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
bool hasAbilityToProduceInLumis () const final
 
bool hasAbilityToProduceInRuns () const final
 

Private Member Functions

int bxNr (int channel)
 
int channel (int wheel, int sector, int bx)
 
edm::InputTag getDTDigiInputTag ()
 
edm::InputTag getDTPHTFInputTag ()
 
int sector (int channel)
 
int wheel (int channel)
 

Private Attributes

edm::EDGetTokenT< L1MuDTChambPhContainerChPh_tok
 
edm::EDGetTokenT< L1MuDTChambThContainerChTh_tok
 
edm::InputTag DTDigiInputTag
 
edm::InputTag DTPHTFInputTag
 
unsigned int eventNum
 
edm::EDGetTokenT< L1MuDTTrackContainerTrk_tok
 

Additional Inherited Members

- Public Types inherited from edm::stream::EDProducer<>
typedef CacheContexts< T... > CacheTypes
 
typedef CacheTypes::GlobalCache GlobalCache
 
typedef AbilityChecker< T... > HasAbility
 
typedef CacheTypes::LuminosityBlockCache LuminosityBlockCache
 
typedef LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCacheLuminosityBlockContext
 
typedef CacheTypes::LuminosityBlockSummaryCache LuminosityBlockSummaryCache
 
typedef CacheTypes::RunCache RunCache
 
typedef RunContextT< RunCache, GlobalCacheRunContext
 
typedef CacheTypes::RunSummaryCache RunSummaryCache
 

Detailed Description

L1 DT Track Finder Digi-to-Raw

J. Troconiz UAM Madrid

Definition at line 29 of file DTTFFEDSim.h.

Constructor & Destructor Documentation

DTTFFEDSim::DTTFFEDSim ( const edm::ParameterSet pset)

Constructor.

Definition at line 26 of file DTTFFEDSim.cc.

References ChPh_tok, ChTh_tok, DTDigiInputTag, DTPHTFInputTag, edm::ParameterSet::getParameter(), and Trk_tok.

26  : eventNum(0) {
27 
28  produces<FEDRawDataCollection>();
29 
30  DTDigiInputTag = pset.getParameter<edm::InputTag>("DTDigi_Source");
31  DTPHTFInputTag = pset.getParameter<edm::InputTag>("DTTracks_Source");
32 
33  ChPh_tok = consumes<L1MuDTChambPhContainer>(DTDigiInputTag);
34  ChTh_tok = consumes<L1MuDTChambThContainer>(DTDigiInputTag);
35  Trk_tok = consumes<L1MuDTTrackContainer>(DTPHTFInputTag);
36 }
T getParameter(std::string const &) const
edm::EDGetTokenT< L1MuDTChambThContainer > ChTh_tok
Definition: DTTFFEDSim.h:63
unsigned int eventNum
Definition: DTTFFEDSim.h:45
edm::InputTag DTDigiInputTag
Definition: DTTFFEDSim.h:47
edm::EDGetTokenT< L1MuDTTrackContainer > Trk_tok
Definition: DTTFFEDSim.h:64
edm::EDGetTokenT< L1MuDTChambPhContainer > ChPh_tok
Definition: DTTFFEDSim.h:62
edm::InputTag DTPHTFInputTag
Definition: DTTFFEDSim.h:48
DTTFFEDSim::~DTTFFEDSim ( )
override

Destructor.

Definition at line 38 of file DTTFFEDSim.cc.

38 {}

Member Function Documentation

int DTTFFEDSim::bxNr ( int  channel)
private

Definition at line 291 of file DTTFFEDSim.cc.

References channel().

291  {
292 
293  int myChannel = channel;
294 
295  if (myChannel > 127)
296  myChannel -= 2;
297 
298  if (myChannel < 0 || myChannel > 251) {
299  return -999;
300  }
301 
302  int myBx = 1 - (myChannel % 3);
303 
304  return myBx;
305 }
int channel(int wheel, int sector, int bx)
Definition: DTTFFEDSim.cc:264
int DTTFFEDSim::channel ( int  wheel,
int  sector,
int  bx 
)
private

Definition at line 264 of file DTTFFEDSim.cc.

References funct::abs().

Referenced by bxNr(), fillRawData(), sector(), and wheel().

264  {
265 
266  // wheel : -3 -2 -1 +1 +2 +3 <=> PHTF's : N2, N1, N0, P0, P1, P2
267  // 0 <=> ETTF
268  // sector : 0 -> 11
269  // bx : -1 -> +1
270 
271  int myChannel = 255;
272 
273  if (abs(bx) > 1) {
274  return myChannel;
275  }
276  if (sector < 0 || sector > 11) {
277  return myChannel;
278  }
279  if (abs(wheel) > 3) {
280  return myChannel;
281  }
282 
283  myChannel = sector * 21 + wheel * 3 - bx + 10;
284 
285  if (myChannel > 125)
286  myChannel += 2;
287 
288  return myChannel;
289 }
int sector(int channel)
Definition: DTTFFEDSim.cc:307
int wheel(int channel)
Definition: DTTFFEDSim.cc:321
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool DTTFFEDSim::fillRawData ( edm::Event e,
FEDRawDataCollection data 
)

Generate and fill FED raw data for a full event.

Definition at line 52 of file DTTFFEDSim.cc.

References L1MuDTChambPhContainer::bxSize(), L1MuDTChambThContainer::bxSize(), L1MuDTTrackContainer::bxSize(), dt_crc::calcCRC(), channel(), ChPh_tok, ChTh_tok, FEDRawData::data(), edm::EventID::event(), eventNum, FEDRawDataCollection::FEDData(), edm::Event::getByToken(), L1MuDTChambThContainer::getContainer(), L1MuDTChambPhContainer::getContainer(), L1MuDTTrackContainer::getContainer(), edm::EventBase::id(), createfilelist::int, groupFilesInBlocks::lines, FEDRawData::resize(), and Trk_tok.

Referenced by produce().

52  {
53 
54  eventNum = e.id().event();
55 
56  int lines = 2;
57 
59  e.getByToken(ChPh_tok, phtrig);
60  lines += phtrig->bxSize(-1, 1);
61 
63  e.getByToken(ChTh_tok, thtrig);
64  lines += thtrig->bxSize(-1, 1);
65 
67  e.getByToken(Trk_tok, trtrig);
68  lines += trtrig->bxSize(-1, 1) * 3;
69 
70  FEDRawData &dttfdata = data.FEDData(0x30C);
71  dttfdata.resize(lines * 8); // size in bytes
72  unsigned char *LineFED = dttfdata.data();
73 
74  int *dataWord1 = new int;
75  int *dataWord2 = new int;
76 
77  //--> Header
78 
79  *dataWord1 = 0x50000000 + (eventNum & 0xFFFFFF);
80  *dataWord2 = 0x00030C00;
81 
82  int newCRC = 0xFFFF;
83  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
84 
85  *((int *)LineFED) = *dataWord2;
86  LineFED += 4;
87  *((int *)LineFED) = *dataWord1;
88 
89  //--> DTTF data
90 
91  int TS1Id[4], TS2Id[4]; // word identifier for TS #1,#2 for stations
92  TS1Id[0] = 0x0E;
93  TS2Id[0] = 0x1E;
94  TS1Id[1] = 0x2E;
95  TS2Id[1] = 0x3E;
96  TS1Id[3] = 0x4E;
97  TS2Id[3] = 0x5E;
98  TS1Id[2] = 0x8FFF8;
99  TS2Id[2] = 0x9FFF8;
100 
101  // Input
103 
104  for (tsphi = phtrig->getContainer()->begin();
105  tsphi != phtrig->getContainer()->end(); tsphi++) {
106  if (tsphi->code() != 7) {
107 
108  int wheelID = tsphi->whNum() + 1;
109  if (wheelID <= 0)
110  wheelID -= 2;
111  int stationID = tsphi->stNum() - 1;
112  int is2nd = tsphi->Ts2Tag();
113 
114  int channelNr = channel(wheelID, tsphi->scNum(), tsphi->bxNum() - is2nd);
115  if (channelNr == 255)
116  continue;
117  int TSId = (is2nd == 0) ? TS1Id[stationID] : TS2Id[stationID];
118 
119  *dataWord1 = ((channelNr & 0xFF) << 24) + 0x00FFFFFF;
120 
121  if (stationID != 2) {
122  *dataWord2 = ((TSId & 0x0FF) << 24) + (~(tsphi->code() + 1) & 0x007) +
123  ((~tsphi->phiB() & 0x3FF) << 3) +
124  ((~tsphi->phi() & 0xFFF) << 13);
125  } else {
126  *dataWord2 = ((TSId & 0xFFFFF) << 12) +
127  (~(tsphi->code() + 1) & 0x00007) +
128  ((~tsphi->phi() & 0x00FFF) << 3);
129  }
130 
131  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
132 
133  LineFED += 4;
134  *((int *)LineFED) = *dataWord2;
135  LineFED += 4;
136  *((int *)LineFED) = *dataWord1;
137  }
138  }
139  // Input
140 
141  // Input
143 
144  for (tsthe = thtrig->getContainer()->begin();
145  tsthe != thtrig->getContainer()->end(); tsthe++) {
146 
147  int wheelTh = tsthe->whNum();
148  int sectorID = tsthe->scNum();
149 
150  int channelNr = channel(0, sectorID, tsthe->bxNum());
151  if (channelNr == 255)
152  continue;
153  int TSId = wheelTh + 2;
154 
155  *dataWord1 = ((channelNr & 0xFF) << 24) + 0x00FFFFFF;
156 
157  *dataWord2 = ((TSId & 0x07) << 28) + 0x0FFFFFFF;
158 
159  int stationID = tsthe->stNum() - 1;
160  for (int bti = 0; bti < 7; bti++)
161  if (wheelTh == -2 || wheelTh == -1 ||
162  (wheelTh == 0 && (sectorID == 0 || sectorID == 3 || sectorID == 4 ||
163  sectorID == 7 || sectorID == 8 || sectorID == 11)))
164  *dataWord2 -= (tsthe->position(bti) & 0x1) << (stationID * 7 + bti);
165  else
166  *dataWord2 -= (tsthe->position(6 - bti) & 0x1) << (stationID * 7 + bti);
167 
168  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
169 
170  LineFED += 4;
171  *((int *)LineFED) = *dataWord2;
172  LineFED += 4;
173  *((int *)LineFED) = *dataWord1;
174  }
175  // Input
176 
177  // Output
179 
180  for (tstrk = trtrig->getContainer()->begin();
181  tstrk != trtrig->getContainer()->end(); tstrk++) {
182 
183  int channelNr = channel(tstrk->whNum(), tstrk->scNum(), tstrk->bx());
184  if (channelNr == 255)
185  continue;
186  int TSId = (tstrk->TrkTag() == 0) ? 0xAFFF : 0xBFFF;
187 
188  *dataWord1 = ((channelNr & 0xFF) << 24) + 0x00FFFFFF;
189 
190  *dataWord2 = ((TSId & 0xFFFF) << 16) + (tstrk->stNum(4) & 0x0000F) +
191  ((tstrk->stNum(3) & 0x0000F) << 4) +
192  ((tstrk->stNum(2) & 0x0000F) << 8) +
193  ((tstrk->stNum(1) & 0x00003) << 12);
194 
195  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
196 
197  LineFED += 4;
198  *((int *)LineFED) = *dataWord2;
199  LineFED += 4;
200  *((int *)LineFED) = *dataWord1;
201 
202  TSId = (tstrk->TrkTag() == 0) ? 0xCFFE : 0xDFFE;
203 
204  *dataWord1 = ((channelNr & 0xFF) << 24) + 0x00FFFFFF;
205 
206  *dataWord2 = ((TSId & 0xFFFE) << 16) + (~tstrk->quality_packed() & 0x0007) +
207  ((tstrk->phi_packed() & 0x00FF) << 3) +
208  ((~tstrk->charge_packed() & 0x0001) << 11) +
209  ((~tstrk->pt_packed() & 0x001F) << 12);
210 
211  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
212 
213  LineFED += 4;
214  *((int *)LineFED) = *dataWord2;
215  LineFED += 4;
216  *((int *)LineFED) = *dataWord1;
217 
218  channelNr = channel(0, tstrk->scNum(), tstrk->bx());
219  if (channelNr == 255)
220  continue;
221  TSId = (tstrk->whNum() + 3) << 16;
222  TSId += (tstrk->whNum() < 0) ? 0x8FFFC : 0x7FFFC;
223 
224  *dataWord1 = ((channelNr & 0xFF) << 24) + 0x00FFFFFF;
225 
226  *dataWord2 = (TSId & 0xFFFFC) << 12;
227 
228  if (tstrk->TrkTag() == 0) {
229  *dataWord2 += 0x3F80 + (tstrk->eta_packed() & 0x003F) +
230  ((~tstrk->finehalo_packed() & 0x0001) << 6);
231  } else {
232  *dataWord2 += 0x007F + ((tstrk->eta_packed() & 0x003F) << 7) +
233  ((~tstrk->finehalo_packed() & 0x0001) << 13);
234  }
235 
236  dt_crc::calcCRC(*dataWord1, *dataWord2, newCRC);
237 
238  LineFED += 4;
239  *((int *)LineFED) = *dataWord2;
240  LineFED += 4;
241  *((int *)LineFED) = *dataWord1;
242  }
243  // Output
244 
245  //--> Trailer
246 
247  *dataWord1 = 0xA0000000 + (lines & 0xFFFFFF);
248  *dataWord2 = 0;
249 
250  dt_crc::calcCRC(*dataWord1, *dataWord2 & 0xFFFF, newCRC);
251 
252  *dataWord2 += (newCRC & 0xFFFF) << 16;
253 
254  LineFED += 4;
255  *((int *)LineFED) = *dataWord2;
256  LineFED += 4;
257  *((int *)LineFED) = *dataWord1;
258 
259  delete dataWord1;
260  delete dataWord2;
261  return true;
262 }
EventNumber_t event() const
Definition: EventID.h:41
The_Container const * getContainer() const
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
int bxSize(int step1, int step2) const
edm::EDGetTokenT< L1MuDTChambThContainer > ChTh_tok
Definition: DTTFFEDSim.h:63
The_Container::const_iterator The_iterator
int bxSize(int step1, int step2) const
unsigned int eventNum
Definition: DTTFFEDSim.h:45
const FEDRawData & FEDData(int fedid) const
retrieve data for fed
void resize(size_t newsize)
Definition: FEDRawData.cc:32
Phi_Container::const_iterator Phi_iterator
int channel(int wheel, int sector, int bx)
Definition: DTTFFEDSim.cc:264
TrackContainer::const_iterator Trackiterator
Phi_Container const * getContainer() const
int bxSize(int step1, int step2) const
edm::EventID id() const
Definition: EventBase.h:59
void calcCRC(long, int &)
Definition: DTCRC.cc:3
edm::EDGetTokenT< L1MuDTTrackContainer > Trk_tok
Definition: DTTFFEDSim.h:64
edm::EDGetTokenT< L1MuDTChambPhContainer > ChPh_tok
Definition: DTTFFEDSim.h:62
const unsigned char * data() const
Return a const pointer to the beginning of the data buffer.
Definition: FEDRawData.cc:28
TrackContainer const * getContainer() const
edm::InputTag DTTFFEDSim::getDTDigiInputTag ( )
inlineprivate

Definition at line 59 of file DTTFFEDSim.h.

References DTDigiInputTag.

59 { return DTDigiInputTag; }
edm::InputTag DTDigiInputTag
Definition: DTTFFEDSim.h:47
edm::InputTag DTTFFEDSim::getDTPHTFInputTag ( )
inlineprivate

Definition at line 60 of file DTTFFEDSim.h.

References DTPHTFInputTag.

60 { return DTPHTFInputTag; }
edm::InputTag DTPHTFInputTag
Definition: DTTFFEDSim.h:48
void DTTFFEDSim::produce ( edm::Event e,
const edm::EventSetup c 
)
override

Produce digis out of raw data.

Definition at line 40 of file DTTFFEDSim.cc.

References data, fillRawData(), eostools::move(), and edm::Event::put().

Referenced by JSONExport.JsonExport::export(), HTMLExport.HTMLExport::export(), and HTMLExport.HTMLExportStatic::export().

40  {
41 
43 
44  if (!fillRawData(e, data))
45  return;
46 
47  unique_ptr<FEDRawDataCollection> fed_product(new FEDRawDataCollection(data));
48 
49  e.put(std::move(fed_product));
50 }
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:125
bool fillRawData(edm::Event &e, FEDRawDataCollection &data)
Generate and fill FED raw data for a full event.
Definition: DTTFFEDSim.cc:52
char data[epos_bytes_allocation]
Definition: EPOS_Wrapper.h:82
def move(src, dest)
Definition: eostools.py:511
int DTTFFEDSim::sector ( int  channel)
private

Definition at line 307 of file DTTFFEDSim.cc.

References channel().

Referenced by geometryXMLparser.DTAlignable::index().

307  {
308 
309  int myChannel = channel;
310 
311  if (myChannel > 127)
312  myChannel -= 2;
313 
314  if (myChannel < 0 || myChannel > 251) {
315  return -999;
316  }
317 
318  return myChannel / 21;
319 }
int channel(int wheel, int sector, int bx)
Definition: DTTFFEDSim.cc:264
int DTTFFEDSim::wheel ( int  channel)
private

Definition at line 321 of file DTTFFEDSim.cc.

References channel().

Referenced by geometryXMLparser.DTAlignable::index().

321  {
322 
323  int myChannel = channel;
324 
325  if (myChannel > 127)
326  myChannel -= 2;
327 
328  if (myChannel < 0 || myChannel > 251) {
329  return -999;
330  }
331 
332  int myWheel = ((myChannel % 21) / 3) - 3;
333 
334  return myWheel;
335 }
int channel(int wheel, int sector, int bx)
Definition: DTTFFEDSim.cc:264

Member Data Documentation

edm::EDGetTokenT<L1MuDTChambPhContainer> DTTFFEDSim::ChPh_tok
private

Definition at line 62 of file DTTFFEDSim.h.

Referenced by DTTFFEDSim(), and fillRawData().

edm::EDGetTokenT<L1MuDTChambThContainer> DTTFFEDSim::ChTh_tok
private

Definition at line 63 of file DTTFFEDSim.h.

Referenced by DTTFFEDSim(), and fillRawData().

edm::InputTag DTTFFEDSim::DTDigiInputTag
private

Definition at line 47 of file DTTFFEDSim.h.

Referenced by DTTFFEDSim(), and getDTDigiInputTag().

edm::InputTag DTTFFEDSim::DTPHTFInputTag
private

Definition at line 48 of file DTTFFEDSim.h.

Referenced by DTTFFEDSim(), and getDTPHTFInputTag().

unsigned int DTTFFEDSim::eventNum
private

Definition at line 45 of file DTTFFEDSim.h.

Referenced by fillRawData().

edm::EDGetTokenT<L1MuDTTrackContainer> DTTFFEDSim::Trk_tok
private

Definition at line 64 of file DTTFFEDSim.h.

Referenced by DTTFFEDSim(), and fillRawData().