Belle II Software development
SVDUnpackerModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <svd/modules/svdUnpacker/SVDUnpackerModule.h>
10
11
12#include <framework/datastore/DataStore.h>
13#include <framework/datastore/StoreObjPtr.h>
14#include <framework/logging/Logger.h>
15
16#include <boost/crc.hpp> // for boost::crc_basic, boost::augmented_crc
17#define CRC16POLYREV 0x8005 // CRC-16 polynomial, normal representation
18
19#include <arpa/inet.h>
20
21#include <sstream>
22#include <iomanip>
23#include <cstring>
24#include <vector>
25#include <set>
26#include <map>
27#include <utility>
28#include <algorithm>
29
30using namespace std;
31using namespace Belle2;
32using namespace Belle2::SVD;
33
34//-----------------------------------------------------------------
35// Register the Module
36//-----------------------------------------------------------------
37REG_MODULE(SVDUnpacker);
38
39//-----------------------------------------------------------------
40// Implementation
41//-----------------------------------------------------------------
42
43std::string Belle2::SVD::SVDUnpackerModule::m_xmlFileName = std::string("SVDChannelMapping.xml");
44
49{
50 //Set module properties
51 setDescription("Produce SVDShaperDigits from RawSVD. NOTE: only zero-suppressed mode is currently supported!");
53
54 addParam("SVDEventInfo", m_svdEventInfoName, "Name of the SVDEventInfo object", string(""));
55 addParam("rawSVDListName", m_rawSVDListName, "Name of the raw SVD List", string(""));
56 addParam("svdShaperDigitListName", m_svdShaperDigitListName, "Name of the SVDShaperDigits list", string(""));
57 addParam("shutUpFTBError", m_shutUpFTBError,
58 "if >0 is the number of reported FTB header ERRORs before quiet operations. If <0 full log produced.", -1);
59 addParam("FADCTriggerNumberOffset", m_FADCTriggerNumberOffset,
60 "number to be added to the FADC trigger number to match the main trigger number", 0);
61 addParam("svdDAQDiagnosticsListName", m_svdDAQDiagnosticsListName, "Name of the DAQDiagnostics list", string(""));
62 addParam("softwarePipelineAddressEmulation", m_emulatePipelineAddress, "Estimate emulated pipeline address", bool(true));
63 addParam("killDigitsFromUpsetAPVs", m_killUpsetDigits, "Delete digits from upset APVs", bool(false));
64 addParam("silentlyAppend", m_silentAppend, "Append digits to a pre-existing non-empty storeArray", bool(false));
65 addParam("badMappingFatal", m_badMappingFatal, "Throw B2FATAL if there's a wrong payload in the database", bool(false));
66 addParam("UnpackerErrorRate", m_errorRate, "Unpacker will print one error every UnpackerErrorRate", int(1000));
67 addParam("PrintRawData", m_printRaw, "Printing Raw data words for debugging", bool(false));
68}
69
73
75{
76 m_eventMetaDataPtr.isRequired();
77 // Don't panic if no SVD data.
78 m_rawSVD.isOptional(m_rawSVDListName);
79
80 // Register default SVDEventInfo for unpacking Raw Data
82
85
88
89}
90
92{
93 if (!m_mapping.isValid())
94 B2FATAL("no valid SVD Channel Mapping. We stop here.");
95
96 m_wrongFTBcrc = 0;
97 if (m_mapping.hasChanged()) { m_map = std::make_unique<SVDOnlineToOfflineMap>(m_mapping->getFileName()); }
98
99 if (! m_map) { //give up
100 B2ERROR("SVD xml map not loaded." << std::endl <<
101 "No SVDShaperDigit will be produced for this run!");
102 return;
103 }
104
105 //number of FADC boards
106 nFADCboards = m_map->getFADCboardsNumber();
107
108 //passing APV<->FADC mapping from SVDOnlineToOfflineMap object
109 APVmap = &(m_map->APVforFADCmap);
110
111 //setting UnpackerErrorRate factor to use it for BadMapping error suppression
112 m_map->setErrorRate(m_errorRate);
113
116 nUpsetAPVsErrors = -1;
117 nSEURecoveryCase = -1;
120 nFADCMatchErrors = -1;
121 nAPVErrors = -1;
122 nFTBFlagsErrors = -1;
124
126
127 //get the relative time shift
128 if (!m_svdGlobalConfig.isValid())
129 B2FATAL("SVDGlobalConfigParameters not valid!!");
130
131 m_relativeTimeShift = m_svdGlobalConfig->getRelativeTimeShift();
132
133}
134
135#ifndef __clang__
136#pragma GCC diagnostic push
137#pragma GCC diagnostic ignored "-Wstack-usage="
138#endif
140{
141 if (!m_rawSVD || !m_rawSVD.getEntries())
142 return;
143
145 B2WARNING("Unpacking SVDShaperDigits to a non-empty pre-existing \n"
146 << "StoreArray. This can lead to undesired behaviour. At least\n"
147 << "remember to use SVDShaperDigitSorter in your path and \n"
148 << "set the silentlyAppend parameter of SVDUnpacker to true.");
149
150 SVDDAQDiagnostic* currentDAQDiagnostic;
151 vector<SVDDAQDiagnostic*> vDiagnostic_ptr;
152
153 map<SVDShaperDigit, SVDDAQDiagnostic*> diagnosticMap;
154 // Store encountered pipeline addresses with APVs in which they were observed
155 map<unsigned short, set<pair<unsigned short, unsigned short> > > apvsByPipeline;
156
157 if (!m_eventMetaDataPtr.isValid()) { // give up...
158 B2ERROR("Missing valid EventMetaData." << std::endl << "No SVDShaperDigit produced for this event!");
159 return;
160 }
161
162 bool nAPVmatch = true;
163 bool badMapping = false;
164 bool badHeader = false;
165 bool badTrailer = false;
166 bool missedHeader = false;
167 bool missedTrailer = false;
168
169 // flag to set SVDEventInfo once per event
170 bool isSetEventInfo = false;
171
172 //flag to set nAPVsamples in SVDEventInfo once per event
173 bool isSetNAPVsamples = false;
174
175 unsigned short nAPVheaders = 999;
176 set<short> seenAPVHeaders = {};
177
178 unsigned short nEntries_rawSVD = m_rawSVD.getEntries();
179 auto eventNo = m_eventMetaDataPtr->getEvent();
180
181 short fadc = 255, apv = 63;
182
183 unsigned short cntFADCboards = 0;
184 for (unsigned int i = 0; i < nEntries_rawSVD; i++) {
185
186 unsigned int numEntries_rawSVD = m_rawSVD[ i ]->GetNumEntries();
187 for (unsigned int j = 0; j < numEntries_rawSVD; j++) {
188
189 const unsigned short maxNumOfCh = m_rawSVD[i]->GetMaxNumOfCh(j);
190
191 std::vector<unsigned short> nWords;
192 nWords.reserve(maxNumOfCh);
193 std::vector<uint32_t*> data32tab(maxNumOfCh); //vector of pointers
194 for (unsigned int k = 0; k < maxNumOfCh; k++) {
195 nWords.push_back(m_rawSVD[i]->GetDetectorNwords(j, k));
196 data32tab[k] = (uint32_t*)m_rawSVD[i]->GetDetectorBuffer(j, k); // points at the beginning of the 1st buffer
197 }
198
199 unsigned short ftbError = 0;
200 unsigned short trgType = 0;
201 unsigned short trgNumber = 0;
202 unsigned short daqMode = -1;
203 unsigned short daqType = 0;
204 unsigned short cmc1;
205 unsigned short cmc2;
206 unsigned short apvErrors;
207 unsigned short pipAddr;
208 unsigned short ftbFlags = 0;
209 unsigned short apvErrorsOR = 0;
210
211 bool is3sampleData = false;
212 bool is6sampleData = false;
213
214 for (unsigned int buf = 0; buf < maxNumOfCh; buf++) { // loop over 4(COPPER) or 48(PCIe40) buffers
215
216 if (data32tab[buf] == nullptr || nWords.at(buf) == 0) {
217 if (data32tab[buf] != nullptr || nWords.at(buf) != 0) {
218 B2WARNING("Invalid combination of buffer pointer and nWords:" <<
219 LogVar("COPPER/PCIe40 ID", i) <<
220 LogVar("COPPER/PCIe40 Entry", j) <<
221 LogVar("COPPER/PCIe40 Channel", buf) <<
222 LogVar("data32tab[buf]", data32tab[buf]) <<
223 LogVar("nWords[buf]", nWords.at(buf)));
224 }
225 continue;
226 }
227 if (m_printRaw) printB2Debug(data32tab[buf], data32tab[buf], &data32tab[buf][nWords.at(buf) - 1], nWords.at(buf));
228
229 cntFADCboards++;
230
231 missedHeader = false;
232 missedTrailer = false;
233
234 uint32_t* data32_it = data32tab[buf];
235 short strip, sample[6];
236 vector<uint32_t> crc16vec;
237
238 //reset value for headers and trailers check
240
241 for (; data32_it != &data32tab[buf][nWords.at(buf)]; data32_it++) {
242 m_data32 = *data32_it; //put current 32-bit frame to union
243
244 if (m_data32 == 0xffaa0000) { // first part of FTB header
245 crc16vec.clear(); // clear the input container for crc16 calculation
246 crc16vec.push_back(m_data32);
247
248 seenHeadersAndTrailers |= 0x1; // we found FTB header
249
250 data32_it++; // go to 2nd part of FTB header
251 crc16vec.push_back(*data32_it);
252
253 m_data32 = *data32_it; //put the second 32-bit frame to union
254
255 ftbError = m_FTBHeader.errorsField;
256
257 if (ftbError != 240) {
259
261 switch (ftbError - 240) {
262 case 3:
263 B2ERROR("FADC Event Number is different from (FTB & TTD) Event Numbers");
264 break;
265 case 5:
266 B2ERROR("TTD Event Number is different from (FTB & FADC) Event Numbers");
267 break;
268 case 6:
269 B2ERROR("FTB Event Number is different from (TTD & FADC) Event Numbers");
270 break;
271 case 7:
272 B2ERROR("(FTB, TTD & FADC) Event Numbers are different from each other");
273 break;
274 default:
275 B2ERROR("Problem with errorsField variable in FTB Header" << LogVar("abnormal value", ftbError));
276 }
277 }
278 }
279
280 if (m_FTBHeader.eventNumber !=
281 (eventNo & 0xFFFFFF)) {
284 m_shutUpFTBError -= 1;
285 B2ERROR("Event number mismatch detected! The event number given by EventMetaData object is different from the one in the FTB Header."
286 << LogVar("Expected event number & 0xFFFFFF",
287 (eventNo & 0xFFFFFF)) << LogVar("Event number in the FTB", m_FTBHeader.eventNumber));
288 }
289 }
290
291 continue;
292 } // is FTB Header
293
294 crc16vec.push_back(m_data32);
295
296 if (m_MainHeader.check == 6) { // FADC header
297
298 seenHeadersAndTrailers |= 0x2; //we found FADC Header
299
300 fadc = m_MainHeader.FADCnum;
301 trgType = m_MainHeader.trgType;
302 trgNumber = m_MainHeader.trgNumber;
303 daqMode = m_MainHeader.DAQMode;
304 daqType = m_MainHeader.DAQType;
305
306 //Let's add run-dependent info: daqMode="11" in case of 3-mixed-6 sample acquisition mode.
307 if (daqType) daqMode = 3;
308
309 nAPVheaders = 0; // start counting APV headers for this FADC
310 nAPVmatch = true; //assume correct # of APV headers
311 badMapping = false; //assume correct mapping
312 badHeader = false;
313 badTrailer = false;
314
315 is3sampleData = false;
316 is6sampleData = false;
317
318 if (daqMode == 0) B2ERROR("SVDDataFormatCheck: the event " << eventNo <<
319 " is apparently taken with 1-sample mode, this is not expected.");
320 if (daqMode == 1) is3sampleData = true;
321 if (daqMode == 2) is6sampleData = true;
322
323 if (
324 m_MainHeader.trgNumber !=
325 ((eventNo - m_FADCTriggerNumberOffset) & 0xFF)) {
326
329 B2ERROR("Event number mismatch detected! The event number given by EventMetaData object is different from the one in the FADC Header. "
330 << LogVar("Event number", eventNo) << LogVar("FADC", fadc) << LogVar("Trigger number LSByte reported by the FADC",
331 m_MainHeader.trgNumber) << LogVar("+ offset", m_FADCTriggerNumberOffset) << LogVar("expected", (eventNo & 0xFF)));
332 badHeader = true;
333 }
334
335 // create SVDModeByte object from MainHeader vars
336 m_SVDModeByte = SVDModeByte(m_MainHeader.runType, 0, daqMode, m_MainHeader.trgTiming);
337
338 // create SVDEventInfo and fill it with SVDModeByte & SVDTriggerType objects
339 if (!isSetEventInfo) {
341 m_svdEventInfoPtr.create();
342 m_svdEventInfoPtr->setModeByte(m_SVDModeByte);
343 m_svdEventInfoPtr->setTriggerType(m_SVDTriggerType);
344 m_svdEventInfoPtr->setAPVClock(m_hwClock);
345
346 //set relative time shift
347 m_svdEventInfoPtr->setRelativeShift(m_relativeTimeShift);
348 // set X-talk info online from Raw Data
349 m_svdEventInfoPtr->setCrossTalk(m_MainHeader.xTalk);
350
351 isSetEventInfo = true;
352 } else { // let's check if the current SVDModeByte and SVDTriggerType are consistent with the one stored in SVDEventInfo
353 if (m_SVDModeByte != m_svdEventInfoPtr->getModeByte()) {m_svdEventInfoPtr->setMatchModeByte(false); badHeader = true; nEventInfoMatchErrors++;}
354 if (trgType != (m_svdEventInfoPtr->getTriggerType()).getType()) { m_svdEventInfoPtr->setMatchTriggerType(false); badHeader = true; nEventInfoMatchErrors++;}
355 }
356 } // is FADC header
357
358 if (m_APVHeader.check == 2) { // APV header
359
360 nAPVheaders++;
361 apv = m_APVHeader.APVnum;
362 seenAPVHeaders.insert(apv);
363
364 cmc1 = m_APVHeader.CMC1;
365 cmc2 = m_APVHeader.CMC2;
366 apvErrors = m_APVHeader.apvErr;
367 pipAddr = m_APVHeader.pipelineAddr;
368
369 if (apvErrors != 0) {
370 nAPVErrors++;
371 if (!(nAPVErrors % m_errorRate)
372 or nAPVErrors < 100) B2ERROR("APV error has been detected." << LogVar("FADC", fadc) << LogVar("APV", apv) << LogVar("Error value",
373 apvErrors));
374 }
375 // temporary SVDDAQDiagnostic object (no info from trailers and APVmatch code)
376 currentDAQDiagnostic = m_storeDAQDiagnostics.appendNew(trgNumber, trgType, pipAddr, cmc1, cmc2, apvErrors, ftbError, true,
377 nAPVmatch,
378 badHeader, missedHeader, missedTrailer,
379 fadc, apv);
380 vDiagnostic_ptr.push_back(currentDAQDiagnostic);
381
382 apvsByPipeline[pipAddr].insert(make_pair(fadc, apv));
383
384 // Let's check if the data frame does not come after APV header with piplAddr = 255 (special SEU recovery pseudo-data)
385 if (pipAddr == 255) {
386 data32_it++;
387 m_data32 = *data32_it;
388 if (m_data_A.check == 0) B2ERROR("The strip data frame is detected for the special SEU recovery mode. " << LogVar("pipline address",
389 pipAddr) << "This is unexpected!");
390 data32_it--;
391 m_data32 = *data32_it;
392 }
393
394 } //is APV Header
395
396 if (m_data_A.check == 0) { // data
397 strip = m_data_A.stripNum;
398
399 sample[0] = m_data_A.sample1;
400 sample[1] = m_data_A.sample2;
401 sample[2] = m_data_A.sample3;
402
403 sample[3] = 0;
404 sample[4] = 0;
405 sample[5] = 0;
406
407 // Let's check the next rawdata word to determine if we acquired 3 or 6 sample
408 data32_it++;
409 m_data32 = *data32_it;
410
411 if (m_data_B.check == 0 && strip == m_data_B.stripNum) { // 2nd data frame with the same strip number -> six samples
412
413 if (!isSetNAPVsamples) {
414 m_svdEventInfoPtr->setNSamples(6);
415 isSetNAPVsamples = true;
416 } else {
417 if (is3sampleData)
418 B2ERROR("DAQMode value (indicating 3-sample acquisition mode) doesn't correspond to the actual number of samples (6) in the data! The data might be corrupted!");
419 }
420
421 crc16vec.push_back(m_data32);
422
423 sample[3] = m_data_B.sample4;
424 sample[4] = m_data_B.sample5;
425 sample[5] = m_data_B.sample6;
426 }
427
428 else { // three samples
429 data32_it--;
430 m_data32 = *data32_it;
431
432 if (!isSetNAPVsamples) {
433 m_svdEventInfoPtr->setNSamples(3);
434 isSetNAPVsamples = true;
435 } else {
436 if (is6sampleData)
437 B2ERROR("DAQMode value (indicating 6-sample acquisition mode) doesn't correspond to the actual number of samples (3) in the data! The data might be corrupted!");
438 }
439 }
440
441 // Generating SVDShaperDigit object
442 SVDShaperDigit* newShaperDigit = m_map->NewShaperDigit(fadc, apv, strip, sample, 0.0);
443 if (newShaperDigit) {
444 diagnosticMap.insert(make_pair(*newShaperDigit, currentDAQDiagnostic));
445 delete newShaperDigit;
446 } else if (m_badMappingFatal) {
447 B2FATAL("Respective FADC/APV combination not found -->> incorrect payload in the database! ");
448 } else {
449 badMapping = true;
450 }
451
452 } //is data frame
453
454
455 if (m_FADCTrailer.check == 14) { // FADC trailer
456
457 seenHeadersAndTrailers |= 0x4; // we found FADC trailer
458
459 //additional check if we have a faulty/fake FADC that is not in the map
460 if (APVmap->find(fadc) == APVmap->end()) badMapping = true;
461
462 //comparing number of APV chips and the number of APV headers, for the current FADC
463 unsigned short nAPVs = APVmap->count(fadc);
464
465 if (nAPVheaders == 0) {
466 currentDAQDiagnostic = m_storeDAQDiagnostics.appendNew(0, 0, 0, 0, 0, 0, ftbError, true, nAPVmatch, badHeader, 0, 0, fadc, 0);
467 vDiagnostic_ptr.push_back(currentDAQDiagnostic);
468 }
469
470 if (nAPVs != nAPVheaders) {
471 // There is an APV missing, detect which it is.
472 for (const auto& fadcApv : *APVmap) {
473 if (fadcApv.first != fadc) continue;
474 if (seenAPVHeaders.find(fadcApv.second) == seenAPVHeaders.end()) {
475 // We have a missing APV. Look if it is a known one.
476 auto missingRec = m_missingAPVs.find(make_pair(fadcApv.first, fadcApv.second));
477 if (missingRec != m_missingAPVs.end()) {
478 // This is known to be missing, so keep quiet and just update event counters
479 if (missingRec->second.first > eventNo)
480 missingRec->second.first = eventNo;
481 if (missingRec->second.second < eventNo)
482 missingRec->second.second = eventNo;
483 } else {
484 // We haven't seen this previously.
486 m_missingAPVs.insert(make_pair(
487 make_pair(fadcApv.first, fadcApv.second),
488 make_pair(eventNo, eventNo)
489 ));
490 if (!(nMissingAPVsErrors % m_errorRate)) B2ERROR("missing APV header! " << LogVar("Event number", eventNo) << LogVar("APV",
491 int(fadcApv.second)) << LogVar("FADC",
492 int(fadcApv.first)));
493 }
494 }
495 }
496 nAPVmatch = false;
497 } // is nAPVs != nAPVheaders
498
499 seenAPVHeaders.clear();
500
501 ftbFlags = m_FADCTrailer.FTBFlags;
502 if ((ftbFlags >> 5) != 0) badTrailer = true;
503 if (ftbFlags != 0) {
505 if (!(nFTBFlagsErrors % m_errorRate) or nFTBFlagsErrors < 100) {
506 B2ERROR(" FTB Flags variable has an active error bit(s)" << LogVar("on FADC number", fadc));
507
508 if (ftbFlags & 16) B2ERROR("----> CRC error has been detected. Data might be corrupted!");
509 if (ftbFlags & 8) B2ERROR("----> Bad Event indication has been detected. Data might be corrupted!");
510 if (ftbFlags & 4) B2ERROR("----> Double Header has been detected. Data might be corrupted!");
511 if (ftbFlags & 2) B2ERROR("----> Time Out has been detected. Data might be corrupted!");
512 if (ftbFlags & 1) B2ERROR("----> Event Too Long! Data might be corrupted!");
513 }
514 }
515
516 apvErrorsOR = m_FADCTrailer.apvErrOR;
517
518
519 }// is FADC trailer
520
521 if (m_FTBTrailer.controlWord == 0xff55) {// FTB trailer
522
523 seenHeadersAndTrailers |= 0x8; // we found FTB trailer
524
525 //check CRC16
526 crc16vec.pop_back();
527 unsigned short iCRC = crc16vec.size();
528 std::vector<uint32_t> crc16input;
529 crc16input.reserve(iCRC);
530
531 for (unsigned short icrc = 0; icrc < iCRC; icrc++)
532 crc16input.push_back(htonl(crc16vec.at(icrc)));
533
534 //verify CRC16
535 boost::crc_basic<16> bcrc(0x8005, 0xffff, 0, false, false);
536 bcrc.process_bytes(crc16input.data(), crc16input.size() * sizeof(uint32_t));
537 unsigned int checkCRC = bcrc.checksum();
538
539 if (checkCRC != m_FTBTrailer.crc16) {
540 B2WARNING("FTB CRC16 checksum DOES NOT MATCH" << LogVar("for FADC no.", fadc));
542 }
543
544 } // is FTB trailer
545
546 } // end loop over 32-bit frames in each buffer
547
548 //Let's check if all the headers and trailers were in place in the last frame
549 if (seenHeadersAndTrailers != 0xf) {
550 if (!(seenHeadersAndTrailers & 1)) {B2ERROR("Missing FTB Header is detected. SVD data might be corrupted!" << LogVar("Event number", eventNo) << LogVar("FADC", fadc)); missedHeader = true;}
551 if (!(seenHeadersAndTrailers & 2)) {B2ERROR("Missing FADC Header is detected -> related FADC number couldn't be retrieved. SVD data might be corrupted! " << LogVar("Event number", eventNo) << LogVar("previous FADC", fadc)); missedHeader = true;}
552 if (!(seenHeadersAndTrailers & 4)) {B2ERROR("Missing FADC Trailer is detected. SVD data might be corrupted!" << LogVar("Event number", eventNo) << LogVar("FADC", fadc)); missedTrailer = true;}
553 if (!(seenHeadersAndTrailers & 8)) {B2ERROR("Missing FTB Trailer is detected. SVD data might be corrupted!" << LogVar("Event number", eventNo) << LogVar("FADC", fadc)); missedTrailer = true;}
554 }
555
556 for (auto p : vDiagnostic_ptr) {
557 // adding remaining info to Diagnostic object
558 p->setFTBFlags(ftbFlags);
559 p->setApvErrorOR(apvErrorsOR);
560 p->setAPVMatch(nAPVmatch);
561 p->setBadMapping(badMapping);
562 p->setBadTrailer(badTrailer);
563 p->setMissedHeader(missedHeader);
564 p->setMissedTrailer(missedTrailer);
565 }
566
567 vDiagnostic_ptr.clear();
568
569 } // end iteration on 4(COPPER)/48(PCIe40) data buffers
570
571 } // end event loop
572
573 }// end loop over RawSVD objects
574
575 // Check the number of FADC boards
576 if (cntFADCboards != nFADCboards) { // nFADCboards=52
578 if (!(nFADCMatchErrors % m_errorRate)) B2ERROR("Number of data objects in rawSVD do not match the number of FADC boards" <<
579 LogVar("# of data objects in rawSVD",
580 cntFADCboards) << LogVar("# of FADCs", nFADCboards) << LogVar("Event number", eventNo));
581
582 // We override all FADCMatch fields in diagnostics and set it to false.
583 for (auto& p : m_storeDAQDiagnostics) {
584 p.setFADCMatch(false);
585 }
586 }
587
588
589 // Check if we have special data for SEU recovery mode (pipAddr = 0xFF)
590 // If so, let's monitor affected FADC/APV's and the event range with seuRecMap
591 auto itPtr = apvsByPipeline.find(255);
592 if (itPtr != apvsByPipeline.end()) {
593 for (const auto& fadcApv : itPtr->second) {
594
595 auto seuRec = m_seuRecMap.find(make_pair(fadcApv.first, fadcApv.second));
596 if (seuRec != m_seuRecMap.end()) {
597 if (seuRec->second.first > eventNo)
598 seuRec->second.first = eventNo;
599 if (seuRec->second.second < eventNo)
600 seuRec->second.second = eventNo;
601 } else {
603 m_seuRecMap.insert(make_pair(
604 make_pair(fadcApv.first, fadcApv.second),
605 make_pair(eventNo, eventNo)
606 ));
608 B2ERROR("Special Recovery Data (Dummy APV Header) found due to the detection of Single Event Upset (SEU)!!!" << LogVar("APV",
609 int(fadcApv.second)) << LogVar("FADC", int(fadcApv.first)) << LogVar("Event number", eventNo));
610
611 }
612 for (auto& pp : m_storeDAQDiagnostics) {
613 if (pp.getFADCNumber() == fadcApv.first and pp.getAPVNumber() == fadcApv.second)
614 pp.setSEURecoData(true);
615 }
616 }
617 }
618
619
620 // Detect upset APVs and report/treat
621 auto majorAPV = max_element(apvsByPipeline.begin(), apvsByPipeline.end(),
622 [](const decltype(apvsByPipeline)::value_type & p1,
623 const decltype(apvsByPipeline)::value_type & p2) -> bool
624 { return p1.second.size() < p2.second.size(); }
625 );
626 // We set emuPipelineAddress fields in diagnostics to this.
628 for (auto& p : m_storeDAQDiagnostics)
629 p.setEmuPipelineAddress(majorAPV->first);
630
631 unsigned short apvsByPipelineSize = apvsByPipeline.size();
632 if (majorAPV->first == 255) apvsByPipelineSize = 1;
633
634 // And report any upset apvs or update records
635 if (apvsByPipelineSize > 1)
636 for (const auto& p : apvsByPipeline) {
637 if (p.first == majorAPV->first or p.first == 255) continue;
638 for (const auto& fadcApv : p.second) {
639 // We have an upset APV. Look if it is a known one.
640 auto upsetRec = m_upsetAPVs.find(make_pair(fadcApv.first, fadcApv.second));
641 if (upsetRec != m_upsetAPVs.end()) {
642 // This is known to be upset, so keep quiet and update event counters
643 if (upsetRec->second.first > eventNo)
644 upsetRec->second.first = eventNo;
645 if (upsetRec->second.second < eventNo)
646 upsetRec->second.second = eventNo;
647 } else {
648 // We haven't seen this one previously.
650 m_upsetAPVs.insert(make_pair(
651 make_pair(fadcApv.first, fadcApv.second),
652 make_pair(eventNo, eventNo)
653 ));
654
655 if (!(nUpsetAPVsErrors % m_errorRate)) B2ERROR("Upset APV detected!!!" << LogVar("APV", int(fadcApv.second)) << LogVar("FADC",
656 int(fadcApv.first)) << LogVar("Event number", eventNo));
657 }
658
659 for (auto& pp : m_storeDAQDiagnostics) {
660
661 if (pp.getFADCNumber() == fadcApv.first and pp.getAPVNumber() == fadcApv.second)
662 pp.setUpsetAPV(true);
663 }
664
665 }
666 }
667
668 // Here we can delete digits coming from upset APVs. We detect them by comparing
669 // actual and emulated pipeline address fields in DAQDiagnostics.
670 for (auto& p : diagnosticMap) {
671
672 if ((m_killUpsetDigits && p.second->getPipelineAddress() != p.second->getEmuPipelineAddress()) || p.second->getFTBError() != 240
673 || p.second->getFTBFlags() || p.second->getAPVError() || !(p.second->getAPVMatch()) || !(p.second->getFADCMatch())
674 || p.second->getBadHeader()
675 || p.second->getBadMapping() || p.second->getUpsetAPV() || p.second->getMissedHeader() || p.second->getMissedTrailer()) continue;
676 m_storeShaperDigits.appendNew(p.first);
677 }
678
679 if (!m_svdEventInfoPtr->getMatchTriggerType()) {if (!(nEventInfoMatchErrors % m_errorRate) or nEventInfoMatchErrors < 200) B2WARNING("Inconsistent SVD Trigger Type value for: " << LogVar("Event number", eventNo));}
680 if (!m_svdEventInfoPtr->getMatchModeByte()) {if (!(nEventInfoMatchErrors % m_errorRate) or nEventInfoMatchErrors < 200) B2WARNING("Inconsistent SVD ModeByte object for: " << LogVar("Event number", eventNo));}
681
682
683} //end event function
684#ifndef __clang__
685#pragma GCC diagnostic pop
686#endif
687
689{
690
691 // Summary report on missing APVs
692 if (m_missingAPVs.size() > 0) {
693 B2WARNING("SVDUnpacker summary 1: Missing APVs");
694 for (const auto& miss : m_missingAPVs)
695 B2WARNING(LogVar("Missing APV", miss.first.second) << LogVar("FADC", miss.first.first) << LogVar("since event",
696 miss.second.first) << LogVar("to event", miss.second.second));
697 }
698 if (m_upsetAPVs.size() > 0) {
699 B2WARNING("SVDUnpacker summary 2: Upset APVs");
700 for (const auto& upst : m_upsetAPVs)
701 B2WARNING(LogVar("Upset APV", upst.first.second) << LogVar("FADC", upst.first.first) <<
702 LogVar("since event", upst.second.first) << LogVar("to event", upst.second.second));
703 }
704
705 if (m_seuRecMap.size() > 0) {
706 B2WARNING("SVDUnpacker summary 2: Special SEU Recovery Data");
707 for (const auto& seu : m_seuRecMap)
708 B2WARNING(LogVar("SEU recovery data: APV", seu.first.second) << LogVar("FADC", seu.first.first) <<
709 LogVar("since event", seu.second.first) << LogVar("to event", seu.second.second));
710 }
711
712}
713
714
715// additional printing function
716void SVDUnpackerModule::printB2Debug(uint32_t* data32, uint32_t* data32_min, uint32_t* data32_max, int nWords)
717{
718
719 uint32_t* min = std::max((data32 - nWords), data32_min);
720 uint32_t* max = std::min((data32 + nWords), data32_max);
721
722 size_t counter{0};
723 std::stringstream os;
724 os << std::hex << std::setfill('0');
725 for (uint32_t* ptr = min; ptr <= max; ++ptr) {
726 os << std::setw(8) << *ptr;
727 if (++counter % 10 == 0) os << std::endl;
728 else os << " ";
729 }
730
731 os << std::endl;
732 B2INFO(os.str());
733 return;
734
735}
@ c_ErrorIfAlreadyRegistered
If the object/array was already registered, produce an error (aborting initialisation).
Definition DataStore.h:72
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
Class to store SVD DAQ diagnostic information.
Class to store SVD mode information.
Definition SVDModeByte.h:69
The SVD ShaperDigit class.
Class to store Trigger Type information.
int m_wrongFTBcrc
FTB CRC no-Match counter.
StoreArray< SVDDAQDiagnostic > m_storeDAQDiagnostics
SVDDAQDiagnostic array.
virtual ~SVDUnpackerModule()
Destructor of the module.
unsigned short seenHeadersAndTrailers
this 4-bits value should be 1111 if no headers/trailers are missing
int m_shutUpFTBError
regulates the number of "Event number mismatch" errors reported
MainHeader m_MainHeader
Implementation of FADC Header.
StoreArray< RawSVD > m_rawSVD
output for RawSVD
FADCTrailer m_FADCTrailer
Implementation of FADC Trailer.
StoreObjPtr< SVDEventInfo > m_svdEventInfoPtr
SVDEventInfo output per event.
virtual void initialize() override
Initializes the Module.
StoreArray< SVDShaperDigit > m_storeShaperDigits
SVDShaperDigit array.
virtual void event() override
event
SVDUnpackerModule()
Constructor of the module.
int nAPVErrors
counter of APV errors
FTBHeader m_FTBHeader
Implementation of FTB Header.
virtual void endRun() override
end run
std::string m_rawSVDListName
RawSVD StoreArray name.
DBObjPtr< PayloadFile > m_mapping
pointer to the payload with the mapping
int nFADCMatchErrors
counter of FADC boards =/= n of RawData objects errors
bool m_killUpsetDigits
Optionally, we can kill digits coming from upset APVs right in the unpacker.
data_B m_data_B
Implementation of 2nd data word.
int m_FADCTriggerNumberOffset
FADC Trigger Offset.
std::map< std::pair< unsigned short, unsigned short >, std::pair< std::size_t, std::size_t > > m_missingAPVs
Map to store a list of missing APVs.
int nMissingAPVsErrors
counter of missing APVs errors
static std::string m_xmlFileName
XML filename.
int nEventMatchErrors
counter of Event match errors
FTBTrailer m_FTBTrailer
Implementation of FTB Trailer.
int m_errorRate
The parameter that indicates what fraction of B2ERRORs messages should be suppressed to not overload ...
APVHeader m_APVHeader
Implementation of APV Header.
int nUpsetAPVsErrors
counter of upset APV errors
void printB2Debug(uint32_t *data32, uint32_t *data32_min, uint32_t *data32_max, int nWords)
additional function that prints raw data words
virtual void beginRun() override
begin run
std::string m_svdShaperDigitListName
SVDShaperDigit StoreArray name.
DBObjPtr< HardwareClockSettings > m_hwClock
system clock
std::string m_svdEventInfoName
SVDEventInfo name.
DBObjPtr< SVDGlobalConfigParameters > m_svdGlobalConfig
SVDGlobal Configuration payload.
bool m_printRaw
Optionally we can get printout of Raw Data words.
int nEventInfoMatchErrors
counter of inconsistencies in SVDEventInfo within an event
bool m_badMappingFatal
Optionally we can stop the unpacking if there is a missing APV/FADC combination in the mapping -> wro...
int nErrorFieldErrors
counter of event mismatch errors in FTB's ErrorField
std::map< std::pair< unsigned short, unsigned short >, std::pair< std::size_t, std::size_t > > m_seuRecMap
Map to store a list of APVs for special data for SEU recovery.
data_A m_data_A
Implementation of 1st data word.
uint32_t m_data32
Input 32-bit data word.
StoreObjPtr< EventMetaData > m_eventMetaDataPtr
Required input for EventMetaData.
bool m_silentAppend
Silently append new SVDShaperDigits to a pre-existing non-empty SVDShaperDigits storeArray.
std::unordered_multimap< unsigned char, unsigned char > * APVmap
pointer to APVforFADCmap filled by mapping procedure
SVDTriggerType m_SVDTriggerType
SVDTriggerType object.
SVDModeByte m_SVDModeByte
instance of SVDModeByte for the event
std::unique_ptr< SVDOnlineToOfflineMap > m_map
Pointer to online-to-offline map.
int m_relativeTimeShift
latency difference between the 3- and 6-sample acquired events in usint of APV clock / 4,...
std::map< std::pair< unsigned short, unsigned short >, std::pair< std::size_t, std::size_t > > m_upsetAPVs
Map to store a list of upset APVs.
int nFTBFlagsErrors
counter of errors in FTBFlags variable
int nSEURecoveryCase
counter of SEU Special Recovery data cases
bool m_emulatePipelineAddress
Software emulation of pipeline address This is a replacement of hardware pipeline address emulation.
std::string m_svdDAQDiagnosticsListName
SVDDAQDiagnostic StoreArray name.
int nTriggerMatchErrors
counters for specific ERRORS produced by the Unpacker
unsigned short nFADCboards
how many FADCs we have
Class to store variables with their name which were sent to the logging service.
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Namespace to encapsulate code needed for simulation and reconstrucion of the SVD.
Abstract base class for different kinds of events.
STL namespace.