CMS 3D CMS Logo

RawDataConverter.h
Go to the documentation of this file.
7 
8 // Forward declarations
9 class TFile;
10 class TTree;
11 
12 class RawDataConverter : public edm::one::EDAnalyzer<edm::one::WatchRuns> {
13 public:
14  explicit RawDataConverter(const edm::ParameterSet&);
15  ~RawDataConverter() override = default;
16 
17 private:
19  void beginJob() override;
20  void beginRun(edm::Run const&, edm::EventSetup const&) override;
21  void analyze(const edm::Event&, const edm::EventSetup&) override;
22  void endRun(edm::Run const&, edm::EventSetup const&) override {}
23  void endJob() override;
24 
25  void fillDetectorId(void);
26  void ClearData(void);
28  const edm::Event&
29  iEvent); // Check what kind of file is being processed and get valid module and instance labels returns the type of Digis that was found
30 
31  template <class T>
32  void GetDigis(const edm::Event&); // Copy the Digis into the local container (theData)
33 
34  std::vector<std::string> theDigiModuleLabels;
35  std::vector<std::string> theProductInstanceLabels;
36 
39 
40  TFile* theOutputFile;
41  TTree* theOutputTree;
43 
44  int latency;
46  int runnumber;
47  int lumiBlock;
49 };
50 
51 // Copy the Digis into the local container (theData)
52 // Currently this has only been implemented and tested for SiStripDigis
53 // SiStripRawDigis and SiStripProcessedRawDigis will need some changes to work (this is the final goal)
54 template <class Digitype>
56  LogDebug("RawDataConverter") << "Fill ZeroSuppressed Digis into the Tree";
57 
58  // Get the DetSetVector for the SiStripDigis
59  // This is a vector with all the modules, each module containing zero or more strips with signal (Digis)
60  edm::Handle<edm::DetSetVector<Digitype> > detSetVector; // Handle for holding the DetSetVector
61  iEvent.getByLabel(CurrentModuleLabel, CurrentInstanceLabel, detSetVector);
62  if (!detSetVector.isValid())
63  throw std::runtime_error("Could not find the Digis");
64 
65  // set everything in the local container to zero
66  ClearData();
67 
68  // Fill the Digis into the Raw Data Container
69 
70  LASGlobalLoop loop; // loop helper
71  int det, ring, beam, disk, pos; // and its variables
72 
73  // loop over TEC+- (internal) modules
74  det = 0;
75  ring = 0;
76  beam = 0;
77  disk = 0;
78  do {
79  // Find the module in the DetSetVector and get a pointer (iterator) to it
81  detSetVector->find(detectorId.GetTECEntry(det, ring, beam, disk));
82 
83  if (theModule != detSetVector->end()) {
84  // loop over all the Digis in this Module
86  for (theDigi = theModule->data.begin(); theDigi != theModule->data.end(); ++theDigi) {
87  // fill the number of adc counts into the local container
88  if (theDigi->channel() < 512)
89  theData.GetTECEntry(det, ring, beam, disk).at(theDigi->channel()) = theDigi->adc();
90  }
91  }
92  } while (loop.TECLoop(det, ring, beam, disk));
93 
94  // loop TIB/TOB
95  det = 2;
96  beam = 0;
97  pos = 0; // <- set det = 2 (TIB)
98  do {
99  // Find the module in the DetSetVector and get a pointer (iterator) to it
101  detSetVector->find(detectorId.GetTIBTOBEntry(det, beam, pos));
102 
103  if (theModule != detSetVector->end()) {
104  // loop over all the Digis in this Module
105  typename edm::DetSet<Digitype>::const_iterator theDigi;
106  for (theDigi = theModule->data.begin(); theDigi != theModule->data.end(); ++theDigi) {
107  // fill the number of adc counts into the local container
108  if (theDigi->channel() < 512)
109  theData.GetTIBTOBEntry(det, beam, pos).at(theDigi->channel()) = theDigi->adc();
110  }
111  }
112  } while (loop.TIBTOBLoop(det, beam, pos));
113 
114  // loop TEC (AT)
115  det = 0;
116  beam = 0;
117  disk = 0;
118  do {
119  // Find the module in the DetSetVector and get a pointer (iterator) to it
121  detSetVector->find(detectorId.GetTEC2TECEntry(det, beam, disk));
122 
123  if (theModule != detSetVector->end()) {
124  // loop over all the Digis in this Module
125  typename edm::DetSet<Digitype>::const_iterator theDigi;
126  for (theDigi = theModule->data.begin(); theDigi != theModule->data.end(); ++theDigi) {
127  // fill the number of adc counts into the local container
128  if (theDigi->channel() < 512)
129  theData.GetTEC2TECEntry(det, beam, disk).at(theDigi->channel()) = theDigi->adc();
130  }
131  }
132  } while (loop.TEC2TECLoop(det, beam, disk));
133 }
void fillDetectorId(void)
DigiType GetValidLabels(const edm::Event &iEvent)
~RawDataConverter() override=default
LASGlobalData< int > detectorId
void endRun(edm::Run const &, edm::EventSetup const &) override
void beginJob() override
void analyze(const edm::Event &, const edm::EventSetup &) override
int iEvent
Definition: GenABIO.cc:224
std::vector< std::string > theDigiModuleLabels
void beginRun(edm::Run const &, edm::EventSetup const &) override
void GetDigis(const edm::Event &)
T & GetTIBTOBEntry(int subdetector, int beam, int tibTobPosition)
iterator end()
Return the off-the-end iterator.
Definition: DetSetVector.h:316
RawDataConverter(const edm::ParameterSet &)
T & GetTEC2TECEntry(int subdetector, int beam, int tecDisk)
void endJob() override
T & GetTECEntry(int subdetector, int tecRing, int beam, int tecDisk)
Definition: LASGlobalData.h:84
bool isValid() const
Definition: HandleBase.h:70
std::string CurrentModuleLabel
LASGlobalData< std::vector< float > > theData
iterator begin()
Return an iterator to the first DetSet.
Definition: DetSetVector.h:305
std::string CurrentInstanceLabel
collection_type::const_iterator const_iterator
Definition: DetSet.h:31
collection_type::const_iterator const_iterator
Definition: DetSetVector.h:102
std::vector< std::string > theProductInstanceLabels
Definition: Run.h:45
#define LogDebug(id)