File indexing completed on 2023-03-17 11:00:00
0001 #include "PhysicsToBitConverter.h"
0002 #include "rctDataBase.h"
0003
0004 #include "DataFormats/L1CaloTrigger/interface/L1CaloEmCand.h"
0005 #include "DataFormats/L1TCalorimeter/interface/CaloEmCand.h"
0006 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegion.h"
0007 #include "DataFormats/L1TCalorimeter/interface/CaloRegion.h"
0008 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegionDetId.h"
0009 #include "DataFormats/L1CaloTrigger/interface/L1CaloCollections.h"
0010
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012
0013 #include "EventFilter/L1TRawToDigi/plugins/UnpackerFactory.h"
0014
0015 #include <iostream>
0016 #include <fstream>
0017
0018 #include "EventFilter/L1TRawToDigi/interface/Block.h"
0019
0020 #include "CaloCollections.h"
0021 #include "RCTEmRegionUnpacker.h"
0022
0023 namespace l1t {
0024 namespace stage1 {
0025 void unpack_em(const Block& block, UnpackerCollections* coll) {
0026 int nBX, firstBX, lastBX;
0027 nBX = int(ceil(block.header().getSize() / 6.));
0028 getBXRange(nBX, firstBX, lastBX);
0029
0030 auto resRCTEMCands_ = static_cast<CaloCollections*>(coll)->getCaloEmCands();
0031
0032
0033 int unsigned i = 0;
0034
0035 for (int bx = firstBX; bx <= lastBX; bx++) {
0036 unsigned int crate;
0037 bool even = false;
0038
0039 std::vector<uint32_t> uint;
0040 uint.reserve(6);
0041
0042 PhysicsToBitConverter converter;
0043 rctDataBase database;
0044 int mp7link = (int)(block.header().getID() / 2);
0045 database.GetLinkRCT(mp7link, crate, even);
0046
0047 uint.push_back(block.payload()[i++]);
0048 uint.push_back(block.payload()[i++]);
0049 uint.push_back(block.payload()[i++]);
0050 uint.push_back(block.payload()[i++]);
0051 uint.push_back(block.payload()[i++]);
0052 uint.push_back(block.payload()[i++]);
0053
0054 LogDebug("L1T") << "--------------- mp7 link =" << mp7link << "RCT crate id=" << crate
0055 << ", RCT crate even=" << even << std::endl;
0056
0057 if (!even) {
0058 for (int i = 0; i < 6; i++)
0059 converter.Set32bitWordLinkOdd(i, uint[i]);
0060 converter.Convert();
0061
0062 for (int j = 0; j < 4; j++) {
0063 unsigned int rank = (unsigned int)converter.GetNEEt(j);
0064 unsigned int reg = (unsigned int)converter.GetNEReg(j);
0065 unsigned int card = (unsigned int)converter.GetNECard(j);
0066
0067 LogDebug("L1T") << "UNPACKER, CRATE" << crate << "NON ISO em rank=" << rank << ", region=" << reg
0068 << ", card=" << card << std::endl;
0069
0070 L1CaloEmCand em = L1CaloEmCand(rank, reg, card, crate, false, j, bx);
0071 resRCTEMCands_->push_back(em);
0072 }
0073
0074 for (int j = 0; j < 4; j++) {
0075 unsigned int rank = converter.GetIEEt(j);
0076 unsigned int reg = converter.GetIEReg(j);
0077 unsigned int card = converter.GetIECard(j);
0078
0079 LogDebug("L1T") << "UNPACKER, CRATE" << crate << "ISO em rank=" << rank << ", region=" << reg
0080 << ", card=" << card << std::endl;
0081 L1CaloEmCand em = L1CaloEmCand(rank, reg, card, crate, true, j, bx);
0082 resRCTEMCands_->push_back(em);
0083 }
0084 }
0085 }
0086 }
0087
0088 void unpack_region(const Block& block, UnpackerCollections* coll) {
0089 int nBX, firstBX, lastBX;
0090 nBX = int(ceil(block.header().getSize() / 6.));
0091 getBXRange(nBX, firstBX, lastBX);
0092
0093 auto resRCTRegions_ = static_cast<CaloCollections*>(coll)->getCaloRegions();
0094
0095
0096 int unsigned i = 0;
0097 std::vector<uint32_t> uint;
0098 uint.reserve(6);
0099
0100 for (int bx = firstBX; bx <= lastBX; bx++) {
0101 unsigned int crate;
0102 bool even = false;
0103
0104 PhysicsToBitConverter converter;
0105 rctDataBase database;
0106 int mp7link = (int)(block.header().getID() / 2);
0107 database.GetLinkRCT(mp7link, crate, even);
0108
0109 uint.push_back(block.payload()[i++]);
0110 uint.push_back(block.payload()[i++]);
0111 uint.push_back(block.payload()[i++]);
0112 uint.push_back(block.payload()[i++]);
0113 uint.push_back(block.payload()[i++]);
0114 uint.push_back(block.payload()[i++]);
0115
0116 LogDebug("L1T") << "--------------- mp7 link =" << mp7link << "RCT crate id=" << crate
0117 << ", RCT crate even=" << even << std::endl;
0118
0119 if (!even) {
0120 for (int i = 0; i < 6; i++)
0121 converter.Set32bitWordLinkOdd(i, uint[i]);
0122 converter.Convert();
0123
0124 for (int j = 0; j < 8; j++) {
0125 unsigned int hfet = (unsigned int)converter.GetHFEt(j);
0126 bool hfgrain = (bool)converter.GetHFFg(j);
0127
0128 LogDebug("L1T") << "UNPACKER, CRATE" << crate << "region=" << j << ", rgnEt=" << hfet << std::endl;
0129 L1CaloRegion rgn = L1CaloRegion(hfet, hfgrain, crate, j);
0130 rgn.setBx(bx);
0131 resRCTRegions_->push_back(rgn);
0132 }
0133 }
0134
0135 else {
0136 for (int i = 0; i < 6; i++)
0137 converter.Set32bitWordLinkEven(i, uint[i]);
0138 converter.Convert();
0139
0140 for (int j = 0; j < 7; j++) {
0141 for (int k = 0; k < 2; k++) {
0142 unsigned int RCet = (unsigned int)converter.GetRCEt(j, k);
0143 bool overflow = (bool)converter.GetRCOf(j, k);
0144 bool tauveto = (bool)converter.GetRCTau(j, k);
0145 bool hadveto = (bool)converter.GetRCHad(j, k);
0146 bool quiet = false;
0147
0148 LogDebug("L1T") << "UNPACKER, CRATE=" << crate << ",region=" << k << ", card=" << j << ", rgnEt=" << RCet
0149 << ", overflow=" << overflow << ", tauveto=" << tauveto << ", hadveto=" << hadveto
0150 << std::endl;
0151 L1CaloRegion rgn = L1CaloRegion(RCet, overflow, tauveto, hadveto, quiet, crate, j, k);
0152 rgn.setBx(bx);
0153 resRCTRegions_->push_back(rgn);
0154 }
0155 }
0156 }
0157 }
0158 }
0159
0160 bool RCTEmRegionUnpacker::unpack(const Block& block, UnpackerCollections* coll) {
0161 if (block.header().getCapID() == 0) {
0162 unpack_region(block, coll);
0163 } else if (block.header().getCapID() == 1) {
0164 unpack_em(block, coll);
0165 } else {
0166 return false;
0167 }
0168 return true;
0169 }
0170 }
0171 }
0172
0173 DEFINE_L1T_UNPACKER(l1t::stage1::RCTEmRegionUnpacker);