Functions
test_TpcToolBasedRawDigitPrepService.cxx File Reference
#include <string>
#include <iostream>
#include <sstream>
#include <fstream>
#include <iomanip>
#include "art/Framework/Services/Registry/ServiceHandle.h"
#include "lardata/DetectorInfoServices/DetectorClocksService.h"
#include "lardataobj/RawData/RawDigit.h"
#include "dunecore/ArtSupport/ArtServiceHelper.h"
#include "dunecore/ArtSupport/DuneToolManager.h"
#include "dunecore/DuneInterface/Data/AdcTypes.h"
#include "dunecore/DuneInterface/Service/RawDigitPrepService.h"
#include "dunecore/DuneInterface/Data/WiredAdcChannelDataMap.h"
#include <cassert>
#include "TH1.h"

Go to the source code of this file.

Functions

bool sigequal (AdcSignal sig1, AdcSignal sig2)
 
bool flagequal (AdcFlag flg1, AdcFlag flg2)
 
int test_TpcToolBasedRawDigitPrepService (bool useExistingFcl=false)
 
int main (int argc, char *argv[])
 

Function Documentation

bool flagequal ( AdcFlag  flg1,
AdcFlag  flg2 
)

Definition at line 52 of file test_TpcToolBasedRawDigitPrepService.cxx.

52  {
53  if ( flg1 == flg2 ) return true;
54  if ( flg1 == AdcInterpolated && flg2 == AdcStuckOff ) return true;
55  if ( flg1 == AdcInterpolated && flg2 == AdcStuckOn ) return true;
56  if ( flg1 == AdcSetFixed ) return true;
57  cout << "flagequal: " << flg1 << " != " << flg2 << endl;
58  return false;
59 }
const AdcFlag AdcSetFixed
Definition: AdcTypes.h:41
const AdcFlag AdcStuckOn
Definition: AdcTypes.h:37
const AdcFlag AdcStuckOff
Definition: AdcTypes.h:36
const AdcFlag AdcInterpolated
Definition: AdcTypes.h:42
QTextStream & endl(QTextStream &s)
int main ( int  argc,
char *  argv[] 
)

Definition at line 320 of file test_TpcToolBasedRawDigitPrepService.cxx.

320  {
321  bool useExistingFcl = false;
322  if ( argc > 1 ) {
323  string sarg(argv[1]);
324  if ( sarg == "-h" ) {
325  cout << "Usage: " << argv[0] << " [UseExisting]" << endl;
326  cout << " If UseExisting = true, existing FCL file is used [false]." << endl;
327  return 0;
328  }
329  useExistingFcl = sarg == "true" || sarg == "1";
330  }
331  TH1::AddDirectory(false);
332  return test_TpcToolBasedRawDigitPrepService(useExistingFcl);
333 }
int test_TpcToolBasedRawDigitPrepService(bool useExistingFcl=false)
QTextStream & endl(QTextStream &s)
bool sigequal ( AdcSignal  sig1,
AdcSignal  sig2 
)

Definition at line 43 of file test_TpcToolBasedRawDigitPrepService.cxx.

43  {
44  AdcSignal sigdiff = sig2 - sig1;
45  if ( sigdiff < -0.5 || sigdiff > 0.5 ) {
46  cout << "sigequal: " << sig1 << " != " << sig2 << endl;
47  return false;
48  }
49  return true;
50 }
float AdcSignal
Definition: AdcTypes.h:21
QTextStream & endl(QTextStream &s)
int test_TpcToolBasedRawDigitPrepService ( bool  useExistingFcl = false)

Definition at line 63 of file test_TpcToolBasedRawDigitPrepService.cxx.

