Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:29:46

0001 #include <memory>
0002 #include <iostream>
0003 #include <vector>
0004 
0005 #include "FWCore/Framework/interface/Frameworkfwd.h"
0006 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Framework/interface/MakerMacros.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 
0011 #include "DataFormats/FTLDigi/interface/FTLDigiCollections.h"
0012 
0013 class MTDDigiDump : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0014 public:
0015   explicit MTDDigiDump(const edm::ParameterSet&);
0016   ~MTDDigiDump() override;
0017 
0018   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0019 
0020 private:
0021   void beginJob() override;
0022   void analyze(const edm::Event&, const edm::EventSetup&) override;
0023   void endJob() override;
0024 
0025   // ----------member data ---------------------------
0026 
0027   edm::EDGetTokenT<BTLDigiCollection> tok_BTL_digi;
0028   edm::EDGetTokenT<ETLDigiCollection> tok_ETL_digi;
0029 };
0030 
0031 MTDDigiDump::MTDDigiDump(const edm::ParameterSet& iConfig)
0032 
0033 {
0034   tok_BTL_digi = consumes<BTLDigiCollection>(edm::InputTag("mix", "FTLBarrel"));
0035   tok_ETL_digi = consumes<ETLDigiCollection>(edm::InputTag("mix", "FTLEndcap"));
0036 }
0037 
0038 MTDDigiDump::~MTDDigiDump() {}
0039 
0040 //
0041 // member functions
0042 //
0043 
0044 // ------------ method called for each event ------------
0045 void MTDDigiDump::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0046   using namespace std;
0047 
0048   edm::Handle<BTLDigiCollection> h_BTL_digi;
0049   iEvent.getByToken(tok_BTL_digi, h_BTL_digi);
0050 
0051   edm::Handle<ETLDigiCollection> h_ETL_digi;
0052   iEvent.getByToken(tok_ETL_digi, h_ETL_digi);
0053 
0054   // --- BTL DIGIs:
0055 
0056   if (!h_BTL_digi->empty()) {
0057     std::cout << " ----------------------------------------" << std::endl;
0058     std::cout << " BTL DIGI collection: \n" << std::endl;
0059 
0060     for (const auto& dataFrame : *h_BTL_digi) {
0061       // --- detector element ID:
0062       std::cout << "\n BTL DIGI:  row = " << dataFrame.row() << " col = " << dataFrame.column()
0063                 << " BTLDetId = " << dataFrame.id();
0064 
0065       // --- loop over the dataFrame samples
0066       for (int isample = 0; isample < dataFrame.size(); ++isample) {
0067         const auto& sample = dataFrame.sample(isample);
0068 
0069         std::cout << "       sample " << isample << ":";
0070         if (sample.data() == 0 && sample.toa() == 0) {
0071           std::cout << std::endl;
0072           continue;
0073         }
0074         std::cout << "  amplitude = " << sample.data() << "  time1 = " << sample.toa() << "  time2 = " << sample.toa2()
0075                   << " r/c = " << (uint32_t)sample.row() << " / " << (uint32_t)sample.column()
0076                   << " th = " << sample.threshold() << " mode = " << sample.mode() << std::endl;
0077 
0078       }  // isaple loop
0079 
0080     }  // digi loop
0081 
0082   }  // if ( h_BTL_digi->size() > 0 )
0083 
0084   // --- ETL DIGIs:
0085 
0086   if (!h_ETL_digi->empty()) {
0087     std::cout << "\n ----------------------------------------" << std::endl;
0088     std::cout << " ETL DIGI collection: \n" << std::endl;
0089 
0090     for (const auto& dataFrame : *h_ETL_digi) {
0091       // --- detector element ID:
0092       std::cout << "\n ETL DIGI:  row = " << dataFrame.row() << " col = " << dataFrame.column()
0093                 << " ETLDetId = " << dataFrame.id();
0094 
0095       // --- loop over the dataFrame samples
0096       for (int isample = 0; isample < dataFrame.size(); ++isample) {
0097         const auto& sample = dataFrame.sample(isample);
0098 
0099         std::cout << "       sample " << isample << ":";
0100         if (sample.data() == 0 && sample.toa() == 0) {
0101           std::cout << std::endl;
0102           continue;
0103         }
0104         std::cout << "  amplitude = " << sample.data() << "  time = " << sample.toa() << " r/c = " << sample.row()
0105                   << " / " << sample.column() << " th = " << sample.threshold() << " mode = " << sample.mode()
0106                   << std::endl;
0107 
0108       }  // isample loop
0109 
0110     }  // digi loop
0111 
0112   }  // if ( h_ETL_digi->size() > 0 )
0113 }
0114 
0115 // ------------ method called once each job just before starting event loop  ------------
0116 void MTDDigiDump::beginJob() {}
0117 
0118 // ------------ method called once each job just after ending the event loop  ------------
0119 void MTDDigiDump::endJob() {}
0120 
0121 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0122 void MTDDigiDump::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0123   //The following says we do not know what parameters are allowed so do no validation
0124   // Please change this to state exactly what you do use, even if it is no parameters
0125   edm::ParameterSetDescription desc;
0126   desc.setUnknown();
0127   descriptions.addDefault(desc);
0128 }
0129 
0130 //define this as a plug-in
0131 DEFINE_FWK_MODULE(MTDDigiDump);