CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HcalRawToDigi.cc
Go to the documentation of this file.
1 using namespace std;
9 #include <iostream>
10 
12  unpacker_(conf.getUntrackedParameter<int>("HcalFirstFED",int(FEDNumbering::MINHCALFEDID)),conf.getParameter<int>("firstSample"),conf.getParameter<int>("lastSample")),
13  filter_(conf.getParameter<bool>("FilterDataQuality"),conf.getParameter<bool>("FilterDataQuality"),
14  false,
15  0, 0,
16  -1),
17  fedUnpackList_(conf.getUntrackedParameter<std::vector<int> >("FEDs", std::vector<int>())),
18  firstFED_(conf.getUntrackedParameter<int>("HcalFirstFED",FEDNumbering::MINHCALFEDID)),
19  unpackCalib_(conf.getUntrackedParameter<bool>("UnpackCalib",false)),
20  unpackZDC_(conf.getUntrackedParameter<bool>("UnpackZDC",false)),
21  unpackTTP_(conf.getUntrackedParameter<bool>("UnpackTTP",false)),
22  silent_(conf.getUntrackedParameter<bool>("silent",true)),
23  complainEmptyData_(conf.getUntrackedParameter<bool>("ComplainEmptyData",false)),
24  unpackerMode_(conf.getUntrackedParameter<int>("UnpackerMode",0)),
25  expectedOrbitMessageTime_(conf.getUntrackedParameter<int>("ExpectedOrbitMessageTime",-1))
26 {
27  tok_data_ = consumes<FEDRawDataCollection>(conf.getParameter<edm::InputTag>("InputLabel"));
28 
29  if (fedUnpackList_.empty()) {
31  fedUnpackList_.push_back(i);
32  }
33 
36  std::ostringstream ss;
37  for (unsigned int i=0; i<fedUnpackList_.size(); i++)
38  ss << fedUnpackList_[i] << " ";
39  edm::LogInfo("HCAL") << "HcalRawToDigi will unpack FEDs ( " << ss.str() << ")";
40 
41  // products produced...
42  produces<HBHEDigiCollection>();
43  produces<HFDigiCollection>();
44  produces<HODigiCollection>();
45  produces<HcalTrigPrimDigiCollection>();
46  produces<HOTrigPrimDigiCollection>();
47  produces<HcalUnpackerReport>();
48  if (unpackCalib_)
49  produces<HcalCalibDigiCollection>();
50  if (unpackZDC_)
51  produces<ZDCDigiCollection>();
52  if (unpackTTP_)
53  produces<HcalTTPDigiCollection>();
54 
55  memset(&stats_,0,sizeof(stats_));
56 
57 }
58 
59 // Virtual destructor needed.
61 
62 // Functions that gets called by framework every event
64 {
65  // Step A: Get Inputs
67  e.getByToken(tok_data_,rawraw);
68  // get the mapping
70  es.get<HcalDbRecord>().get( pSetup );
71  const HcalElectronicsMap* readoutMap=pSetup->getHcalMapping();
72 
73  // Step B: Create empty output : three vectors for three classes...
74  std::vector<HBHEDataFrame> hbhe;
75  std::vector<HODataFrame> ho;
76  std::vector<HFDataFrame> hf;
77  std::vector<HcalTriggerPrimitiveDigi> htp;
78  std::vector<HcalCalibDataFrame> hc;
79  std::vector<ZDCDataFrame> zdc;
80  std::vector<HcalTTPDigi> ttp;
81  std::vector<HOTriggerPrimitiveDigi> hotp;
82  std::auto_ptr<HcalUnpackerReport> report(new HcalUnpackerReport);
83 
84  // Heuristics: use ave+(max-ave)/8
85  if (stats_.max_hbhe>0) hbhe.reserve(stats_.ave_hbhe+(stats_.max_hbhe-stats_.ave_hbhe)/8);
86  if (stats_.max_ho>0) ho.reserve(stats_.ave_ho+(stats_.max_ho-stats_.ave_ho)/8);
87  if (stats_.max_hf>0) hf.reserve(stats_.ave_hf+(stats_.max_hf-stats_.ave_hf)/8);
89  if (stats_.max_tp>0) htp.reserve(stats_.ave_tp+(stats_.max_tp-stats_.ave_tp)/8);
90  if (stats_.max_tpho>0) hotp.reserve(stats_.ave_tpho+(stats_.max_tpho-stats_.ave_tpho)/8);
91 
92  if (unpackZDC_) zdc.reserve(24);
93 
94 
96  colls.hbheCont=&hbhe;
97  colls.hoCont=&ho;
98  colls.hfCont=&hf;
99  colls.tpCont=&htp;
100  colls.tphoCont=&hotp;
101  colls.calibCont=&hc;
102  colls.zdcCont=&zdc;
103  if (unpackTTP_) colls.ttp=&ttp;
104 
105  // Step C: unpack all requested FEDs
106  for (std::vector<int>::const_iterator i=fedUnpackList_.begin(); i!=fedUnpackList_.end(); i++) {
107  const FEDRawData& fed = rawraw->FEDData(*i);
108  if (fed.size()==0) {
109  if (complainEmptyData_) {
110  if (!silent_) edm::LogWarning("EmptyData") << "No data for FED " << *i;
111  report->addError(*i);
112  }
113  } else if (fed.size()<8*3) {
114  if (!silent_) edm::LogWarning("EmptyData") << "Tiny data " << fed.size() << " for FED " << *i;
115  report->addError(*i);
116  } else {
117  try {
118  unpacker_.unpack(fed,*readoutMap,colls, *report,silent_);
119  report->addUnpacked(*i);
120  } catch (cms::Exception& e) {
121  if (!silent_) edm::LogWarning("Unpacking error") << e.what();
122  report->addError(*i);
123  } catch (...) {
124  if (!silent_) edm::LogWarning("Unpacking exception");
125  report->addError(*i);
126  }
127  }
128  }
129 
130 
131  // gather statistics
132  stats_.max_hbhe=std::max(stats_.max_hbhe,(int)hbhe.size());
133  stats_.ave_hbhe=(stats_.ave_hbhe*stats_.n+hbhe.size())/(stats_.n+1);
134  stats_.max_ho=std::max(stats_.max_ho,(int)ho.size());
135  stats_.ave_ho=(stats_.ave_ho*stats_.n+ho.size())/(stats_.n+1);
136  stats_.max_hf=std::max(stats_.max_hf,(int)hf.size());
137  stats_.ave_hf=(stats_.ave_hf*stats_.n+hf.size())/(stats_.n+1);
138  stats_.max_tp=std::max(stats_.max_tp,(int)htp.size());
139  stats_.ave_tp=(stats_.ave_tp*stats_.n+htp.size())/(stats_.n+1);
140  stats_.max_tpho=std::max(stats_.max_tpho,(int)hotp.size());
141  stats_.ave_tpho=(stats_.ave_tpho*stats_.n+hotp.size())/(stats_.n+1);
142  stats_.max_calib=std::max(stats_.max_calib,(int)hc.size());
143  stats_.ave_calib=(stats_.ave_calib*stats_.n+hc.size())/(stats_.n+1);
144 
145 
146  stats_.n++;
147 
148  // Step B: encapsulate vectors in actual collections
149  std::auto_ptr<HBHEDigiCollection> hbhe_prod(new HBHEDigiCollection());
150  std::auto_ptr<HFDigiCollection> hf_prod(new HFDigiCollection());
151  std::auto_ptr<HODigiCollection> ho_prod(new HODigiCollection());
152  std::auto_ptr<HcalTrigPrimDigiCollection> htp_prod(new HcalTrigPrimDigiCollection());
153  std::auto_ptr<HOTrigPrimDigiCollection> hotp_prod(new HOTrigPrimDigiCollection());
154 
155  hbhe_prod->swap_contents(hbhe);
156  hf_prod->swap_contents(hf);
157  ho_prod->swap_contents(ho);
158  htp_prod->swap_contents(htp);
159  hotp_prod->swap_contents(hotp);
160 
161  // Step C2: filter FEDs, if required
162  if (filter_.active()) {
163  HBHEDigiCollection filtered_hbhe=filter_.filter(*hbhe_prod,*report);
164  HODigiCollection filtered_ho=filter_.filter(*ho_prod,*report);
165  HFDigiCollection filtered_hf=filter_.filter(*hf_prod,*report);
166 
167  hbhe_prod->swap(filtered_hbhe);
168  ho_prod->swap(filtered_ho);
169  hf_prod->swap(filtered_hf);
170  }
171 
172 
173  // Step D: Put outputs into event
174  // just until the sorting is proven
175  hbhe_prod->sort();
176  ho_prod->sort();
177  hf_prod->sort();
178  htp_prod->sort();
179  hotp_prod->sort();
180 
181  e.put(hbhe_prod);
182  e.put(ho_prod);
183  e.put(hf_prod);
184  e.put(htp_prod);
185  e.put(hotp_prod);
186 
188  if (unpackCalib_) {
189  std::auto_ptr<HcalCalibDigiCollection> hc_prod(new HcalCalibDigiCollection());
190  hc_prod->swap_contents(hc);
191 
192  if (filter_.active()) {
193  HcalCalibDigiCollection filtered_calib=filter_.filter(*hc_prod,*report);
194  hc_prod->swap(filtered_calib);
195  }
196 
197  hc_prod->sort();
198  e.put(hc_prod);
199  }
200 
202  if (unpackZDC_) {
203  std::auto_ptr<ZDCDigiCollection> prod(new ZDCDigiCollection());
204  prod->swap_contents(zdc);
205 
206  if (filter_.active()) {
207  ZDCDigiCollection filtered_zdc=filter_.filter(*prod,*report);
208  prod->swap(filtered_zdc);
209  }
210 
211  prod->sort();
212  e.put(prod);
213  }
214 
215  if (unpackTTP_) {
216  std::auto_ptr<HcalTTPDigiCollection> prod(new HcalTTPDigiCollection());
217  prod->swap_contents(ttp);
218 
219  prod->sort();
220  e.put(prod);
221  }
222  e.put(report);
223 
224 
225 }
226 
227 
virtual char const * what() const
Definition: Exception.cc:141
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
HcalDataFrameFilter filter_
Definition: HcalRawToDigi.h:38
std::vector< HcalTTPDigi > * ttp
Definition: HcalUnpacker.h:31
edm::SortedCollection< HcalTriggerPrimitiveDigi > HcalTrigPrimDigiCollection
std::vector< HFDataFrame > * hfCont
Definition: HcalUnpacker.h:26
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:434
void swap(SortedCollection &other)
virtual void produce(edm::Event &e, const edm::EventSetup &c)
std::vector< HBHEDataFrame > * hbheCont
Definition: HcalUnpacker.h:24
std::vector< HOTriggerPrimitiveDigi > * tphoCont
Definition: HcalUnpacker.h:30
edm::SortedCollection< HOTriggerPrimitiveDigi > HOTrigPrimDigiCollection
size_t size() const
Lenght of the data buffer in bytes.
Definition: FEDRawData.h:47
edm::SortedCollection< ZDCDataFrame > ZDCDigiCollection
tuple report
Definition: zeeHLT_cff.py:9
edm::SortedCollection< HODataFrame > HODigiCollection
bool complainEmptyData_
Definition: HcalRawToDigi.h:42
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:116
virtual ~HcalRawToDigi()
std::vector< HcalTriggerPrimitiveDigi > * tpCont
Definition: HcalUnpacker.h:29
struct HcalRawToDigi::Statistics stats_
int expectedOrbitMessageTime_
Definition: HcalRawToDigi.h:43
HcalRawToDigi(const edm::ParameterSet &ps)
tuple conf
Definition: dbtoconf.py:185
void setExpectedOrbitMessageTime(int time)
Definition: HcalUnpacker.h:38
void setMode(int mode)
Definition: HcalUnpacker.h:47
void unpack(const FEDRawData &raw, const HcalElectronicsMap &emap, std::vector< HcalHistogramDigi > &histoDigis)
edm::SortedCollection< HcalCalibDataFrame > HcalCalibDigiCollection
std::vector< HcalCalibDataFrame > * calibCont
Definition: HcalUnpacker.h:27
std::vector< int > fedUnpackList_
Definition: HcalRawToDigi.h:39
const T & get() const
Definition: EventSetup.h:55
bool active() const
whether any filters are on
std::vector< HODataFrame > * hoCont
Definition: HcalUnpacker.h:25
HBHEDigiCollection filter(const HBHEDigiCollection &incol, HcalUnpackerReport &r)
filter HB/HE data frames
HcalUnpacker unpacker_
Definition: HcalRawToDigi.h:37
edm::SortedCollection< HcalTTPDigi > HcalTTPDigiCollection
susybsm::HSCParticleCollection hc
Definition: classes.h:25
edm::SortedCollection< HFDataFrame > HFDigiCollection
std::vector< ZDCDataFrame > * zdcCont
Definition: HcalUnpacker.h:28
edm::EDGetTokenT< FEDRawDataCollection > tok_data_
Definition: HcalRawToDigi.h:36
volatile std::atomic< bool > shutdown_flag false
edm::SortedCollection< HBHEDataFrame > HBHEDigiCollection