63  {
64  const string myname = "test_TpcToolBasedRawDigitPrepService: ";
65 #ifdef NDEBUG
66  cout << myname << "NDEBUG must be off." << endl;
67  abort();
68 #endif
69  string line = "-----------------------------";
70 
71  cout << myname << line << endl;
72  cout << myname << "Create top-level FCL." << endl;
73  std::string const fclfile{"test_TpcToolBasedRawDigitPrepService.fcl"};
74  if (!useExistingFcl) {
75  std::ofstream fout{fclfile};
76  fout << "#include \"services_dune.fcl\"" << endl;
77  fout << "services: @local::dune35t_services_legacy" << endl; // Need geometry for wire building.
78  fout << "#include \"tools_dune.fcl\"" << endl;
79  fout << "services.AdcWireBuildingService: {" << endl;
80  fout << " service_provider: StandardAdcWireBuildingService" << endl;
81  fout << " LogLevel: 1" << endl;
82  fout << "}" << endl;
83  fout << "services.RawDigitPrepService: {" << endl;
84  fout << " service_provider: TpcToolBasedRawDigitPrepService" << endl;
85  fout << " LogLevel: 3" << endl;
86  fout << " DoWires: true" << endl;
87  fout << " ToolNames: [" << endl;
88  fout << " \"digitReader\"," << endl;
89  // fout << " \"adcChannelDumper\"," << endl;
90  fout << " \"rawAdcPlotter\"," << endl;
91  fout << " \"adcSampleFiller\"," << endl;
92  fout << " \"preparedAdcPlotter\"," << endl;
93  fout << " \"adcThresholdSignalFinder\"" << endl;
94  // fout << ", \"adcRoiViewer\"" << endl;
95  fout << " ]" << endl;
96  fout << " CallgrindToolNames: []" << endl;
97  fout << "}" << endl;
98  fout.close();
99  }
100  assert( DuneToolManager::fclFilename(fclfile) == fclfile );
101  std::ifstream config{fclfile};
103 
104  cout << myname << line << endl;
105  cout << myname << "Fetching tool manager." << endl;
107  assert ( ptm != nullptr );
108  DuneToolManager& tm = *ptm;
109  tm.print();
110  assert( tm.toolNames().size() > 10 );
111 
112  AdcCount lowbits = 63;
113  AdcCount highbits = 4095 - 63;
114 
115  cout << myname << line << endl;
116  cout << myname << "Create raw digits." << endl;
117  AdcChannel nchan = 8;
118  unsigned int nsig = 64;
119  AdcSignalVectorVector sigsin(nchan);
120  float fac = 25.0;
121  unsigned int isig_stucklo = 15;
122  unsigned int isig_stuckhi = 25;
123  bool doSticky = false;
124  // nchan is 8
125  AdcSignal peds[8] = {2000.2, 2010.1, 2020.3, 1990.4, 1979.6, 1979.2, 1995.0, 2001.3};
126  vector<RawDigit> digs;
127  map<AdcChannel, AdcCountVector> adcsmap;
128  map<AdcChannel, AdcFlagVector> expflagsmap;
129  // Each channel is filled with a bipolar signal offset by one from preceding channel.
130  for ( AdcChannel chan=0; chan<nchan; ++chan ) {
131  unsigned int isig1 = 10 + chan;
132  for ( unsigned int isig=0; isig<isig1; ++isig ) sigsin[chan].push_back(0);
133  for ( unsigned int i=0; i<10; ++i ) sigsin[chan].push_back(fac*i);
134  for ( unsigned int i=10; i<1000; --i ) sigsin[chan].push_back(fac*i);
135  for ( unsigned int i=19; i<1000; --i ) sigsin[chan].push_back(-sigsin[chan][i+isig1]);
136  for ( unsigned int isig=sigsin[chan].size();
137  isig<nsig; ++isig ) sigsin[chan].push_back(0);
138  assert(sigsin[chan].size() == nsig);
139  AdcCountVector adcsin;
140  for ( unsigned int isig=0; isig<nsig; ++isig) {
141  AdcSignal sig = sigsin[chan][isig] + peds[chan];
142  AdcCount adc = 0;
143  if ( sig > 0.0 ) adc = int(sig+0.5);
144  if ( adc > 4095 ) adc = 4095;
145  AdcCount adchigh = adc & highbits;
146  if ( doSticky ) {
147  if ( isig == isig_stucklo ) adc = adchigh; // Stuck low bits to zero.
148  if ( isig == isig_stuckhi ) adc = adchigh + lowbits; // Stuck low bits to one.
149  }
150  adcsin.push_back(adc);
151  }
152  assert(adcsin.size() == nsig);
153  adcsmap[chan] = adcsin;
154  // Create RawDigit.
155  RawDigit& dig = *digs.emplace(digs.end(), chan, nsig, adcsin, raw::kNone);
156  dig.SetPedestal(peds[chan]);
157  cout << myname << " Compressed size: " << dig.NADC() << endl;
158  cout << myname << " Uncompressed size: " << dig.Samples() << endl;
159  cout << myname << " Pedestal: " << dig.GetPedestal() << endl;
160  cout << myname << " Channel: " << dig.Channel() << endl;
161  assert(dig.Samples() == nsig);
162  assert(dig.Channel() == chan);
163  assert(dig.GetPedestal() == peds[chan]);
164  }
165  assert( adcsmap.size() == nchan );
166 
167  cout << myname << line << endl;
168  cout << myname << "Create the expected flag vector." << endl;
169  for ( AdcChannel chan=0; chan<nchan; ++chan ) {
170  AdcFlagVector expflags(nsig, AdcGood);
171  for ( unsigned int isig=0; isig<nsig; ++isig) {
172  AdcCount adc = adcsmap[chan][isig];
173  AdcCount adclow = adc & lowbits;
174  if ( adc <= 0 ) expflags[isig] = AdcUnderflow;
175  else if ( adc >= 4095 ) expflags[isig] = AdcOverflow;
176  else if ( doSticky ) {
177  if ( adclow == 0 ) expflags[isig] = AdcStuckOff;
178  else if ( adclow == lowbits ) expflags[isig] = AdcStuckOn;
179  }
180  }
181  expflagsmap[chan] = expflags;
182  }
183 
184  cout << myname << line << endl;
185  cout << myname << "Fetch raw digit prep service." << endl;
187  hrdp->print();
188 
189  cout << myname << line << endl;
190  cout << myname << "Prep data from digits." << endl;
191  AdcSignalVector sigs;
192  AdcFlagVector flags;
193  AdcChannelDataMap prepdigs;
194  for ( unsigned int idig=0; idig<digs.size(); ++idig ) {
195  const RawDigit& dig = digs[idig];
196  assert( prepdigs.find(dig.Channel()) == prepdigs.end() );
197  AdcChannelData& data = prepdigs[dig.Channel()];
198  data.setChannelInfo(dig.Channel());
199  data.digitIndex = idig;
200  data.digit = &dig;
201  data.setEventInfo(123, 1);
202  }
203  std::vector<recob::Wire> wires;
204  wires.reserve(nchan);
205  vector<string> snames;
206  WiredAdcChannelDataMap intStates(snames, nchan);
207  assert( intStates.dataMaps.size() == snames.size() );
208  assert( intStates.wires.size() == snames.size() );
209  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataForJob();
210  assert( hrdp->prepare(clockData, prepdigs, &wires, &intStates) == 0 );
211  cout << myname << " # prepared digit channels: " << prepdigs.size() << endl;
212  cout << myname << " # wire channels: " << wires.size() << endl;
213  cout << myname << " # intermediate state channels: " << intStates.dataMaps.size() << endl;
214  for ( const auto& namedadm : intStates.dataMaps ) {
215  string sname = namedadm.first;
216  const AdcChannelDataMap& adm = namedadm.second;
217  auto iwco = intStates.wires.find(sname);
218  const vector<Wire>* pwires = nullptr;
219  if ( iwco == intStates.wires.end() ) {
220  cout << myname << " Wires not found for intermediate state " << sname << "." << endl;
221  assert( iwco != intStates.wires.end() );
222  } else {
223  pwires = iwco->second;
224  }
225  assert( pwires != nullptr );
226  cout << myname << " State " << sname << " has " << adm.size() << " ADC channels";
227  if ( pwires != nullptr ) cout << " and " << pwires->size() << " wires";
228  cout <<"." << endl;
229  assert( pwires->size() == adm.size() );
230  }
231  cout << myname << " # intermediate wires: " << intStates.wires.size() << endl;
232  for ( AdcChannelDataMap::const_iterator ichdat=prepdigs.begin(); ichdat!=prepdigs.end(); ++ichdat ) {
233  AdcChannel chan = ichdat->first;
234  const AdcChannelData& acd = ichdat->second;
235  const AdcSignalVector& sigs = acd.samples;
236  const AdcFlagVector& flags = acd.flags;
237  const raw::RawDigit* pdig = acd.digit;
238  const Wire* pwire = acd.wire;
239  AdcSignal ped = acd.pedestal;
240  cout << myname << "----- Channel " << chan << endl;
241  cout << myname << " Final signal tick count: " << sigs.size() << endl;
242  cout << myname << " Final flag tick count: " << flags.size() << endl;
243  cout << myname << " Final flag tick count: " << flags.size() << endl;
244  cout << myname << " Pedestal: " << ped << endl;
245  cout << myname << " samples[0]: " << sigs[0] << endl;
246  cout << myname << "Check final data." << endl;
247  assert( pwire != nullptr );
248  assert( pwire->SignalROI().size() > 0 );
249  assert( sigs.size() == nsig );
250  assert( flags.size() == nsig );
251  assert( pdig != nullptr );
252  assert( pdig == &digs[chan] );
253  assert( pdig->Channel() == chan );
254  assert( ichdat->second.digitIndex == chan );
255  const AdcFlagVector& expflags = expflagsmap[chan];
256  assert( expflagsmap[chan].size() == nsig );
257  // Fetch intermediate data.
258  cout << myname << "Fetch intermediate data." << endl;
259  vector<const AdcSignalVector*> intSigs;
260  vector<const AdcFlagVector*> intFlags;
261  cout << myname << " ...bad" << endl;
262  auto iacd = intStates.dataMaps.find("bad");
263  assert( iacd == intStates.dataMaps.end() );
264  string header = " ch-tk raw";
265  unsigned int nintexp = 0;
266  for ( string sname : snames ) {
267  cout << myname << " ..." << sname << endl;
268  iacd = intStates.dataMaps.find(sname);
269  assert( iacd != intStates.dataMaps.end() );
270  const AdcChannelData& intAcd = iacd->second[chan];
271  intSigs.push_back(&intAcd.samples);
272  intFlags.push_back(&intAcd.flags);
273  assert( intAcd.wire != nullptr );
274  for ( unsigned int i=sname.size(); i<12; ++i ) header += " ";
275  header += sname;
276  ++nintexp;
277  cout << myname << " Checking sample and flag tick counts." << endl;
278  assert( intSigs.back() != nullptr );
279  assert( intSigs.back()->size() != 0 );
280  assert( intSigs.back()->size() == nsig );
281  assert( intFlags.back() != nullptr );
282  assert( intFlags.back()->size() != 0 );
283  assert( intFlags.back()->size() == nsig );
284  cout << myname << " Wire ROI count: " << intAcd.wire->SignalROI().n_ranges() << endl;
285  cout << myname << " Wire Tick count: " << intAcd.wire->SignalROI().size() << endl;
286  assert( intAcd.wire->SignalROI().n_ranges() == 1 );
287  assert( intAcd.wire->SignalROI().size() == nsig );
288  }
289  header += " final";
290  assert( intStates.dataMaps.size() == nintexp );
291  //assert( intStates.wires.size() == nintexp );
292  assert( intSigs.size() == nintexp );
293  assert( intFlags.size() == nintexp );
294  // Display results.
295  cout << myname << "Display intermediate and final samples." << endl;
296  cout << myname << header << endl;
297  for ( unsigned int isig=0; isig<nsig; ++isig ) {
298  cout << myname;
299  cout << setw(4) << chan << "-"
300  << setw(2) << isig << ": " << setw(4) << adcsmap[chan][isig];
301  for ( unsigned int ista=0; ista<intSigs.size(); ++ista ) {
302  cout << fixed << setprecision(1) << setw(8) << intSigs[ista]->at(isig);
303  cout << " [" << intFlags[ista]->at(isig) << "]";
304  }
305  cout << fixed << setprecision(1) << setw(8) << sigs[isig]
306  << " [" << flags[isig] << "]" << endl;
307  assert( adcsmap[chan][isig] == acd.raw[isig] );
308  if ( flags[isig] == AdcGood ) assert( sigequal(sigs[isig], sigsin[chan][isig]) );
309  assert( flagequal(flags[isig], expflags[isig]) );
310  }
311  }
312 
313  cout << myname << line << endl;
314  cout << "Done." << endl;
315  return 0;
316 }
float GetPedestal() const
Definition: RawDigit.h:214
bool flagequal(AdcFlag flg1, AdcFlag flg2)
std::vector< AdcCount > AdcCountVector
Definition: AdcTypes.h:19
ULong64_t Samples() const
Number of samples in the uncompressed ADC data.
Definition: RawDigit.h:213
Collection of charge vs time digitized from a single readout channel.
Definition: RawDigit.h:69
size_type size() const
Returns the size of the vector.
const std::vector< std::string > & toolNames() const
std::string string
Definition: nybbler.cc:12
float AdcSignal
Definition: AdcTypes.h:21
std::vector< AdcFlag > AdcFlagVector
Definition: AdcTypes.h:30
ChannelID_t Channel() const
DAQ channel this raw data was read from.
Definition: RawDigit.h:212
int16_t adc
Definition: CRTFragment.hh:202
const AdcFlag AdcUnderflow
Definition: AdcTypes.h:33
intermediate_table::const_iterator const_iterator
void print() const
const recob::Wire * wire
const raw::RawDigit * digit
virtual int prepare(detinfo::DetectorClocksData const &clockData, AdcChannelDataMap &prepdigs, std::vector< recob::Wire > *pwires=nullptr, WiredAdcChannelDataMap *pwiredData=nullptr) const =0
no compression
Definition: RawTypes.h:9
static void load_services(std::string const &config)
tm
Definition: demo.py:21
Q_EXPORT QTSManip setprecision(int p)
Definition: qtextstream.h:343
const AdcFlag AdcGood
Definition: AdcTypes.h:32
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:92
fInnerVessel push_back(Point(-578.400000, 0.000000, 0.000000))
const AdcFlag AdcOverflow
Definition: AdcTypes.h:34
static Config * config
Definition: config.cpp:1054
size_t NADC() const
Number of elements in the compressed ADC sample vector.
Definition: RawDigit.h:207
const AdcFlag AdcStuckOn
Definition: AdcTypes.h:37
const RegionsOfInterest_t & SignalROI() const
Returns the list of regions of interest.
Definition: Wire.h:228
static std::string fclFilename(std::string setName="", int dbg=1)
AdcCountVector raw
Q_EXPORT QTSManip setw(int w)
Definition: qtextstream.h:331
const AdcFlag AdcStuckOff
Definition: AdcTypes.h:36
AdcSignal pedestal
std::vector< AdcSignalVector > AdcSignalVectorVector
Definition: AdcTypes.h:23
unsigned int AdcChannel
Definition: AdcTypes.h:50
void SetPedestal(float ped, float sigma=1.)
Set pedestal and its RMS (the latter is 0 by default)
Definition: RawDigit.cxx:68
void line(double t, double *p, double &x, double &y, double &z)
Class holding the regions of interest of signal from a channel.
Definition: Wire.h:118
std::vector< AdcSignal > AdcSignalVector
Definition: AdcTypes.h:22
virtual std::ostream & print(std::ostream &out=std::cout, std::string prefix="") const =0
std::map< AdcChannel, AdcChannelData > AdcChannelDataMap
short AdcCount
Definition: AdcTypes.h:18
bool sigequal(AdcSignal sig1, AdcSignal sig2)
static DuneToolManager * instance(std::string fclname="", int dbg=1)
AdcSignalVector samples
QTextStream & endl(QTextStream &s)
AdcFlagVector flags