LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
DrawRawHist_tool.cc
Go to the documentation of this file.
1 
12 #include "larevt/CalibrationDBI/Interface/DetPedestalProvider.h"
13 #include "larevt/CalibrationDBI/Interface/DetPedestalService.h"
14 
16 
22 
23 #include "TH1F.h"
24 
25 namespace evdb_tool {
26 
27  class DrawRawHist : public IWaveformDrawer {
28  public:
29  explicit DrawRawHist(const fhicl::ParameterSet& pset);
30 
31  ~DrawRawHist();
32 
33  void configure(const fhicl::ParameterSet& pset) override;
34  void Fill(evdb::View2D&, raw::ChannelID_t&, float, float) override;
35  void Draw(const std::string&, float, float) override;
36 
37  float getMaximum() const override { return fMaximum; };
38  float getMinimum() const override { return fMinimum; };
39 
40  private:
41  void BookHistogram(raw::ChannelID_t&, float, float);
42 
43  float fMaximum;
44  float fMinimum;
45 
46  std::unique_ptr<TH1F> fRawDigitHist;
47  };
48 
49  //----------------------------------------------------------------------
50  // Constructor.
52  {
53  configure(pset);
54  }
55 
57 
59  {
60  return;
61  }
62 
64  raw::ChannelID_t& channel,
65  float lowBin,
66  float numTicks)
67  {
69 
70  //grab the singleton with the event
72  if (!event) return;
73 
74  // Handle histograms
75  BookHistogram(channel, lowBin, numTicks);
76 
77  fMinimum = std::numeric_limits<float>::max();
78  fMaximum = std::numeric_limits<float>::lowest();
79 
80  // Loop over the possible producers of RawDigits
81  for (const auto& rawDataLabel : rawOpt->fRawDataLabels) {
82  art::Handle<std::vector<raw::RawDigit>> rawDigitVecHandle;
83  event->getByLabel(rawDataLabel, rawDigitVecHandle);
84 
85  if (!rawDigitVecHandle.isValid()) continue;
86 
87  for (size_t rawDigitIdx = 0; rawDigitIdx < rawDigitVecHandle->size(); rawDigitIdx++) {
88  art::Ptr<raw::RawDigit> rawDigit(rawDigitVecHandle, rawDigitIdx);
89 
90  if (rawDigit->Channel() != channel) continue;
91 
92  // We will need the pedestal service...
93  const lariov::DetPedestalProvider& pedestalRetrievalAlg =
95 
96  // recover the pedestal
97  float pedestal = 0;
98 
99  if (rawOpt->fPedestalOption == 0) { pedestal = pedestalRetrievalAlg.PedMean(channel); }
100  else if (rawOpt->fPedestalOption == 1) {
101  pedestal = rawDigit->GetPedestal();
102  }
103  else if (rawOpt->fPedestalOption == 2) {
104  pedestal = 0;
105  }
106  else {
107  mf::LogWarning("DrawRawHist")
108  << " PedestalOption is not understood: " << rawOpt->fPedestalOption
109  << ". Pedestals not subtracted.";
110  }
111 
112  std::vector<short> uncompressed(rawDigit->Samples());
113  raw::Uncompress(rawDigit->ADCs(), uncompressed, rawDigit->Compression());
114 
115  TH1F* histPtr = fRawDigitHist.get();
116 
117  for (size_t idx = 0; idx < uncompressed.size(); idx++) {
118  float signalVal = float(uncompressed[idx]) - pedestal;
119 
120  histPtr->Fill(float(idx) + 0.5, signalVal);
121  }
122 
123  short minimumVal = *std::min_element(uncompressed.begin(), uncompressed.end());
124  short maximumVal = *std::max_element(uncompressed.begin(), uncompressed.end());
125 
126  fMinimum = float(minimumVal) - pedestal;
127  fMaximum = float(maximumVal) - pedestal;
128 
129  histPtr->SetLineColor(kBlack);
130 
131  // There is only one channel displayed so if here we are done
132  break;
133  }
134  }
135 
136  return;
137  }
138 
139  void DrawRawHist::Draw(const std::string& options, float maxLowVal, float maxHiVal)
140  {
141  TH1F* histPtr = fRawDigitHist.get();
142 
143  // Do we have valid limits to set?
144  histPtr->SetMaximum(maxHiVal);
145  histPtr->SetMinimum(maxLowVal);
146 
147  histPtr->Draw(options.c_str());
148 
149  return;
150  }
151 
152  //......................................................................
153  void DrawRawHist::BookHistogram(raw::ChannelID_t& channel, float startTick, float numTicks)
154  {
157 
158  // Get rid of the previous histograms
159  if (fRawDigitHist.get()) fRawDigitHist.reset();
160 
161  // figure out the signal type for this plane, assume that
162  // plane n in each TPC/cryostat has the same type
163  geo::SigType_t sigType = geo->SignalType(channel);
164  int numBins = numTicks;
165 
166  fRawDigitHist = std::make_unique<TH1F>(
167  "fRAWQHisto", ";t [ticks];q [ADC]", numBins, startTick, startTick + numTicks);
168 
169  TH1F* histPtr = fRawDigitHist.get();
170 
171  histPtr->SetMaximum(cst->fRawQHigh[(size_t)sigType]);
172  histPtr->SetMinimum(cst->fRawQLow[(size_t)sigType]);
173 
174  histPtr->SetLineColor(kBlack);
175  histPtr->SetLineWidth(1);
176 
177  histPtr->GetXaxis()->SetLabelSize(0.10); // was 0.15
178  histPtr->GetXaxis()->SetLabelOffset(0.01); // was 0.00
179  histPtr->GetXaxis()->SetTitleSize(0.10); // was 0.15
180  histPtr->GetXaxis()->SetTitleOffset(0.60); // was 0.80
181 
182  histPtr->GetYaxis()->SetLabelSize(0.10); // was 0.15
183  histPtr->GetYaxis()->SetLabelOffset(0.002); // was 0.00
184  histPtr->GetYaxis()->SetTitleSize(0.10); // was 0.15
185  histPtr->GetYaxis()->SetTitleOffset(0.16); // was 0.80
186  }
187 
189 }
float GetPedestal() const
Definition: RawDigit.h:221
const ADCvector_t & ADCs() const
Reference to the compressed ADC count vector.
Definition: RawDigit.h:209
const art::Event * GetEvent() const
Definition: EventHolder.cxx:45
std::vector< double > fRawQLow
low edge of ADC values for drawing raw digits
#define DEFINE_ART_CLASS_TOOL(tool)
Definition: ToolMacros.h:42
ULong64_t Samples() const
Number of samples in the uncompressed ADC data.
Definition: RawDigit.h:217
void Draw(const std::string &, float, float) override
void Fill(evdb::View2D &, raw::ChannelID_t &, float, float) override
ChannelID_t Channel() const
DAQ channel this raw data was read from.
Definition: RawDigit.h:213
Definition of basic raw digits.
Singleton to hold the current art::Event for the event display.
bool isValid() const noexcept
Definition: Handle.h:203
float getMinimum() const override
The color scales used by the event display.
enum geo::_plane_sigtype SigType_t
Enumerate the possible plane projections.
static EventHolder * Instance()
Definition: EventHolder.cxx:15
Collect all the RawData header files together.
This provides an interface for tools which are tasked with drawing the "wire" data (deconvolved wavef...
float getMaximum() const override
DrawRawHist(const fhicl::ParameterSet &pset)
std::vector< double > fRawQHigh
high edge of ADC values for drawing raw digits
raw::Compress_t Compression() const
Compression algorithm used to store the ADC counts.
Definition: RawDigit.h:229
SigType_t SignalType(PlaneID const &pid) const
Returns the type of signal on the channels of specified TPC plane.
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
std::vector< art::InputTag > fRawDataLabels
module label that made the raw digits, default is daq
void configure(const fhicl::ParameterSet &pset) override
unsigned int ChannelID_t
Type representing the ID of a readout channel.
Definition: RawTypes.h:28
void Uncompress(const std::vector< short > &adc, std::vector< short > &uncompressed, raw::Compress_t compress)
Uncompresses a raw data buffer.
Definition: raw.cxx:744
int fPedestalOption
0: use DetPedestalService; 1: Use pedestal in raw::RawDigt; 2: no ped subtraction ...
Namespace collecting geometry-related classes utilities.
std::unique_ptr< TH1F > fRawDigitHist
void BookHistogram(raw::ChannelID_t &, float, float)
art framework interface to geometry description
Event finding and building.