10#include <ecl/modules/eclDQMExtended/eclDQMExtended.h>
13#include <ecl/dataobjects/ECLElementNumbers.h>
14#include <ecl/utility/ECLDspEmulator.h>
15#include <ecl/utility/ECLDspUtilities.h>
18#include <framework/core/HistoModule.h>
19#include <framework/logging/Logger.h>
22#include <TDirectory.h>
27#include <boost/format.hpp>
67 "histogram directory in ROOT file", std::string(
"ECL"));
69 "How to initialize DSP coeffs: ''DB'' or ''File'' are acceptable ", std::string(
"DB"));
71 "directory for DSP coeffs", std::string(
"/hsm/belle2/bdata/users/dmitry/dsp/"));
73 "Name of run with DSP files", std::string(
"run0000/"));
75 "(ampdiff_{cellid,shaper}, timediff_{cellid,shaper} histograms) for "
76 "failed fits.",
false);
78 "time determination in ShaperDSP emulator",
true);
80 std::vector<int> default_skip = {
87 "source (from TRGSummary)", default_skip);
100 TDirectory* oldDir = gDirectory;
110 std::vector<TH1F*> h_amp_qualityfail_raw, h_time_qualityfail_raw;
112 for (
int i = 0; i < 4; i++) {
113 std::string h_name_a, h_title_a, h_name_t, h_title_t;
114 h_name_a = str(boost::format(
"amp_timefail_q%1%") % i);
115 h_title_a = str(boost::format(
"Amp for time mismatches w/ FPGA fit qual=%1%") % i);
116 h_name_t = str(boost::format(
"time_ampfail_q%1%") % i);
117 h_title_t = str(boost::format(
"Time for amp mismatches w/ FPGA fit qual=%1%") % i);
118 TH1F* h_a =
new TH1F(h_name_a.c_str(), h_title_a.c_str(), 1200, 0, 262144);
119 TH1F* h_t =
new TH1F(h_name_t.c_str(), h_title_t.c_str(), 240, -2050, 2050);
124 for (
int i = 0; i < 4; i++) {
125 for (
int j = 0; j < 4; j++) {
127 std::string h_name_a, h_title_a, h_name_t, h_title_t;
128 h_name_a = str(boost::format(
"amp_qf%1%_qd%2%") % i % j);
129 h_title_a = str(boost::format(
"Amp for C++ fit qual=%1% and FPGA fit qual=%2%") % i % j);
130 h_name_t = str(boost::format(
"time_qf%1%_qd%2%") % i % j);
131 h_title_t = str(boost::format(
"Time for C++ fit qual=%1% and FPGA fit qual=%2%") % i % j);
132 TH1F* h_a =
new TH1F(h_name_a.c_str(), h_title_a.c_str(), 1200, 0, 262144);
133 TH1F* h_t =
new TH1F(h_name_t.c_str(), h_title_t.c_str(), 240, -2050, 2050);
134 h_amp_qualityfail_raw.push_back(h_a);
135 h_time_qualityfail_raw.push_back(h_t);
140 h_amp_qualityfail_raw.clear();
141 h_time_qualityfail_raw.clear();
144 h_ampfail_quality =
new TH1F(
"ampfail_quality",
"Number of FPGA <-> C++ fitter #bf{amp} inconsistencies vs fit qual", 5, -1, 4);
146 h_ampfail_quality->GetXaxis()->SetTitle(
"FPGA fit qual. -1-all evts,0-good,1-int overflow,2-low amp,3-bad chi2");
149 h_timefail_quality =
new TH1F(
"timefail_quality",
"Number of FPGA <-> C++ fitter #bf{time} inconsistencies vs fit qual", 5, -1, 4);
151 h_timefail_quality->GetXaxis()->SetTitle(
"FPGA fit qual. -1-all evts,0-good,1-int overflow,2-low amp,3-bad chi2");
153 h_ampfail_cellid =
new TH1F(
"ampfail_cellid",
"Cell IDs w/ amp inconsistencies",
157 h_timefail_cellid =
new TH1F(
"timefail_cellid",
"Cell IDs w/ time inconsistencies",
161 h_amptimefail_cellid =
new TH1F(
"amptimefail_cellid",
"Cell IDs w/ time and amp inconsistencies",
165 h_ampfail_shaperid =
new TH1F(
"ampfail_shaperid",
"Shaper IDs w/ amp inconsistencies", 624, 1, 625);
168 h_timefail_shaperid =
new TH1F(
"timefail_shaperid",
"Shaper IDs w/ time inconsistencies", 624, 1, 625);
171 h_amptimefail_shaperid =
new TH1F(
"amptimefail_shaperid",
"Shaper IDs w/ time and amp inconsistencies", 624, 1, 625);
174 h_ampfail_crateid =
new TH1F(
"ampfail_crateid",
"Crate IDs w/ amp inconsistencies", 52, 1, 53);
177 h_timefail_crateid =
new TH1F(
"timefail_crateid",
"Crate IDs w/ time inconsistencies", 52, 1, 53);
180 h_amptimefail_crateid =
new TH1F(
"amptimefail_crateid",
"Crate IDs w/ time and amp inconsistencies", 52, 1, 53);
187 h_qualityfail_shaperid =
new TH1F(
"qualityfail_shaperid",
"Shaper IDs w/ fit qual inconsistencies", 624, 1, 625);
190 h_qualityfail_crateid =
new TH1F(
"qualityfail_crateid",
"Crate IDs w/ fit qual inconsistencies", 52, 1, 53);
193 h_fail_shaperid =
new TH1F(
"fail_shaperid",
"Shaper IDs w/ inconsistencies", 624, 1, 625);
196 h_fail_crateid =
new TH1F(
"fail_crateid",
"Crate IDs w/ inconsistencies", 52, 1, 53);
197 h_fail_crateid->GetXaxis()->SetTitle(
"Crate ID (same as ECLCollector ID)");
203 h_ampdiff_cellid =
new TH2F(
"ampdiff_cellid",
"Amp. diff. (Emulator-Data) for amp inconsistencies",
208 h_timediff_cellid =
new TH2F(
"timediff_cellid",
"Time diff. (Emulator-Data) for time inconsistencies",
214 "for amp inconsistencies vs Shaper Id",
215 624, 1, 625, 239, -262143, 262143);
220 "for time inconsistencies vs Shaper Id",
221 624, 1, 625, 239, -4095, 4095);
226 h_ampdiff_quality =
new TH2F(
"ampdiff_quality",
"Amp. diff. (Emulator-Data) for amp. inconsistencies", 4, 0, 4, 239,
228 h_ampdiff_quality->GetXaxis()->SetTitle(
"FPGA fit quality. 0-good, 1-int overflow, 2-low amp, 3-bad chi2");
231 h_timediff_quality =
new TH2F(
"timediff_quality",
"Time diff. (Emulator-Data) for time inconsistencies", 4, 0, 4, 239,
233 h_timediff_quality->GetXaxis()->SetTitle(
"FPGA fit quality. 0-good, 1-int overflow, 2-low amp, 3-bad chi2");
236 h_quality_fit_data =
new TH2F(
"quality_fit_data",
"C++ fitter vs FPGA, fit quality inconsistencies", 4, 0, 4, 4, 0, 4);
237 h_quality_fit_data->GetXaxis()->SetTitle(
"C++ fit qual. 0-good,1-int overflow,2-low amp,3-bad chi2");
238 h_quality_fit_data->GetYaxis()->SetTitle(
"FPGA fit qual. 0-good,1-int overflow,2-low amp,3-bad chi2");
240 h_ampflag_qualityfail =
new TH2F(
"ampflag_qualityfail",
"Amp flag (0/1) for fit qual inconsistencies", 4, 0, 4, 4, -1,
242 h_ampflag_qualityfail->GetXaxis()->SetTitle(
"FPGA fit quality. 0-good,1-int overflow,2-low amp,3-bad chi2");
245 h_timeflag_qualityfail =
new TH2F(
"timeflag_qualityfail",
"Time flag (0/1) for fit qual inconsistencies", 4, 0, 4, 4,
247 h_timeflag_qualityfail->GetXaxis()->SetTitle(
"FPGA fit quality. 0-good,1-int overflow,2-low amp,3-bad chi2");
251 h_missing_ecldigits->SetTitle(
"Fit quality flag for missing ECLDigits (missing fit results)");
252 h_missing_ecldigits->GetXaxis()->SetTitle(
"C++ fit quality. 0-good,1-int overflow,2-low amplitude,3-bad chi2");
266 if (!
mapper.initFromDB()) B2FATAL(
"ECL DQM logic test FATAL:: Can't initialize eclChannelMapper");
270 else B2FATAL(
"ECL DQM logic test FATAL: No way to initialize DSP coeffs!!! Please choose InitKey = DB or InitKey = File");
275 const std::vector<float> intermediate = cal->getCalibVector();
276 constants.resize(intermediate.size());
277 for (
size_t i = 0; i < constants.size(); i++) constants[i] = (
short int)intermediate[i];
282 std::map<std::string, std::vector<short int>>& map1,
283 std::map<std::string, short int>& map2)
286 dspdata->
getF(map1[
"F"]);
287 dspdata->
getF1(map1[
"F1"]);
288 dspdata->
getF31(map1[
"F31"]);
289 dspdata->
getF32(map1[
"F32"]);
290 dspdata->
getF33(map1[
"F33"]);
291 dspdata->
getF41(map1[
"F41"]);
292 dspdata->
getF43(map1[
"F43"]);
294 map2[
"k_a"] = (
short int)dspdata->
getka();
295 map2[
"k_b"] = (
short int)dspdata->
getkb();
296 map2[
"k_c"] = (
short int)dspdata->
getkc();
298 map2[
"k_1"] = (
short int)dspdata->
getk1();
299 map2[
"k_2"] = (
short int)dspdata->
getk2();
306 int iCrate =
mapper.getCrateID(cellID);
307 int iShaperPosition =
mapper.getShaperPosition(cellID);
308 return (iCrate - 1) * 12 + iShaperPosition;
313 size_t size = vectorFrom.size();
314 if (size % 16) B2ERROR(
"ECL DQM logic test error: Split is impossible!" <<
LogVar(
"Vector size", size));
315 return (vectorFrom.data() + (size / 16) * (iChannel - 1));
342 if (iShaper - (iShaper - 1) / 12 * 12 > 10)
continue;
344 if (iShaper - (iShaper - 1) / 12 * 12 > 8)
continue;
354 const std::filesystem::path RunSubDir(
m_RunName);
355 const std::regex
Filter(
".*(crate)([0-9]{2})/.*(dsp)([0-9]{2})(.dat)");
356 if (!exists(MainDir / RunSubDir)) B2FATAL(
"ECL DQM logic test FATAL: Directory w/ DSP files don't exist" <<
LogVar(
"Directory",
357 MainDir / RunSubDir));
358 for (
const std::filesystem::directory_entry& x : std::filesystem::recursive_directory_iterator(MainDir / RunSubDir)) {
359 if (!std::regex_match(x.path().string(),
Filter) || !std::filesystem::is_regular_file(x.path()))
continue;
360 int iCrate = atoi(std::regex_replace(x.path().string(),
Filter,
"$2").c_str());
361 int iShaperPosition = atoi(std::regex_replace(x.path().string(),
Filter,
"$4").c_str());
362 int iShaper = (iCrate - 1) * 12 + iShaperPosition;
363 if (iCrate > 36 && iCrate < 45) {
364 if (iShaperPosition > 10)
continue;
365 }
else if (iCrate > 44) {
366 if (iShaperPosition > 8)
continue;
380 int iChannelPosition =
mapper.getShaperChannel(cellID);
383 const short int* f =
vectorsplit(map_vec.at(
"F"), iChannelPosition);
384 const short int* f1 =
vectorsplit(map_vec.at(
"F1"), iChannelPosition);
385 const short int* fg31 =
vectorsplit(map_vec.at(
"F31"), iChannelPosition);
386 const short int* fg32 =
vectorsplit(map_vec.at(
"F32"), iChannelPosition);
387 const short int* fg33 =
vectorsplit(map_vec.at(
"F33"), iChannelPosition);
388 const short int* fg41 =
vectorsplit(map_vec.at(
"F41"), iChannelPosition);
389 const short int* fg43 =
vectorsplit(map_vec.at(
"F43"), iChannelPosition);
392 int k_a = map_coef.at(
"k_a");
393 int k_b = map_coef.at(
"k_b");
394 int k_c = map_coef.at(
"k_c");
395 int k_1 = map_coef.at(
"k_1");
396 int k_2 = map_coef.at(
"k_2");
397 int k_16 = map_coef.at(
"k_16");
398 int chi_thres = map_coef.at(
"chi_thres");
404 int* y = adc_data.data();
405 int ttrig2 = trigger_time - 2 * (trigger_time / 8);
407 auto result = lftda_(f, f1, fg41, fg43, fg31, fg32, fg33, y, ttrig2, A0,
408 Ahard, Askip, k_a, k_b, k_c, k_16, k_1, k_2, chi_thres,
414 if (result.skip_thr || result.hit_thr)
m_QualityFit += 4;
419 for (
int i = 0; i < 4; i++) {
423 for (
int i = 0; i < 4; i++) {
424 for (
int j = 0; j < 3; j++) {
469 if (timing_type == skipped_timing_type)
return;
473 int iAmpflag_qualityfail = 0;
474 int iTimeflag_qualityfail = 0;
476 if (!
m_ECLTrigs.isValid()) B2FATAL(
"ECL DQM logic test FATAL: Trigger time information is not available");
480 std::vector<int> DspArray = aECLDsp.getDspA();
488 B2ERROR(
"ECL DQM logic test error: ECL Digit does not exist for A_emulator > Thr_skip"
532 for (
int i = 0; i < 4; i++) {
533 for (
int j = 0; j < 3; j++) {
535 if (i == 0) k = j + 1;
536 if (i == 1) {
if (j == 0) k = j;
else k = j + 1; }
537 if (i == 2) {
if (j == 2) k = j + 1;
else k = j; }
Class for accessing objects in the database.
TH2F * h_timediff_cellid
Histogram: Time diff.
virtual ~ECLDQMEXTENDEDModule() override
Destructor.
TH1F * h_amptimefail_shaperid
Histogram: ShaperIDs w/ failed amplitudes and times.
int m_TimeData
Signal time from ECL data.
void initDspfromDB()
Get DSP coeffs and auxiliary constants from DB.
TH2F * h_timediff_shaperid
Histogram: Time diff.
TH1F * h_ampfail_crateid
Histogram: CrateIDs w/ failed amplitudes.
TH1F * h_fail_crateid
Histogram: CrateIDs w/ failed logic.
bool m_SaveDetailedFitData
Save detailed fit data for failed fits.
std::vector< short int > v_totalthrA0
Vector to store Low amplitude thresholds.
TH1F * h_timefail_quality
Histogram: Time control flags in bins of QualityData.
bool m_adjusted_timing
Use modified time determination algorithm in emulator, same as in ShaperDSP version >= 1....
int m_AmpFit
Signal amplitude obtaining from DSP emulator.
int m_TrigTime
Trigger time value.
int m_CellId
Cell ID number.
StoreArray< ECLDsp > m_ECLDsps
ECL DSP data.
virtual void initialize() override
Initialize the module.
TH1F * h_timefail_shaperid
Histogram: ShaperIDs w/ failed times.
std::vector< short int > v_totalthrAhard
Vector to store hit thresholds.
int m_QualityData
Quality flag from ECL data.
TH1F * h_qualityfail_shaperid
Histogram: ShaperIDs w/ failed qualities.
virtual void event() override
Event processor.
TH1F * h_fail_shaperid
Histogram: ShaperIDs w/ failed logic.
void emulator(int, int, std::vector< int >)
Call for DSP emulator.
TH2F * h_timediff_quality
Histogram: Time diff.
virtual void endRun() override
Call when a run ends.
TH1F * h_qualityfail_crateid
Histogram: CrateIDs w/ failed qualities.
TH1F * h_amptimefail_crateid
Histogram: CrateIDs w/ failed amplitudes and times.
virtual void terminate() override
Terminate.
int m_AmpData
Signal amplitude from ECL data.
void callbackCalibration(DBObjPtr< ECLCrystalCalib > &, std::vector< short int > &)
Read calibration values for thresholds from DBObject.
std::vector< short int > v_totalthrAskip
Vector to store skip amplitude threshold.
std::map< int, std::map< std::string, std::vector< short int > > > map_container_vec
Map to store DSP coeffs.
std::vector< int > m_skipEvents
Skip events with specified type of timing source (see TRGSummary class, ETimingType enum)
StoreArray< ECLDigit > m_ECLDigits
ECL digits.
DBArray< ECLDspData > m_ECLDspDataArray2
DBArray for payload 'ECLDSPPars2'.
std::vector< TH1F * > h_amp_timefail
Histogram vector: Amplitude for time mismacthes (AmplitudeFit == AmplitudeData) w/ various QualityDat...
TH1F * h_amptimefail_cellid
Histogram: CellIDs w/ failed amplitudes and times.
TH2F * h_ampflag_qualityfail
Histogram: Amp flag (0/1) w/ failed qualities in bins of QualityData.
DBObjPtr< ECLCrystalCalib > m_calibrationThrAhard
Hit threshold.
TH1F * h_timefail_crateid
Histogram: CrateIDs w/ failed times.
const short int * vectorsplit(const std::vector< short int > &, int)
Select from vector of DSP coeffs a subvector corresponding to accurate channel number.
std::string m_InitKey
Key to initialize DSP coeffs.
TH2F * h_quality_fit_data
Histogram: QualityFit vs QualityData for quality fails.
DBObjPtr< ECLCrystalCalib > m_calibrationThrAskip
Skip amplitude threshold.
int conversion(int)
Convert a CellID number to the global Shaper number.
virtual void beginRun() override
Call when a run begins.
std::string m_histogramDirectoryName
Histogram directory in ROOT file.
ECL::ECLChannelMapper mapper
ECL channel mapper.
ECLDQMEXTENDEDModule()
< derived from HistoModule class.
TH2F * h_ampdiff_quality
Histogram: Amp.
DBArray< ECLDspData > m_ECLDspDataArray1
DBArray for payload 'ECLDSPPars1'.
TH1F * h_ampfail_cellid
Histogram: CellIDs w/ failed amplitudes.
TH2F * h_ampdiff_shaperid
Histogram: Amp.
TH1F * h_ampfail_quality
Histogram: Amp.
int m_QualityFit
Quality flag obtaining from DSP emulator.
std::string m_RunName
Run with valid DSP coeffs.
DBArray< ECLDspData > m_ECLDspDataArray0
DBArray for payload 'ECLDSPPars0'.
TH2F * h_ampdiff_cellid
Histogram: Amplitude diff.
TH1F * h_ampfail_shaperid
Histogram: ShaperIDs w/ failed amplitudes.
TH1F * h_timefail_cellid
Histogram: CellIDs w/ failed times.
StoreArray< ECLTrig > m_ECLTrigs
ECL trigger data.
DBObjPtr< ECLCrystalCalib > m_calibrationThrA0
Low amplitude threshold.
TH1F * h_qualityfail_cellid
Histogram: CellIDs w/ failed qualities.
std::map< int, std::map< std::string, short int > > map_container_coef
Map to store auxiliary constants for all shapers.
std::vector< std::vector< TH1F * > > h_amp_qualityfail
Histogram vector: AmplitudeData for quality mismathes w/ various QualityFit (raw) and QualityData (co...
StoreObjPtr< TRGSummary > m_TRGSummary
StoreObjPtr TRGSummary.
std::vector< TH1F * > h_time_ampfail
Histogram vector: Time for amplitude mismathces (TimeFit == TimeData) w/ various QualityData values.
std::string m_DSPDirectoryName
Directory name consisting of DSP coeffs.
int m_TimeFit
Signal time obtaining from DSP emulator.
void initDspfromFile()
Get DSP coeffs and auxiliary constants from Files.
TH2F * h_timeflag_qualityfail
Histogram: Time flag (0/1) w/ failed qualities in bins of Quality Data.
TH1F * h_missing_ecldigits
Histogram: Quality flag (from emulator) for the waveforms where fit result (ECLDigit) should have bee...
std::vector< std::vector< TH1F * > > h_time_qualityfail
Histogram vector: TimeData for quality mismacthes w/ various QualitFit (raw) and QualityData (column)...
virtual void defineHisto() override
Function to define histograms.
Class to store ECL digitized hits (output of ECLDigi) relation to ECLHit filled in ecl/modules/eclDig...
int getAmp() const
Get Fitting Amplitude.
int getQuality() const
Get Fitting Quality.
static ECLDigit * getByCellID(int cid)
Find ECLDigit by Cell ID using linear search.
int getTimeFit() const
Get Fitting Time.
This object contains ECL DSP coefs – electromagnetic calorimeter digital signal processing coefficien...
unsigned char getkb() const
Number of bits for FG32.
short int getchiThresh() const
chi2 threshold for fit quality flag
void getF43(std::vector< short int > &dst) const
Alternative for FG33 for signals with small amplitude.
unsigned char gety0Startr() const
start point for pedestal calculation
void getF(std::vector< short int > &dst) const
Array with tabulated signal waveform.
void getF31(std::vector< short int > &dst) const
Array FG31, used to estimate signal amplitude.
unsigned char getk2() const
multipliers power of 2 for chi2 calculation
void getF33(std::vector< short int > &dst) const
Array FG33, used to estimate pedestal height in signal.
unsigned char getka() const
Number of bits for FG31, FG41.
unsigned char getk1() const
multipliers power of 2 for f, f1
void getF32(std::vector< short int > &dst) const
Array FG32, used to estimate A * delta_t.
void getF41(std::vector< short int > &dst) const
Alternative for FG31 for signals with small amplitude.
void getF1(std::vector< short int > &dst) const
Array with tabulated derivative of signal waveform.
unsigned char getkc() const
Number of bits for FG33, FG43.
double getTimeTrig() const
Get Trig Time.
static ECLTrig * getByCellID(int cid)
Find ECLTrig by Cell ID using linear search.
static ECLDspData * readEclDsp(const char *raw_file, int boardNumber)
Convert ECLDspData from *.dat file to Root object.
This class is used to select pairs, triplets... of objects.
HistoModule()
Constructor.
void setDescription(const std::string &description)
Sets the description of the module.
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
@ TTYP_DPHY
delayed physics events for background
@ TTYP_POIS
poisson random trigger
@ TTYP_RAND
random trigger events
Class to store variables with their name which were sent to the logging service.
void addParam(const std::string &name, T ¶mVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
const int c_NCrystals
Number of crystals.
Abstract base class for different kinds of events.