Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 /**
0002  * \file EcalPnGraphDumperModule.h 
0003  * module dumping TGraph with 50 data frames from Pn Diodes
0004  *   
0005  * 
0006  * \author N. Amapane - S. Argiro'
0007  * \author G. Franzoni
0008  * \author J. Haupt
0009  *
0010  */
0011 
0012 #include <FWCore/Framework/interface/one/EDAnalyzer.h>
0013 #include <FWCore/Framework/interface/Event.h>
0014 #include <FWCore/Framework/interface/MakerMacros.h>
0015 #include <DataFormats/EcalDigi/interface/EcalDigiCollections.h>
0016 #include <DataFormats/EcalDetId/interface/EcalDetIdCollections.h>
0017 
0018 #include <DataFormats/EcalDigi/interface/EcalTriggerPrimitiveDigi.h>
0019 #include <DataFormats/EcalDigi/interface/EcalTriggerPrimitiveSample.h>
0020 
0021 #include <iostream>
0022 #include <vector>
0023 
0024 #include "TFile.h"
0025 #include "TGraph.h"
0026 
0027 class EcalPnGraphDumperModule : public edm::one::EDAnalyzer<> {
0028 public:
0029   EcalPnGraphDumperModule(const edm::ParameterSet& ps);
0030   ~EcalPnGraphDumperModule();
0031 
0032   std::string intToString(int num);
0033 
0034   void analyze(const edm::Event& e, const edm::EventSetup& c);
0035 
0036 protected:
0037   int verbosity;
0038   int eventCounter;
0039 
0040   int ieb_id;
0041   int first_Pn;
0042 
0043   bool inputIsOk;
0044 
0045   std::string fileName;
0046 
0047   std::vector<int> listChannels;
0048   std::vector<int> listAllChannels;
0049   std::vector<int> listPns;
0050   std::vector<int> listAllPns;
0051 
0052   int numPn;
0053 
0054   int abscissa[50];
0055   int ordinate[50];
0056 
0057   std::vector<TGraph> graphs;
0058 
0059   TFile* root_file;
0060 };
0061 
0062 EcalPnGraphDumperModule::EcalPnGraphDumperModule(const edm::ParameterSet& ps) {
0063   fileName = ps.getUntrackedParameter<std::string>("fileName", std::string("toto"));
0064 
0065   ieb_id = ps.getUntrackedParameter<int>("ieb_id", 1);
0066   first_Pn = 0;
0067 
0068   listPns = ps.getUntrackedParameter<std::vector<int> >("listPns", std::vector<int>());
0069 
0070   numPn = ps.getUntrackedParameter<int>("numPn", 9);
0071 
0072   // consistency checks checks
0073   inputIsOk = true;
0074 
0075   std::vector<int>::iterator intIter;
0076 
0077   for (intIter = listPns.begin(); intIter != listPns.end(); intIter++) {
0078     if (((*intIter) < 1) || (10 < (*intIter))) {
0079       std::cout << "[EcalPnGraphDumperModule] pn number : " << (*intIter) << " found in listPns. "
0080                 << " Valid range is 1-10. Returning." << std::endl;
0081       inputIsOk = false;
0082       return;
0083     }
0084     if (!first_Pn)
0085       first_Pn = (*intIter);
0086   }
0087 
0088   // setting the abcissa array once for all
0089   for (int i = 0; i < 50; i++)
0090     abscissa[i] = i;
0091 
0092   // local event counter (in general different from LV1)
0093   eventCounter = 0;
0094 }
0095 
0096 EcalPnGraphDumperModule::~EcalPnGraphDumperModule() {
0097   fileName += (std::string("_iEB") + intToString(ieb_id));
0098   fileName += (std::string("_Pn") + intToString(first_Pn));
0099   fileName += ".graph.root";
0100 
0101   root_file = new TFile(fileName.c_str(), "RECREATE");
0102   std::vector<TGraph>::iterator gr_it;
0103   for (gr_it = graphs.begin(); gr_it != graphs.end(); gr_it++)
0104     (*gr_it).Write();
0105   root_file->Close();
0106 }
0107 
0108 std::string EcalPnGraphDumperModule::intToString(int num) {
0109   //
0110   // outputs the number into the string stream and then flushes
0111   // the buffer (makes sure the output is put into the stream)
0112   //
0113   std::ostringstream myStream;  //creates an ostringstream object
0114   myStream << num << std::flush;
0115 
0116   return (myStream.str());  //returns the string form of the stringstream object
0117 }
0118 
0119 void EcalPnGraphDumperModule::analyze(const edm::Event& e, const edm::EventSetup& c) {
0120   eventCounter++;
0121   if (!inputIsOk)
0122     return;
0123 
0124   // retrieving crystal PN diodes from Event
0125   edm::Handle<EcalPnDiodeDigiCollection> PNs;
0126   e.getByLabel("ecalEBunpacker", PNs);
0127 
0128   // getting the list of all the Pns which will be dumped on TGraph
0129   // - listPns is the list as given by the user
0130   // - numPn is number of Pn's for which graphs are required
0131   std::vector<int>::iterator pn_it;
0132   for (pn_it = listPns.begin(); pn_it != listPns.end(); pn_it++) {
0133     int ipn = (*pn_it);
0134     int hpn = (numPn / 2);
0135 
0136     for (int u = (-hpn); u <= hpn; u++) {
0137       int ipn_c = ipn + u;
0138       if (ipn_c < 1 || ipn_c > 10)
0139         continue;
0140       listAllPns.push_back(ipn_c);
0141     }
0142   }
0143 
0144   // loop over all the available PN digis and make graphs for those which have been chosen by the user
0145   for (EcalPnDiodeDigiCollection::const_iterator pnItr = PNs->begin(); pnItr != PNs->end(); ++pnItr) {
0146     {
0147       int ipn = (*pnItr).id().iPnId();
0148       int ieb = EcalPnDiodeDetId((*pnItr).id()).iDCCId();
0149 
0150       // selecting based on DCCId
0151       if (ieb != ieb_id)
0152         return;
0153 
0154       // selecting desired Pns only
0155       std::vector<int>::iterator iPnIter;
0156       iPnIter = find(listAllPns.begin(), listAllPns.end(), ipn);
0157       if (iPnIter == listAllPns.end()) {
0158         continue;
0159       }
0160 
0161       for (int i = 0; i < (*pnItr).size() && i < 50; ++i) {
0162         ordinate[i] = (*pnItr).sample(i).adc();
0163       }
0164 
0165       TGraph oneGraph(50, abscissa, ordinate);
0166       std::string title;
0167       title = "Graph_ev" + intToString(eventCounter) + "_ipn" + intToString(ipn);
0168       oneGraph.SetTitle(title.c_str());
0169       oneGraph.SetName(title.c_str());
0170       graphs.push_back(oneGraph);
0171 
0172     }  // loop in crystals
0173   }
0174 }
0175 
0176 DEFINE_FWK_MODULE(EcalPnGraphDumperModule);