Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-05 22:26:05

0001 /*
0002  * adapting to CTPPS pixel detector March 2017 - F.Ferro
0003  */
0004 
0005 #include "EventFilter/CTPPSRawToDigi/interface/CTPPSPixelRawToDigi.h"
0006 
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 
0011 #include "DataFormats/CTPPSDigi/interface/CTPPSPixelDigi.h"
0012 
0013 #include "DataFormats/Common/interface/DetSetVector.h"
0014 #include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h"
0015 #include "DataFormats/FEDRawData/interface/FEDRawData.h"
0016 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0017 
0018 #include "EventFilter/CTPPSRawToDigi/interface/CTPPSPixelDataFormatter.h"
0019 
0020 #include "FWCore/Framework/interface/ConsumesCollector.h"
0021 
0022 using namespace std;
0023 
0024 CTPPSPixelRawToDigi::CTPPSPixelRawToDigi(const edm::ParameterSet& conf)
0025     : config_(conf)
0026 
0027 {
0028   FEDRawDataCollection_ = consumes<FEDRawDataCollection>(config_.getParameter<edm::InputTag>("inputLabel"));
0029   CTPPSPixelDAQMapping_ = esConsumes<CTPPSPixelDAQMapping, CTPPSPixelDAQMappingRcd>();
0030 
0031   produces<edm::DetSetVector<CTPPSPixelDigi>>();
0032 
0033   isRun3_ = config_.getParameter<bool>("isRun3");
0034   includeErrors_ = config_.getParameter<bool>("includeErrors");
0035   mappingLabel_ = config_.getParameter<std::string>("mappingLabel");
0036 
0037   if (includeErrors_) {
0038     produces<edm::DetSetVector<CTPPSPixelDataError>>();
0039   }
0040 }
0041 
0042 CTPPSPixelRawToDigi::~CTPPSPixelRawToDigi() {
0043   edm::LogInfo("CTPPSPixelRawToDigi") << " CTPPSPixelRawToDigi destructor!";
0044 }
0045 
0046 void CTPPSPixelRawToDigi::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0047   edm::ParameterSetDescription desc;
0048   desc.add<bool>("isRun3", true);
0049   desc.add<bool>("includeErrors", true);
0050   desc.add<edm::InputTag>("inputLabel", edm::InputTag("rawDataCollector"));
0051   desc.add<std::string>("mappingLabel", "RPix");
0052   descriptions.add("ctppsPixelDigis", desc);
0053 }
0054 
0055 void CTPPSPixelRawToDigi::produce(edm::Event& ev, const edm::EventSetup& es) {
0056   edm::Handle<FEDRawDataCollection> buffers;
0057   ev.getByToken(FEDRawDataCollection_, buffers);
0058 
0059   edm::ESHandle<CTPPSPixelDAQMapping> mapping;
0060 
0061   bool data_exist = false;
0062   for (int fed = FEDNumbering::MINCTPPSPixelsFEDID; fed <= FEDNumbering::MAXCTPPSPixelsFEDID; fed++) {
0063     const FEDRawData& tempRawData = buffers->FEDData(fed);
0064     if (tempRawData.size() != 0) {
0065       data_exist = true;
0066       break;
0067     }
0068   }
0069   /// create product (digis & errors)
0070   auto collection = std::make_unique<edm::DetSetVector<CTPPSPixelDigi>>();
0071 
0072   auto errorcollection = std::make_unique<edm::DetSetVector<CTPPSPixelDataError>>();
0073 
0074   if (data_exist) {
0075     mapping = es.getHandle(CTPPSPixelDAQMapping_);
0076 
0077     fedIds_ = mapping->fedIds();
0078 
0079     CTPPSPixelDataFormatter formatter(mapping->ROCMapping);
0080     formatter.setErrorStatus(includeErrors_);
0081 
0082     bool errorsInEvent = false;
0083     CTPPSPixelDataFormatter::DetErrors nodeterrors;
0084 
0085     for (auto aFed = fedIds_.begin(); aFed != fedIds_.end(); ++aFed) {
0086       int fedId = *aFed;
0087 
0088       edm::LogInfo("CTPPSPixelRawToDigi") << " PRODUCE DIGI FOR FED: " << dec << fedId << endl;
0089 
0090       CTPPSPixelDataFormatter::Errors errors;
0091       /// get event data for this fed
0092       const FEDRawData& fedRawData = buffers->FEDData(fedId);
0093 
0094       formatter.interpretRawData(isRun3_, errorsInEvent, fedId, fedRawData, *collection, errors);
0095 
0096       if (includeErrors_) {
0097         for (auto const& is : errors) {
0098           uint32_t errordetid = is.first;
0099           /// errors given dummy detId must be sorted by Fed
0100           if (errordetid == RPixErrorChecker::dummyDetId) {
0101             nodeterrors.insert(nodeterrors.end(), errors[errordetid].begin(), errors[errordetid].end());
0102           } else {
0103             edm::DetSet<CTPPSPixelDataError>& errorDetSet = errorcollection->find_or_insert(errordetid);
0104             errorDetSet.data.insert(errorDetSet.data.end(), is.second.begin(), is.second.end());
0105           }
0106         }
0107       }
0108     }
0109 
0110     if (includeErrors_) {
0111       errorcollection->find_or_insert(RPixErrorChecker::dummyDetId).data = nodeterrors;
0112     }
0113     if (errorsInEvent)
0114       LogDebug("CTPPSPixelRawToDigi") << "Error words were stored in this event";
0115   }
0116   ///send digis and errors back to framework
0117   ev.put(std::move(collection));
0118 
0119   if (includeErrors_) {
0120     ev.put(std::move(errorcollection));
0121   }
0122 }
0123 
0124 DEFINE_FWK_MODULE(CTPPSPixelRawToDigi);