Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:10:57

0001 #include "EventFilter/RctRawToDigi/plugins/RctRawToDigi.h"
0002 
0003 // System headers
0004 #include <vector>
0005 #include <sstream>
0006 #include <iostream>
0007 
0008 // Framework headers
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/Framework/interface/MakerMacros.h"
0013 
0014 // Raw data collection headers
0015 #include "DataFormats/FEDRawData/interface/FEDHeader.h"
0016 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0017 #include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h"
0018 #include "DataFormats/FEDRawData/interface/FEDTrailer.h"
0019 
0020 #include "EventFilter/L1TRawToDigi/interface/AMCSpec.h"
0021 #include "EventFilter/L1TRawToDigi/interface/Block.h"
0022 #include "EventFilter/L1TRawToDigi/interface/PackingSetup.h"
0023 
0024 // Namespace resolution
0025 using edm::LogError;
0026 using edm::LogWarning;
0027 using std::cout;
0028 using std::dec;
0029 using std::endl;
0030 using std::hex;
0031 using std::string;
0032 using std::vector;
0033 
0034 RctRawToDigi::RctRawToDigi(const edm::ParameterSet& iConfig)
0035     : inputLabel_(iConfig.getParameter<edm::InputTag>("inputLabel")),
0036       fedId_(iConfig.getUntrackedParameter<int>("rctFedId", FEDNumbering::MINRCTFEDID)),
0037       verbose_(iConfig.getUntrackedParameter<bool>("verbose", false)) {
0038   LogDebug("RCT") << "RctRawToDigi will unpack FED Id " << fedId_;
0039   /** Register Products **/
0040   // RCT collections
0041   produces<L1CaloEmCollection>();
0042   produces<L1CaloRegionCollection>();
0043 
0044   // Error collection
0045   consumes<FEDRawDataCollection>(inputLabel_);
0046 }
0047 
0048 RctRawToDigi::~RctRawToDigi() {
0049   // do anything here that needs to be done at destruction time
0050   // (e.g. close files, deallocate resources etc.)
0051   //delete formatTranslator_;
0052 }
0053 
0054 // ------------ method called to produce the data  ------------
0055 void RctRawToDigi::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0056   using namespace edm;
0057 
0058   // Instantiate all the collections the unpacker needs; puts them in event when this object goes out of scope.
0059   std::unique_ptr<RctUnpackCollections> colls(new RctUnpackCollections(iEvent));
0060 
0061   // get raw data collection
0062   edm::Handle<FEDRawDataCollection> feds;
0063   iEvent.getByLabel(inputLabel_, feds);
0064 
0065   // if raw data collection is present, check the headers and do the unpacking
0066   if (feds.isValid()) {
0067     const FEDRawData& rctRcd = feds->FEDData(fedId_);
0068 
0069     //Check FED size
0070     LogDebug("RCT") << "Upacking FEDRawData of size " << std::dec << rctRcd.size();
0071 
0072     //check header size
0073     if (rctRcd.size() < sLinkHeaderSize_ + sLinkTrailerSize_ + amc13HeaderSize_ + amc13TrailerSize_ + MIN_DATA) {
0074       if (rctRcd.size() > 0) {
0075         LogError("L1T") << "Cannot unpack: empty/invalid L1T raw data (size = " << rctRcd.size() << ") for ID "
0076                         << fedId_ << ". Returning empty collections!";
0077       } else {
0078         // there's no MC packer for this payload, so totally expected that sometimes it will be absent... no warning issued.
0079       }
0080       //continue;
0081       return;
0082     }
0083 
0084     // do the unpacking
0085     unpack(rctRcd, iEvent, colls.get());
0086 
0087   } else {
0088     LogError("L1T") << "Cannot unpack: no collection found";
0089 
0090     return;
0091   }
0092 }
0093 
0094 void RctRawToDigi::unpack(const FEDRawData& d, edm::Event& e, RctUnpackCollections* const colls) {
0095   const unsigned char* data = d.data();  // The 8-bit wide raw-data array.
0096   // Data offset - starts at 16 as there is a 64-bit S-Link header followed
0097   // by a 64-bit software-controlled header (for pipeline format version
0098   // info that is not yet used).
0099 
0100   FEDHeader header(data);
0101 
0102   if (header.check()) {
0103     LogDebug("L1T") << "Found SLink header:"
0104                     << " Trigger type " << header.triggerType() << " L1 event ID " << header.lvl1ID() << " BX Number "
0105                     << header.bxID() << " FED source " << header.sourceID() << " FED version " << header.version();
0106   } else {
0107     LogWarning("L1T") << "Did not find a valid SLink header!";
0108   }
0109 
0110   FEDTrailer trailer(data + (d.size() - sLinkTrailerSize_));
0111 
0112   if (trailer.check()) {
0113     LogDebug("L1T") << "Found SLink trailer:"
0114                     << " Length " << trailer.fragmentLength() << " CRC " << trailer.crc() << " Status "
0115                     << trailer.evtStatus() << " Throttling bits " << trailer.ttsBits();
0116   } else {
0117     LogWarning("L1T") << "Did not find a SLink trailer!";
0118   }
0119 
0120   unpackCTP7(reinterpret_cast<const uint32_t*>(data), 0, sizeof(data), colls);
0121 }
0122 
0123 void RctRawToDigi::unpackCTP7(const uint32_t* data,
0124                               const unsigned block_id,
0125                               const unsigned size,
0126                               RctUnpackCollections* const colls) {
0127   //offset from 6 link header words
0128   uint32_t of = 6;
0129   LogDebug("L1T") << "Block ID  = " << block_id << " size = " << size;
0130 
0131   CTP7Format ctp7Format;
0132   RctDataDecoder rctDataDecoder;
0133   uint32_t nBXTemp = 0;
0134   uint32_t ctp7FWVersion;
0135   uint32_t L1ID, L1aBCID;
0136   std::vector<RCTInfo> allCrateRCTInfo[5];
0137 
0138   L1ID = data[1 + of];                          // extract the L1 ID number
0139   L1aBCID = data[5 + of] & 0x00000FFF;          // extract the BCID number of L1A
0140   nBXTemp = (data[5 + of] & 0x00FF0000) >> 16;  // extract number of BXs readout per L1A
0141   ctp7FWVersion = data[4 + of];
0142 
0143   if (nBXTemp != 1 && nBXTemp != 3 && nBXTemp != 5)
0144     nBXTemp = 1;
0145 
0146   const uint32_t nBX = nBXTemp;
0147 
0148   LogDebug("L1T") << "CTP7 L1ID = " << L1ID << " L1A BCID = " << L1aBCID << " BXs in capture = " << nBX
0149                   << " CTP7 DAQ FW = " << ctp7FWVersion;
0150 
0151   struct link_data {
0152     bool even;
0153     unsigned int crateID;
0154     unsigned int ctp7LinkNumber;
0155     std::vector<unsigned int> uint;
0156   };
0157 
0158   //nBX max 5, nLinks max 36 [nBX][nLinks]
0159   link_data allLinks[5][36];
0160   const uint32_t NLinks = ctp7Format.NLINKS;
0161   assert(NLinks <= 36);
0162 
0163   //change this implementation
0164   uint32_t iDAQBuffer = 0;
0165 
0166   //Step 1: Grab all data from ctp7 buffer and put into link format
0167   for (unsigned int iLink = 0; iLink < NLinks; iLink++) {
0168     iDAQBuffer = of + ctp7Format.EVENT_HEADER_WORDS +
0169                  iLink * (ctp7Format.CHANNEL_HEADER_WORDS + nBX * ctp7Format.CHANNEL_DATA_WORDS_PER_BX);
0170 
0171     //first decode linkID
0172     uint32_t linkID = data[iDAQBuffer++];
0173     uint32_t tmp = data[iDAQBuffer++];
0174     uint32_t CRCErrCnt = tmp & 0x0000FFFF;
0175     //uint32_t linkStatus = (tmp & 0xFFFF0000) >> 16;
0176 
0177     uint32_t crateID = 0;
0178     uint32_t expectedCrateID = 0;
0179     bool even = false;
0180     bool expectedEven = false;
0181 
0182     //getExpected Link ID
0183     rctDataDecoder.getExpectedLinkID(iLink, expectedCrateID, expectedEven);
0184     //getDecodedLink ID
0185     rctDataDecoder.decodeLinkID(linkID, crateID, even);
0186 
0187     //Perform a check to see if the link ID is as expected, if not then report an error but continue unpacking
0188     if (expectedCrateID != crateID || even != expectedEven) {
0189       LogError("L1T") << "Expected Crate ID " << expectedCrateID << " expectedEven " << expectedEven
0190                       << "does not match actual Crate ID " << crateID << " even " << even;
0191     }
0192 
0193     if (CRCErrCnt != 0)
0194       LogError("L1T") << "WARNING CRC ErrorFound linkID " << linkID << " expected crateID " << expectedCrateID;
0195 
0196     // Loop over multiple BX
0197     for (uint32_t iBX = 0; iBX < nBX; iBX++) {
0198       allLinks[iBX][iLink].uint.reserve(6);
0199       allLinks[iBX][iLink].ctp7LinkNumber = iLink;
0200       allLinks[iBX][iLink].crateID = expectedCrateID;
0201       allLinks[iBX][iLink].even = expectedEven;
0202 
0203       //Notice 6 words per BX
0204       for (unsigned int iWord = 0; iWord < 6; iWord++) {
0205         allLinks[iBX][iLink].uint.push_back(data[iDAQBuffer + iWord + iBX * 6]);
0206       }
0207     }
0208   }
0209 
0210   //Step 2: Dynamically match links and create RCTInfo Objects
0211   uint32_t nCratesFound = 0;
0212   for (unsigned int iCrate = 0; iCrate < 18; iCrate++) {
0213     bool foundEven = false, foundOdd = false;
0214     link_data even[5];
0215     link_data odd[5];
0216 
0217     for (unsigned int iLink = 0; iLink < NLinks; iLink++) {
0218       if ((allLinks[0][iLink].crateID == iCrate) && (allLinks[0][iLink].even == true)) {
0219         foundEven = true;
0220         for (unsigned int iBX = 0; iBX < nBX; iBX++)
0221           even[iBX] = allLinks[iBX][iLink];
0222       } else if ((allLinks[0][iLink].crateID == iCrate) && (allLinks[0][iLink].even == false)) {
0223         foundOdd = true;
0224         for (unsigned int iBX = 0; iBX < nBX; iBX++)
0225           odd[iBX] = allLinks[iBX][iLink];
0226       }
0227 
0228       //if success then create RCTInfoVector and fill output object
0229       if (foundEven && foundOdd) {
0230         nCratesFound++;
0231 
0232         //fill rctInfoVector for all BX read out
0233         for (unsigned int iBX = 0; iBX < nBX; iBX++) {
0234           //RCTInfoFactory rctInfoFactory;
0235           std::vector<RCTInfo> rctInfoData;
0236           rctDataDecoder.decodeLinks(even[iBX].uint, odd[iBX].uint, rctInfoData);
0237           rctDataDecoder.setRCTInfoCrateID(rctInfoData, iCrate);
0238           allCrateRCTInfo[iBX].push_back(rctInfoData.at(0));
0239         }
0240         break;
0241       }
0242     }
0243   }
0244 
0245   if (nCratesFound != 18)
0246     LogError("L1T") << "Warning -- only found " << nCratesFound << " valid crates";
0247 
0248   //start assuming 1 BX readout
0249   int32_t startBX = 0;
0250   if (nBX == 1)
0251     startBX = 0;
0252   else if (nBX == 3)
0253     startBX = -1;
0254   else if (nBX == 5)
0255     startBX = -2;
0256 
0257   //Step 3: Create Collections from RCTInfo Objects
0258   //Notice, iBX used for grabbing data from array, startBX used for storing in Collection
0259   for (uint32_t iBX = 0; iBX < nBX; iBX++, startBX++) {
0260     for (unsigned int iCrate = 0; iCrate < nCratesFound; iCrate++) {
0261       RCTInfo rctInfo = allCrateRCTInfo[iBX].at(iCrate);
0262       //Use Crate ID to identify eta/phi of candidate
0263       for (int j = 0; j < 4; j++) {
0264         L1CaloEmCand em = L1CaloEmCand(rctInfo.neRank[j], rctInfo.neRegn[j], rctInfo.neCard[j], rctInfo.crateID, false);
0265         em.setBx(startBX);
0266         colls->rctEm()->push_back(em);
0267       }
0268 
0269       for (int j = 0; j < 4; j++) {
0270         L1CaloEmCand em = L1CaloEmCand(rctInfo.ieRank[j], rctInfo.ieRegn[j], rctInfo.ieCard[j], rctInfo.crateID, true);
0271         em.setBx(startBX);
0272         colls->rctEm()->push_back(em);
0273       }
0274 
0275       for (int j = 0; j < 7; j++) {
0276         for (int k = 0; k < 2; k++) {
0277           bool o = (((rctInfo.oBits >> (j * 2 + k)) & 0x1) == 0x1);
0278           bool t = (((rctInfo.tBits >> (j * 2 + k)) & 0x1) == 0x1);
0279           bool m = (((rctInfo.mBits >> (j * 2 + k)) & 0x1) == 0x1);
0280           bool q = (((rctInfo.qBits >> (j * 2 + k)) & 0x1) == 0x1);
0281           L1CaloRegion rgn = L1CaloRegion(rctInfo.rgnEt[j][k], o, t, m, q, rctInfo.crateID, j, k);
0282           rgn.setBx(startBX);
0283           colls->rctCalo()->push_back(rgn);
0284         }
0285       }
0286 
0287       for (int k = 0; k < 4; k++) {
0288         for (int j = 0; j < 2; j++) {
0289           // 0 1 4 5 2 3 6 7
0290           uint32_t offset = j * 2 + k % 2 + (k / 2) * 4;
0291           bool fg = (((rctInfo.hfQBits >> offset) & 0x1) == 0x1);
0292           L1CaloRegion rgn = L1CaloRegion(rctInfo.hfEt[j][k], fg, rctInfo.crateID, (j * 4 + k));
0293           rgn.setBx(startBX);
0294           colls->rctCalo()->push_back(rgn);
0295         }
0296       }
0297     }
0298   }
0299 }
0300 
0301 bool RctRawToDigi::printAll(const unsigned char* data, const unsigned size) {
0302   for (unsigned i = 0; i < size; i++) {
0303     std::cout << data[i] << " ";
0304     if (i % 6 == 5)
0305       std::cout << std::endl;
0306   }
0307   return true;
0308 }
0309 
0310 // ------------ method called once each job just after ending the event loop  ------------
0311 void RctRawToDigi::endJob() {
0312   unsigned total = 0;
0313   std::ostringstream os;
0314 
0315   for (unsigned i = 0; i <= MAX_ERR_CODE; ++i) {
0316     total += errorCounters_.at(i);
0317     os << "Error " << i << " (" << errorCounters_.at(i) << ")";
0318     if (i < MAX_ERR_CODE) {
0319       os << ", ";
0320     }
0321   }
0322 
0323   if (total > 0 && verbose_) {
0324     edm::LogError("RCT") << "Encountered " << total << " unpacking errors: " << os.str();
0325   }
0326 }
0327 
0328 /// make this a plugin
0329 DEFINE_FWK_MODULE(RctRawToDigi);