Bug Summary

File:plugins/monitoring/DAQ_online/JEventProcessor_DAQ_online.cc
Location:line 495, column 3
Description:Value stored to 'evio_buffsize' is never read

Annotated Source Code

1// $Id$
2//
3// File: JEventProcessor_DAQ_online.cc
4// Created: Thu Aug 7 09:37:03 EDT 2014
5// Creator: dalton (on Linux gluon05.jlab.org 2.6.32-358.18.1.el6.x86_64 x86_64)
6//
7
8#include <stdint.h>
9#include <vector>
10
11#include "JEventProcessor_DAQ_online.h"
12#include <JANA/JApplication.h>
13#include <JANA/JFactory.h>
14
15using namespace std;
16using namespace jana;
17
18#include <DAQ/DF1TDCHit.h>
19#include <DAQ/Df250PulseIntegral.h>
20#include <DAQ/JEventSource_EVIO.h>
21#include <TTAB/DTranslationTable.h>
22
23#include <TDirectory.h>
24#include <TH2.h>
25#include <TH1.h>
26#include <TProfile.h>
27#include <TProfile2D.h>
28#include <TROOT.h>
29
30
31static const int highcratenum=100;
32// root hist pointers
33static TH2I *daq_occ_crates[highcratenum];
34static TProfile2D *daq_ped_crates[highcratenum];
35static TProfile2D *daq_TDClocked_crates[highcratenum];
36static TProfile2D *daq_TDCovr_crates[highcratenum];
37static TProfile *daq_hits_per_event;
38static TProfile *daq_words_per_event;
39static TH1D *daq_event_size;
40static TH1D *daq_event_tdiff;
41static TH1D *daq_words_by_type;
42static bool ttab_labels_set = false;
43
44// Routine used to create our JEventProcessor
45extern "C"{
46 void InitPlugin(JApplication *app){
47 InitJANAPlugin(app);
48 app->AddProcessor(new JEventProcessor_DAQ_online());
49 }
50} // "C"
51
52
53//------------------
54// JEventProcessor_DAQ_online (Constructor)
55//------------------
56JEventProcessor_DAQ_online::JEventProcessor_DAQ_online()
57{
58
59}
60
61//------------------
62// ~JEventProcessor_DAQ_online (Destructor)
63//------------------
64JEventProcessor_DAQ_online::~JEventProcessor_DAQ_online()
65{
66
67}
68
69//------------------
70// init
71//------------------
72jerror_t JEventProcessor_DAQ_online::init(void)
73{
74 printf("JEventProcessor_DAQ_online::init()\n");
75
76 // create root folder for DAQ and cd to it, store main dir
77 maindir = gDirectory(TDirectory::CurrentDirectory());
78 daqdir = maindir->mkdir("DAQ");
79 daqdir->cd();
80
81 // Initialise histograms and variables
82 for (int i=0; i<highcratenum; i++) {
83 daq_occ_crates[i] = NULL__null;
84 daq_ped_crates[i] = NULL__null;
85 daq_TDClocked_crates[i] = NULL__null;
86 daq_TDCovr_crates[i] = NULL__null;
87 }
88
89 daq_hits_per_event = new TProfile("daq_hits_per_event", "Hits/event vs. rocid", 100, 0.5, 100.5);
90 daq_words_per_event = new TProfile("daq_words_per_event", "words/event vs. rocid", 100, 0.5, 100.5);
91 daq_event_size = new TH1D("daq_event_size", "Event size in kB", 1000, 0.0, 1.0E3);
92 daq_event_tdiff = new TH1D("daq_event_tdiff", "Time between events", 10000, 0.0, 1.0E2);
93 daq_words_by_type = new TH1D("daq_words_by_type", "Number of words in EVIO file by type", kNEVIOWordTypes, 0, (double)kNEVIOWordTypes);
94
95 daq_words_per_event->GetXaxis()->SetBinLabel(1 ,"Trigger Bank");
96 daq_words_per_event->GetXaxis()->SetBinLabel(99 ,"Residual");
97
98 daq_event_size->SetXTitle("Total event size (kB)");
99 daq_event_tdiff->SetXTitle("#deltat between events (ms)");
100
101 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kUnknown, "unknown");
102 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kEVIOHeader, "EVIO len. & header");
103 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kEVIOEventNumber, "Event Number Word");
104 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kEVIOTimestamp, "Timestamp");
105
106 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kBORData, "BOR record");
107
108 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250BlockHeader, "f250 Block Header");
109 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250BlockTrailer, "f250 Block Trailer");
110 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250EventHeader, "f250 Event Header");
111 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250TriggerTime, "f250 Trigger Time");
112 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250WindowRawData, "f250 Window Raw Data");
113 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250WindowSum, "f250 Window Sum");
114 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250PulseRawData, "f250 Pulse Raw Data");
115 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250PulseIntegral, "f250 Pulse Integral");
116 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250PulseTime, "f250 Pulse Time");
117 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250PulsePedestal, "f250 Pulse Pedestal");
118 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250EventTrailer, "f250 Event Trailer");
119 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250DataNotValid, "f250 Data Not Valid");
120 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf250Filler, "f250 Filler Word");
121
122 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125BlockHeader, "f125 Block Header");
123 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125BlockTrailer, "f125 Block Trailer");
124 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125EventHeader, "f125 Event Header");
125 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125TriggerTime, "f125 Trigger Time");
126 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125WindowRawData, "f125 Window Raw Data");
127 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125CDCPulse, "f125 CDC Pulse");
128 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125FDCPulse6, "f125 FDC Pulse (integral)");
129 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125FDCPulse9, "f125 FDC Pulse (peak)");
130 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125PulseIntegral, "f125 Pulse Integral");
131 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125PulseTime, "f125 Pulse Time");
132 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125PulsePedestal, "f125 Pulse Pedestal");
133 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125EventTrailer, "f125 Event Trailer");
134 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125DataNotValid, "f125 Data Not Valid");
135 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kf125Filler, "f125 Filler Word");
136
137 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2BlockHeader, "F1v2 Block Header");
138 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2BLockTrailer, "F1v2 Block Trailer");
139 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2EventHeader, "F1v2 Event Header");
140 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2TriggerTime, "F1v2 Trigger Time");
141 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2ChipHeader, "F1v2 Chip Header");
142 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2Data, "F1v2 Data");
143 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2Filler, "F1v2 Filler");
144 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v2BreakWord, "F1v2 Break Word");
145
146 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3BlockHeader, "F1v3 Block Header");
147 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3BLockTrailer, "F1v3 Block Trailer");
148 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3EventHeader, "F1v3 Event Header");
149 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3TriggerTime, "F1v3 Trigger Time");
150 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3ChipHeader, "F1v3 Chip Header");
151 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3Data, "F1v3 Data");
152 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3Filler, "F1v3 Filler");
153 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF1v3BreakWord, "F1v3 Break Word");
154
155 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190GlobalHeader, "CAEN1190 GLobal Header");
156 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190GlobalTrailer, "CAEN1190 Global Trailer");
157 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190GlobalTriggerTime, "CAEN1190 Trigger Time");
158 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190TDCHeader, "CAEN1190 TDC Header");
159 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190TDCData, "CAEN1190 TDC Data");
160 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190TDCError, "CAEN1190 TDC Error");
161 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190TDCTrailer, "CAEN1190 TDC Trailer");
162 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kCAEN1190Filler, "CAEN1190 Filler");
163
164 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kConfig, "DAQ Config");
165 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kConfigf250, "DAQ Config f250");
166 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kConfigf125, "DAQ Config f125");
167 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kConfigF1, "DAQ Config F1");
168 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kConfigCAEN1190, "DAQ Config CAEN1190");
169
170 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kEPICSheader, "EPICS header");
171 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kEPICSdata, "EPICS data");
172
173 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kF800FAFA, "0xf800fafa");
174 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kD00DD00D, "0xd00dd00d");
175
176 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kTotWords, "Total words in all events");
177 daq_words_by_type->GetXaxis()->SetBinLabel(1 + kNevents, "Number of events");
178
179 // back to main dir
180 maindir->cd();
181
182 return NOERROR;
183}
184
185//------------------
186// AddROCIDLabels
187//------------------
188void JEventProcessor_DAQ_online::AddROCIDLabels(JEventLoop *loop)
189{
190 /// This is called just once to set the x-axis labels
191 /// of histograms whose x-axis is the rocid so that we
192 /// can label them by detector.
193
194 // NO lock: called in init()
195 const DTranslationTable *ttab = NULL__null;
196 loop->GetSingle(ttab);
197
198 // Loop over all rocid values
199 for(uint32_t rocid=2; rocid<99; rocid++){
200 // We don't actually know what slot/channel combos are defined
201 // for this so we loop until we find one.
202 bool found_chan = false;
203 daq_hits_per_event->GetXaxis()->SetBinLabel(rocid, "");
204 daq_words_per_event->GetXaxis()->SetBinLabel(rocid, "");
205 for(uint32_t slot=2; slot<24; slot++){
206 for(uint32_t channel=0; channel<3; channel++){
207 try{
208 DTranslationTable::csc_t csc = {rocid, slot, channel};
209 const DTranslationTable::DChannelInfo &chinfo = ttab->GetDetectorIndex(csc);
210 daq_hits_per_event->GetXaxis()->SetBinLabel(rocid, ttab->DetectorName(chinfo.det_sys).c_str());
211 daq_words_per_event->GetXaxis()->SetBinLabel(rocid, ttab->DetectorName(chinfo.det_sys).c_str());
212 found_chan = true;
213 break;
214 }catch(JException &e){
215 // Do nothing
216 }
217 }
218 if(found_chan) break;
219 }
220 }
221}
222
223//------------------
224// brun
225//------------------
226jerror_t JEventProcessor_DAQ_online::brun(JEventLoop *eventLoop, int32_t runnumber)
227{
228 // This is called whenever the run number changes
229 return NOERROR;
230}
231
232//------------------
233// evnt
234//------------------
235jerror_t JEventProcessor_DAQ_online::evnt(JEventLoop *loop, uint64_t eventnumber)
236{
237 // This is called for every event. Use of common resources like writing
238 // to a file or filling a histogram should be mutex protected. Using
239 // loop->Get(...) to get reconstructed objects (and thereby activating the
240 // reconstruction algorithm) should be done outside of any mutex lock
241 // since multiple threads may call this method at the same time.
242 // Here's an example:
243 vector<const DF1TDCHit*> f1tdchits;
244 vector<const Df250PulseIntegral*> f250PIs;
245 vector<const Df125PulseIntegral*> f125PIs;
246 vector<const Df125CDCPulse*> f125CDCs;
247 vector<const Df125FDCPulse*> f125FDCs;
248 vector<const DCAEN1290TDCHit*> caen1290hits;
249
250 loop->Get(f1tdchits);
251 loop->Get(f250PIs);
252 loop->Get(f125PIs);
253 loop->Get(f125CDCs);
254 loop->Get(f125FDCs);
255 loop->Get(caen1290hits);
256
257 ParseEventSize(loop->GetJEvent());
258
259 // Set rocid histogram labels based on detector if needed
260 if(!ttab_labels_set){
261 ttab_labels_set = true;
262 AddROCIDLabels(loop);
263 }
264
265 // Initialize counters looking at num. hits per crate
266 uint32_t Nhits_rocid[101];
267 for(uint32_t rocid=0; rocid<101; rocid++) Nhits_rocid[rocid] = 0;
268
269 // Although we are only filling objects local to this plugin, The directory changes: Global ROOT lock
270 japp->RootWriteLock(); //ACQUIRE ROOT LOCK
271
272 if (daqdir!=NULL__null) daqdir->cd();
273
274
275 // Access TDC from DF1TDCHit object
276 for(unsigned int i=0; i<f1tdchits.size(); i++) {
277 const DF1TDCHit *hit = f1tdchits[i];
278 int rocid = hit->rocid;
279 int slot = hit->slot;
280 int channel = hit->channel;
281 int data_word = hit->data_word;
282
283 if(rocid>=0 && rocid<=100) Nhits_rocid[rocid]++;
284
285 if (daq_occ_crates[rocid]==NULL__null) {
286 printf("JEventProcessor_DAQ_online::evnt creating occupancy histogram for crate %i\n",rocid);
287 char cratename[255],title[255];
288 sprintf(cratename,"daq_occ_crate%i",rocid);
289 sprintf(title,"Crate %i occupancy (TDC);Slot;Channel",rocid);
290 daq_occ_crates[rocid] = new TH2I(cratename,title,21,0.5,21.5,32,-0.5,31.5);
291 daq_occ_crates[rocid]->SetStats(0);
292 sprintf(cratename,"daq_TDClocked_crate%i",rocid);
293 sprintf(title,"Crate %i TDC lock status (TDC);Slot;Channel",rocid);
294 daq_TDClocked_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,32,-0.5,31.5);
295 daq_TDClocked_crates[rocid]->SetStats(0);
296 sprintf(cratename,"daq_TDCovr_crate%i",rocid);
297 sprintf(title,"Crate %i TDC overflow status (TDC);Slot;Channel",rocid);
298 daq_TDCovr_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,32,-0.5,31.5);
299 daq_TDCovr_crates[rocid]->SetStats(0);
300
301 }
302 daq_occ_crates[rocid]->Fill(slot,channel);
303 daq_TDClocked_crates[rocid]->Fill(slot,channel,(data_word>>26)&(1));
304 daq_TDCovr_crates[rocid]->Fill(slot,channel,(data_word>>25)&(1));
305 daq_TDCovr_crates[rocid]->Fill(slot,channel,(data_word>>24)&(1));
306 }
307
308 // Access F250 from Df250PulseIntegral object
309 for(unsigned int i=0; i<f250PIs.size(); i++) {
310 const Df250PulseIntegral *hit = f250PIs[i];
311 int rocid = hit->rocid;
312 int slot = hit->slot;
313 int channel = hit->channel;
314
315 if(rocid>=0 && rocid<=100) {
316 Nhits_rocid[rocid]++;
317
318 if (daq_occ_crates[rocid]==NULL__null) {
319 printf("JEventProcessor_DAQ_online::evnt creating occupancy histogram for crate %i\n",rocid);
320 char cratename[255],title[255];
321 sprintf(cratename,"daq_occ_crate%i",rocid);
322 sprintf(title,"Crate %i occupancy (F250);Slot;Channel",rocid);
323 daq_occ_crates[rocid] = new TH2I(cratename,title,21,0.5,21.5,16,-0.5,15.5);
324 daq_occ_crates[rocid]->SetStats(0);
325 }
326 daq_occ_crates[rocid]->Fill(slot,channel);
327
328 if (daq_ped_crates[rocid]==NULL__null) {
329 printf("JEventProcessor_DAQ_online::evnt creating pedestal histogram for crate %i\n",rocid);
330 char cratename[255],title[255];
331 sprintf(cratename,"daq_ped_crate%i",rocid);
332 sprintf(title,"Crate %i Average Pedestal (F250);Slot;Channel",rocid);
333 daq_ped_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,16,-0.5,15.5);
334 daq_ped_crates[rocid]->SetStats(0);
335 }
336 if (hit->pedestal > 0) {
337 daq_ped_crates[rocid]->Fill(slot,channel,hit->pedestal);
338 }
339 }
340 }
341
342 // Access F125 from Df125PulseIntegral object
343 for(unsigned int i=0; i<f125PIs.size(); i++) {
344 const Df125PulseIntegral *hit = f125PIs[i];
345 int rocid = hit->rocid;
346 int slot = hit->slot;
347 int channel = hit->channel;
348
349 if(hit->emulated) continue; // ignore emulated hits
350
351 if(rocid>=0 && rocid<=100) {
352 Nhits_rocid[rocid]++;
353
354 if (daq_occ_crates[rocid]==NULL__null) {
355 printf("JEventProcessor_DAQ_online::evnt creating occupancy histogram for crate %i\n",rocid);
356 char cratename[255],title[255];
357 sprintf(cratename,"daq_occ_crate%i",rocid);
358 sprintf(title,"Crate %i occupancy (F125);Slot;Channel",rocid);
359 daq_occ_crates[rocid] = new TH2I(cratename,title,21,0.5,21.5,16,-0.5,15.5);
360 daq_occ_crates[rocid]->SetStats(0);
361 }
362 daq_occ_crates[rocid]->Fill(slot,channel);
363
364 if (daq_ped_crates[rocid]==NULL__null) {
365 printf("JEventProcessor_DAQ_online::evnt creating pedestal histogram for crate %i\n",rocid);
366 char cratename[255],title[255];
367 sprintf(cratename,"daq_ped_crate%i",rocid);
368 sprintf(title,"Crate %i Average Pedestal (F125);Slot;Channel",rocid);
369 daq_ped_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,16,-0.5,15.5);
370 daq_ped_crates[rocid]->SetStats(0);
371 }
372 if (hit->pedestal > 0) {
373 daq_ped_crates[rocid]->Fill(slot,channel,hit->pedestal);
374 }
375 }
376 }
377
378 // Access F125 from Df125CDCPulse object
379 for(unsigned int i=0; i<f125CDCs.size(); i++) {
380 const Df125CDCPulse *hit = f125CDCs[i];
381 int rocid = hit->rocid;
382 int slot = hit->slot;
383 int channel = hit->channel;
384
385 if(rocid>=0 && rocid<=100) {
386 Nhits_rocid[rocid]++;
387
388 if (daq_occ_crates[rocid]==NULL__null) {
389 printf("JEventProcessor_DAQ_online::evnt creating occupancy histogram for crate %i\n",rocid);
390 char cratename[255],title[255];
391 sprintf(cratename,"daq_occ_crate%i",rocid);
392 sprintf(title,"Crate %i occupancy (F125);Slot;Channel",rocid);
393 daq_occ_crates[rocid] = new TH2I(cratename,title,21,0.5,21.5,16,-0.5,15.5);
394 daq_occ_crates[rocid]->SetStats(0);
395 }
396 daq_occ_crates[rocid]->Fill(slot,channel);
397
398 if (daq_ped_crates[rocid]==NULL__null) {
399 printf("JEventProcessor_DAQ_online::evnt creating pedestal histogram for crate %i\n",rocid);
400 char cratename[255],title[255];
401 sprintf(cratename,"daq_ped_crate%i",rocid);
402 sprintf(title,"Crate %i Average Pedestal (F125);Slot;Channel",rocid);
403 daq_ped_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,16,-0.5,15.5);
404 daq_ped_crates[rocid]->SetStats(0);
405 }
406 if (hit->pedestal > 0) {
407 daq_ped_crates[rocid]->Fill(slot,channel,hit->pedestal);
408 }
409 }
410 }
411
412 // Access F125 from Df125FDCPulse object
413 for(unsigned int i=0; i<f125FDCs.size(); i++) {
414 const Df125FDCPulse *hit = f125FDCs[i];
415 int rocid = hit->rocid;
416 int slot = hit->slot;
417 int channel = hit->channel;
418
419 if(rocid>=0 && rocid<=100) {
420 Nhits_rocid[rocid]++;
421
422 if (daq_occ_crates[rocid]==NULL__null) {
423 printf("JEventProcessor_DAQ_online::evnt creating occupancy histogram for crate %i\n",rocid);
424 char cratename[255],title[255];
425 sprintf(cratename,"daq_occ_crate%i",rocid);
426 sprintf(title,"Crate %i occupancy (F125);Slot;Channel",rocid);
427 daq_occ_crates[rocid] = new TH2I(cratename,title,21,0.5,21.5,16,-0.5,15.5);
428 daq_occ_crates[rocid]->SetStats(0);
429 }
430 daq_occ_crates[rocid]->Fill(slot,channel);
431
432 if (daq_ped_crates[rocid]==NULL__null) {
433 printf("JEventProcessor_DAQ_online::evnt creating pedestal histogram for crate %i\n",rocid);
434 char cratename[255],title[255];
435 sprintf(cratename,"daq_ped_crate%i",rocid);
436 sprintf(title,"Crate %i Average Pedestal (F125);Slot;Channel",rocid);
437 daq_ped_crates[rocid] = new TProfile2D(cratename,title,21,0.5,21.5,16,-0.5,15.5);
438 daq_ped_crates[rocid]->SetStats(0);
439 }
440 if (hit->pedestal > 0) {
441 daq_ped_crates[rocid]->Fill(slot,channel,hit->pedestal);
442 }
443 }
444 }
445
446 // Access CAEN1290 TDC hits
447 for(unsigned int i=0; i<caen1290hits.size(); i++) {
448 const DCAEN1290TDCHit *hit = caen1290hits[i];
449 int rocid = hit->rocid;
450 //int slot = hit->slot;
451 //int channel = hit->channel;
452
453 if(rocid>=0 && rocid<=100) Nhits_rocid[rocid]++;
454 }
455
456 // Fill in hits by crate
457 for(uint32_t rocid=0; rocid<101; rocid++) daq_hits_per_event->Fill(rocid, Nhits_rocid[rocid]);
458
459 maindir->cd();
460 // Unlock ROOT
461 japp->RootUnLock(); //RELEASE ROOT LOCK
462
463 return NOERROR;
464}
465
466//------------------
467// ParseEventSize
468//------------------
469void JEventProcessor_DAQ_online::ParseEventSize(JEvent &event)
470{
471 /// This ugliness is needed to get at the true banks for each event by rocid.
472
473 // Bombproof
474 if(event.GetJEventSource()->className() != string("JEventSource_EVIO")){
475 static bool warned = false;
476 if(!warned){
477 cout << "WARNING: This is not an event source of type JEventSource_EVIO!" << endl;
478 cout << " Event size statistics filling unavailable!" << endl;
479 warned = true;
480 }
481 return;
482 }
483
484 void *ref = event.GetRef();
485 if(!ref) return;
486 uint32_t *istart = JEventSource_EVIO::GetEVIOBufferFromRef(ref);
487 uint32_t evio_buffsize = JEventSource_EVIO::GetEVIOBufferSizeFromRef(ref);
488 uint32_t evio_buffwords = evio_buffsize/sizeof(uint32_t);
489 uint32_t *iend = &istart[evio_buffwords];
490
491 if( istart==NULL__null ) return;
492 if( (evio_buffwords>=10) && (istart[7]==0xc0da0100) ){
493 // NTH is first 8 words so skip them
494 istart= &istart[8];
495 evio_buffsize -= 8*sizeof(uint32_t);
Value stored to 'evio_buffsize' is never read
496 evio_buffwords -= 8;
497 }
498
499 // Check if this is BOR data
500 if( evio_buffwords >= 4 ){
501 uint32_t mask = (0x70<<16) | (0x01);
502 if( (istart[1]&mask) == mask ){
503
504 // FILL HISTOGRAMS
505 // Since we are filling histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
506 japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
507 daq_words_by_type->Fill(kBORData, istart[0]/sizeof(uint32_t));
508 japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK
509 return; // no further parsing needed
510 }
511 }
512
513 // Check if this is EPICS data
514 if( evio_buffwords >= 4 ){
515 if( istart[1] == (0x60<<16) + (0xD<<8) + (0x1<<0) ){
516 if( istart[2] == (0x61<<24) + (0x1<<16) + (0x1<<0) ){
517
518 // FILL HISTOGRAMS
519 // Since we are filling histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
520 japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
521 daq_words_by_type->Fill(kEPICSheader, 3.0); // EVIO outer and segment headers + timestamp
522 daq_words_by_type->Fill(kEPICSdata, istart[0]/sizeof(uint32_t) - 3);
523 japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK
524 return; // no further parsing needed
525 }
526 }
527 }
528
529 if( evio_buffwords < 4 ){
530 cout << "Too few words in event (" << evio_buffwords << ") skipping..." << endl;
531 return;
532 }
533
534 // Physics event length
535 uint32_t physics_event_len = istart[0];
536 if( (istart[1] & 0xFF001000) != 0xFF001000 ) return; // not a physics event
537 if( physics_event_len+1 > evio_buffwords ){
538 cout << "Too many words in physics event: " << physics_event_len+1 << " > " << evio_buffwords << endl;
539 return;
540 }
541
542 // Trigger bank event length
543 uint32_t trigger_bank_len = istart[2];
544 if( (istart[3] & 0xFF202000) != 0xFF202000 ) return; // not a trigger bank
545 uint64_t tlo = istart[2+5];
546 uint64_t thi = istart[2+6];
547 uint64_t timestamp = (thi<<32) + (tlo<<0);
548 if( trigger_bank_len+2 > evio_buffwords ){
549 cout << "Too many words in trigger bank " << trigger_bank_len << " > " << evio_buffwords-2 << endl;
550 return;
551 }
552
553 // Allocate memory to hold stats data
554 uint32_t Nwords[100]; // total data words for each ROC (includes event length words)
555 uint32_t word_stats[kNEVIOWordTypes]; // obtained from parsing event
556 for(uint32_t rocid=0; rocid<100; rocid++) Nwords[rocid] = 0;
557 for(uint32_t i=0; i<kNEVIOWordTypes; i++) word_stats[i] = 0;
558
559 word_stats[kNevents]++;
560 word_stats[kTotWords] += evio_buffwords;
561
562 word_stats[kEVIOHeader] += 4; // physics event and built trigger bank length and header words
563
564 // Loop over data banks
565 uint32_t *iptr = &istart[3+trigger_bank_len];
566 while(iptr < iend){
567
568 uint32_t len = *iptr;
569 uint32_t rocid = (iptr[1]>>16) & 0XFF;
570
571 if(rocid<100) Nwords[rocid] += len+1;
572
573 word_stats[kEVIOHeader] += 2; // ROC data bank length and header words
574
575 uint32_t *imyend = &iptr[len+1];
576 if(imyend > iend) imyend = iend;
577
578 uint64_t Nwords = ((uint64_t)imyend - (uint64_t)iptr)/sizeof(uint32_t);
579 if(Nwords<2){
580 static int Nwarnings = 0;
581 if(Nwarnings<10){
582 cout << "Nwords<2 (?)" << endl;
583 cout << " evio_buffwords = " << evio_buffwords << endl;
584 cout << " physics_event_len = " << physics_event_len << endl;
585 cout << " trigger_bank_len = " << trigger_bank_len << endl;
586 event.Print();
587 if(++Nwarnings == 10) cout << "Last warning!" << endl;
588 }
589 break;
590 }
591
592 DataWordStats(iptr, imyend, word_stats);
593
594 iptr = &iptr[len +1];
595 }
596
597 // FILL HISTOGRAMS
598 // Since we are filling histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
599 japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
600
601 // Calculating time between events is tricky when using multiple-threads.
602 // We need the timestamp of two sequential events, but the order in which
603 // they are processed here will likely not be in event order. Thus, we keep
604 // a running list of the last 128 timestamps and event numbers seen by this
605 // routine. We can then search this for the event prior to this one and if
606 // found, use it.
607 uint32_t event_num = event.GetEventNumber();
608 static uint32_t recent_event_nums[128];
609 static uint64_t recent_timestamps[128];
610 static uint32_t ievent = 0;
611 for(uint32_t i=0; i<ievent; i++){
612 if(i>=128) break;
613 if(recent_event_nums[i] == (event_num-1)){
614 double tdiff = (double)(timestamp - recent_timestamps[i])/250.0E6; // convert to seconds
615 daq_event_tdiff->Fill(tdiff*1000.0); // ms
616 break;
617 }
618 }
619
620 // Record this timestamp/event number in the ring buffer
621 uint32_t idx = ievent%128;
622 recent_event_nums[idx] = event_num;
623 recent_timestamps[idx] = timestamp;
624 ievent++;
625
626 // Fill event size histos
627 double physics_event_len_kB = (double)((physics_event_len+1)*sizeof(uint32_t))/1024.0;
628 daq_event_size->Fill(physics_event_len_kB);
629 uint32_t TotalWords = 0;
630 for(uint32_t rocid=0; rocid<100; rocid++){
631 daq_words_per_event->Fill(rocid, Nwords[rocid]);
632 TotalWords += Nwords[rocid];
633 }
634
635 daq_words_per_event->Fill(1, trigger_bank_len+1);
636 daq_words_per_event->Fill(99, physics_event_len - trigger_bank_len - TotalWords);
637
638 for(uint32_t i=0; i<kNEVIOWordTypes; i++){
639 daq_words_by_type->Fill(i, (double)word_stats[i]);
640 }
641
642 japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK
643}
644
645//------------------
646// DataWordStats
647//------------------
648void JEventProcessor_DAQ_online::DataWordStats(uint32_t *iptr, uint32_t *iend, uint32_t *word_stats)
649{
650 // Upon entry, the iptr will point to the start of the "Physics Event's Data Bank".
651 // It will loop over all sub-banks, tallying the word count as it goes up to
652 // but not including iend.
653
654 iptr++; // advance past length word
655 uint32_t rocid = (*iptr++)>>16 & 0x0FFF;
656 while(iptr < iend){
657 uint32_t data_block_bank_len = *iptr++;
658 uint32_t *iendbank = &iptr[data_block_bank_len];
659 uint32_t det_id = ((*iptr) >> 16) & 0x0FFF;
660
661 if(iendbank > iend) iendbank = iend;
662
663 word_stats[kEVIOHeader] += 2; // data block bank length and header words
664
665 iptr++; // advance to first raw data word
666
667 uint32_t Ntoprocess = data_block_bank_len - 1; // 1 for bank header
668
669#if 0 // I don't know if these words are actually implmented ??
670 word_stats[kEVIOEventNumber]++; // starting event number
671 word_stats[kEVIOTimestamp] += 2; // 48-bit timestamp
672 iptr++; // starting event number
673 iptr++; // 48-bit timestamp
674 iptr++; // 48-bit timestamp
675 Ntoprocess -= 3;
676#endif
677 uint32_t *irawdata = iptr;
678
679 switch(det_id){
680 case 0:
681 case 1:
682 case 3:
683 case 6: // flash 250 module, MMD 2014/2/4
684 case 16: // flash 125 module (CDC), DL 2014/6/19
685 case 26: // F1 TDC module (BCAL), MMD 2014-07-31
686 ParseJLabModuleData(rocid, iptr, iendbank, word_stats);
687 break;
688
689 case 20:
690 ParseCAEN1190(rocid, iptr, iendbank, word_stats);
691 break;
692
693 case 0x55:
694 ParseModuleConfiguration(rocid, iptr, iendbank, word_stats);
695 break;
696
697 default:
698 break;
699 }
700
701 uint32_t Nprocessed = (uint32_t)((uint64_t)iptr - (uint64_t)irawdata)/sizeof(uint32_t);
702 if(Nprocessed < Ntoprocess) word_stats[kUnknown] += Ntoprocess - Nprocessed;
703 iptr = iendbank;
704 }
705
706
707}
708
709//------------------
710// ParseJLabModuleData
711//------------------
712void JEventProcessor_DAQ_online::ParseJLabModuleData(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
713{
714 while(iptr < iend){
715 if(*iptr != 0xf800fafa) break;
716 word_stats[kF800FAFA]++;
717 iptr++;
718 }
719
720 uint32_t mod_id = ((*iptr) >> 18) & 0x000F;
721 switch(mod_id){
722 case DModuleType::FADC250: Parsef250Bank(rocid, iptr, iend, word_stats); break;
723 case DModuleType::FADC125: Parsef125Bank(rocid, iptr, iend, word_stats); break;
724 case DModuleType::F1TDC32: ParseF1v2TDCBank(rocid, iptr, iend, word_stats); break;
725 case DModuleType::F1TDC48: ParseF1v3TDCBank(rocid, iptr, iend, word_stats); break;
726 //case DModuleType::JLAB_TS: ParseTSBank(rocid, iptr, iend, word_stats); break;
727 //case DModuleType::TID: ParseTIBank(rocid, iptr, iend, word_stats); break;
728 }
729}
730
731//------------------
732// Parsef250Bank
733//------------------
734void JEventProcessor_DAQ_online::Parsef250Bank(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
735{
736 while(iptr<iend){
737
738 if(((*iptr>>31) & 0x1) == 0) { word_stats[kUnknown]++ ; iptr++; continue;}
739
740 uint32_t window_width;
741 uint32_t window_words;
742 uint32_t data_type = (*iptr>>27) & 0x0F;
743 switch(data_type){
744 case 0: word_stats[kf250BlockHeader]++; iptr++; break;
745 case 1: word_stats[kf250BlockTrailer]++; iptr++; break;
746 case 2: word_stats[kf250EventHeader]++; iptr++; break;
747 case 3: // Trigger time
748 word_stats[kf250TriggerTime]++;
749 iptr++;
750 if(((*iptr>>31) & 0x1) == 0){ word_stats[kf250TriggerTime]++; iptr++; }
751 break;
752 case 4: // Window Raw Data
753 window_width = (*iptr>>0) & 0x0FFF;
754 window_words = 1 + ((window_width+1)/2); // 1 is for header word + 2 sample per word
755 word_stats[kf250WindowRawData] += window_words;
756 iptr = &iptr[window_words];
757 break;
758 case 7: word_stats[kf250PulseIntegral]++; iptr++; break;
759 case 8: word_stats[kf250PulseTime]++; iptr++; break;
760 case 10: word_stats[kf250PulsePedestal]++; iptr++; break;
761 case 13: word_stats[kf250EventTrailer]++; iptr++; break;
762 case 14: word_stats[kf250DataNotValid]++; iptr++; break;
763 case 15: word_stats[kf250Filler]++; iptr++; break;
764
765 default: word_stats[kUnknown]++; iptr++; break;
766 }
767 }
768}
769
770//------------------
771// Parsef125Bank
772//------------------
773void JEventProcessor_DAQ_online::Parsef125Bank(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
774{
775 while(iptr<iend){
776
777 if(((*iptr>>31) & 0x1) == 0) { word_stats[kUnknown]++ ; iptr++; continue;}
778
779 uint32_t window_width;
780 uint32_t window_words;
781 uint32_t data_type = (*iptr>>27) & 0x0F;
782 switch(data_type){
783 case 0: word_stats[kf125BlockHeader]++; iptr++; break;
784 case 1: word_stats[kf125BlockTrailer]++; iptr++; break;
785 case 2: word_stats[kf125EventHeader]++; iptr++; break;
786 case 3: // Trigger time
787 word_stats[kf125TriggerTime]++;
788 iptr++;
789 if(((*iptr>>31) & 0x1) == 0){ word_stats[kf125TriggerTime]++; iptr++; }
790 break;
791 case 4: // Window Raw Data
792 window_width = (*iptr>>0) & 0x0FFF;
793 window_words = 1 + ((window_width+1)/2); // 1 is for header word + 2 sample per word
794 word_stats[kf125WindowRawData] += window_words;
795 iptr = &iptr[window_words];
796 break;
797 case 5: word_stats[kf125CDCPulse]++;
798 iptr++;
799 if(((*iptr>>31) & 0x1) == 0){ word_stats[kf125CDCPulse]++; iptr++; }
800 break;
801 case 6: word_stats[kf125FDCPulse6]++;
802 iptr++;
803 if(((*iptr>>31) & 0x1) == 0){ word_stats[kf125FDCPulse6]++; iptr++; }
804 break;
805 case 7: word_stats[kf125PulseIntegral]++; iptr++; break;
806 case 8: word_stats[kf125PulseTime]++; iptr++; break;
807 case 9: word_stats[kf125FDCPulse9]++;
808 iptr++;
809 if(((*iptr>>31) & 0x1) == 0){ word_stats[kf125FDCPulse9]++; iptr++; }
810 break;
811 case 10: word_stats[kf125PulsePedestal]++; iptr++; break;
812 case 13: word_stats[kf125EventTrailer]++; iptr++; break;
813 case 14: word_stats[kf125DataNotValid]++; iptr++; break;
814 case 15: word_stats[kf125Filler]++; iptr++; break;
815
816 default: word_stats[kUnknown]++; iptr++; break;
817 }
818 }
819}
820
821//------------------
822// ParseF1v2TDCBank
823//------------------
824void JEventProcessor_DAQ_online::ParseF1v2TDCBank(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
825{
826 while(iptr<iend){
827 switch( (*iptr++) & 0xF8000000 ){
828 case 0xC0000000: word_stats[kF1v2ChipHeader]++; break;
829 case 0xB8000000: word_stats[kF1v2Data]++; break;
830 case 0xF8000000: word_stats[kF1v2Filler]++; break;
831 case 0x80000000: word_stats[kF1v2BlockHeader]++; break;
832 case 0x88000000: word_stats[kF1v2BLockTrailer]++; break;
833 case 0x90000000: word_stats[kF1v2EventHeader]++; break;
834 case 0x98000000: word_stats[kF1v2TriggerTime]++; break;
835 case 0xF0000000: word_stats[kF1v2BreakWord]++; break;
836 default: word_stats[kUnknown]++; break;
837 }
838 }
839}
840
841//------------------
842// ParseF1v3TDCBank
843//------------------
844void JEventProcessor_DAQ_online::ParseF1v3TDCBank(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
845{
846 while(iptr<iend){
847 switch( (*iptr++) & 0xF8000000 ){
848 case 0xC0000000: word_stats[kF1v3ChipHeader]++; break;
849 case 0xB8000000: word_stats[kF1v3Data]++; break;
850 case 0xF8000000: word_stats[kF1v3Filler]++; break;
851 case 0x80000000: word_stats[kF1v3BlockHeader]++; break;
852 case 0x88000000: word_stats[kF1v3BLockTrailer]++; break;
853 case 0x90000000: word_stats[kF1v3EventHeader]++; break;
854 case 0x98000000: word_stats[kF1v3TriggerTime]++; break;
855 case 0xF0000000: word_stats[kF1v3BreakWord]++; break;
856 default: word_stats[kUnknown]++; break;
857 }
858 }
859}
860
861//------------------
862// ParseCAEN1190
863//------------------
864void JEventProcessor_DAQ_online::ParseCAEN1190(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
865{
866 while(iptr<iend){
867
868 // This word appears to be appended to the data.
869 // Probably in the ROL. Ignore it if found.
870 if(*iptr == 0xd00dd00d) {
871 word_stats[kD00DD00D]++;
872 iptr++;
873 continue;
874 }
875
876 uint32_t type = (*iptr++) >> 27;
877 switch(type){
878 case 0b01000: word_stats[kCAEN1190GlobalHeader]++; break;
879 case 0b10000: word_stats[kCAEN1190GlobalTrailer]++; break;
880 case 0b10001: word_stats[kCAEN1190GlobalTriggerTime]++; break;
881 case 0b00001: word_stats[kCAEN1190TDCHeader]++; break;
882 case 0b00000: word_stats[kCAEN1190TDCData]++; break;
883 case 0b00100: word_stats[kCAEN1190TDCError]++; break;
884 case 0b00011: word_stats[kCAEN1190TDCTrailer]++; break;
885 case 0b11000: word_stats[kCAEN1190Filler]++; break;
886 default: word_stats[kUnknown]++; break;
887 }
888 }
889}
890
891//------------------
892// ParseModuleConfiguration
893//------------------
894void JEventProcessor_DAQ_online::ParseModuleConfiguration(uint32_t rocid, uint32_t *&iptr, uint32_t *iend, uint32_t *word_stats)
895{
896 while(iptr < iend){
897
898 word_stats[kConfig]++; // Count headers as generic
899 uint32_t Nvals = ((*iptr++) >> 24) & 0xFF;
900
901 // Loop over all parameters in this section
902 for(uint32_t i=0; i< Nvals; i++){
903
904 switch((*iptr++)>>24){
905 case 0x05: word_stats[kConfigf250]++; break;
906 case 0x0F: word_stats[kConfigf125]++; break;
907 case 0x06: word_stats[kConfigF1]++; break;
908 case 0x10: word_stats[kConfigCAEN1190]++; break;
909 default: word_stats[kConfig]++; break;
910 }
911 }
912 }
913}
914
915//------------------
916// erun
917//------------------
918jerror_t JEventProcessor_DAQ_online::erun(void)
919{
920 // This is called whenever the run number changes, before it is
921 // changed to give you a chance to clean up before processing
922 // events from the next run number.
923
924 // FILL HISTOGRAMS
925 // Since we are filling histograms local to this plugin, it will not interfere with other ROOT operations: can use plugin-wide ROOT fill lock
926 japp->RootFillLock(this); //ACQUIRE ROOT FILL LOCK
927
928 for (int i=0; i<highcratenum; i++) {
929 if (daq_occ_crates[i] != NULL__null) {
930 daq_occ_crates[i]->SetMinimum(daq_occ_crates[i]->GetMinimum(0.001));
931 }
932 if (daq_ped_crates[i] != NULL__null) {
933 daq_ped_crates[i]->SetMinimum(daq_ped_crates[i]->GetMinimum(0.001));
934 }
935 if (daq_TDClocked_crates[i] != NULL__null) {
936 daq_TDClocked_crates[i]->SetMinimum(daq_TDClocked_crates[i]->GetMinimum(0.001));
937 }
938 if (daq_TDCovr_crates[i] != NULL__null) {
939 daq_TDCovr_crates[i]->SetMinimum(daq_TDCovr_crates[i]->GetMinimum(0.001));
940 }
941 }
942
943 japp->RootFillUnLock(this); //RELEASE ROOT FILL LOCK
944
945 return NOERROR;
946}
947
948//------------------
949// fini
950//------------------
951jerror_t JEventProcessor_DAQ_online::fini(void)
952{
953 // Called before program exit after event processing is finished.
954 return NOERROR;
955}
956