14 #include "art_root_io/TFileService.h" 15 #include "cetlib_except/exception.h" 100 if (fOutputHistograms) {
111 dir.make<TH1F>(Form(
"DStopStart_%1zu",
fPlane),
";Delta Stop/Start;", 200, 0., 200.);
113 dir.make<TH1F>(Form(
"DMaxTMinT_%1zu",
fPlane),
";Delta Max/Min Tick;", 200, 0., 200.);
115 dir.make<TH1F>(Form(
"DMaxDMinD_%1zu",
fPlane),
";Delta Max/Min Deriv;", 200, 0., 200.);
122 const recob::Wire::RegionsOfInterest_t::datarange_t& dataRange,
123 const size_t roiStartTick,
124 const size_t channel,
125 const size_t eventCount,
134 const Waveform& waveform = dataRange.data();
137 fWaveformTool->triangleSmooth(rawDerivativeVec, derivativeVec);
140 size_t plane = wids[0].Plane;
141 size_t cryo = wids[0].Cryostat;
142 size_t tpc = wids[0].TPC;
143 size_t wire = wids[0].Wire;
146 hitCandidateVec.clear();
156 if (hitCandidateVec.empty()) {
158 std::cout <<
"** C/T/P: " << cryo <<
"/" << tpc <<
"/" << plane <<
", wire: " << wire
159 <<
" has not hits with input size: " << waveform.size() << std::endl;
164 for (
auto& hitCandidate : hitCandidateVec) {
165 size_t centerIdx = hitCandidate.hitCenter;
167 hitCandidate.hitHeight = waveform.at(centerIdx);
177 Form(
"HitPlane_%1zu/ev%04zu/c%1zut%1zuwire_%05zu", plane, eventCount, cryo, tpc, wire));
179 size_t waveformSize = waveform.size();
180 int waveStart = roiStartTick;
181 int waveStop = waveStart + waveformSize;
183 TProfile* waveHist = dir.make<TProfile>(
184 Form(
"HWfm_%03zu_ctw%01zu-%01zu-%01zu-%05zu", channelCnt, cryo, tpc, plane, wire),
191 TProfile* derivHist = dir.make<TProfile>(
192 Form(
"HDer_%03zu_ctw%01zu-%01zu-%01zu-%05zu", channelCnt, cryo, tpc, plane, wire),
199 TProfile* candHitHist = dir.make<TProfile>(
200 Form(
"HCan_%03zu_ctw%01zu-%01zu-%01zu-%05zu", channelCnt, cryo, tpc, plane, wire),
207 TProfile* maxDerivHist = dir.make<TProfile>(
208 Form(
"HMax_%03zu_ctw%01zu-%01zu-%01zu-%05zu", channelCnt, cryo, tpc, plane, wire),
217 for (
size_t idx = 0; idx < waveform.size(); idx++) {
218 waveHist->Fill(roiStartTick + idx, waveform.at(idx));
219 derivHist->Fill(roiStartTick + idx, derivativeVec.at(idx));
223 for (
const auto& hitCandidate : hitCandidateVec) {
224 candHitHist->Fill(hitCandidate.hitCenter, hitCandidate.hitHeight);
225 maxDerivHist->Fill(hitCandidate.maxTick, hitCandidate.maxDerivative);
226 maxDerivHist->Fill(hitCandidate.minTick, hitCandidate.minDerivative);
228 fDStopStartHist->Fill(hitCandidate.stopTick - hitCandidate.startTick, 1.);
239 const size_t roiStartTick,
241 float dPeakThreshold,
248 std::pair<Waveform::const_iterator, Waveform::const_iterator> minMaxPair =
249 std::minmax_element(startItr, stopItr);
255 if (std::fabs(*maxItr) > std::fabs(*minItr))
260 int deltaTicks = std::distance(maxItr, minItr);
261 float range = *maxItr - *minItr;
264 if (deltaTicks >= dTicksThreshold && range > dPeakThreshold) {
269 int startTick = std::distance(startItr, newEndItr);
275 int stopTick = std::distance(startItr, newStartItr);
278 if (startTick > dTicksThreshold) {
280 if (*(newEndItr - 1) > 0.) {
290 startItr, newEndItr + 1, roiStartTick, dTicksThreshold, dPeakThreshold, hitCandidateVec);
297 std::min_element(maxItr, minItr, [](
const auto&
left,
const auto&
right) {
298 return std::fabs(
left) < std::fabs(
right);
302 if (2 * std::distance(peakItr, minItr) < std::distance(maxItr, peakItr))
304 else if (2 * std::distance(maxItr, peakItr) < std::distance(peakItr, minItr))
307 hitCandidate.
startTick = roiStartTick + startTick;
308 hitCandidate.
stopTick = roiStartTick + stopTick;
309 hitCandidate.
maxTick = roiStartTick + std::distance(startItr, maxItr);
310 hitCandidate.
minTick = roiStartTick + std::distance(startItr, minItr);
311 hitCandidate.
maxDerivative = maxItr != stopItr ? *maxItr : 0.;
312 hitCandidate.
minDerivative = minItr != stopItr ? *minItr : 0.;
313 hitCandidate.
hitCenter = roiStartTick + std::distance(startItr, peakItr) + 0.5;
318 hitCandidateVec.push_back(hitCandidate);
321 if (std::distance(newStartItr, stopItr) > dTicksThreshold) {
323 if (*(newStartItr + 1) < 0.) {
334 roiStartTick + stopTick,
349 if (hitCandidateVec.empty())
return;
357 size_t lastStopTick = hitCandidateVec.front().stopTick;
360 for (
const auto& hitCandidate : hitCandidateVec) {
365 !groupedHitVec.empty()) {
366 mergedHitsVec.emplace_back(groupedHitVec);
368 groupedHitVec.clear();
372 groupedHitVec.emplace_back(hitCandidate);
374 lastStopTick = hitCandidate.stopTick;
379 if (!groupedHitVec.empty()) mergedHitsVec.emplace_back(groupedHitVec);
393 while ((lastItr + 1) != stopItr) {
394 if (*(lastItr + 1) > *lastItr)
break;
411 if (std::distance(startItr, minItr) > 0) {
413 while ((lastItr - 1) != startItr) {
414 if (*(lastItr - 1) < *lastItr)
break;
430 if (std::distance(startItr, lastItr) > 0) {
438 while (loopItr != startItr) {
440 if (*loopItr < 0. || !(*loopItr < *lastItr))
break;
458 if (std::distance(minItr, stopItr) > 1) {
462 while (loopItr != stopItr) {
464 if (*loopItr > 0. || !(*loopItr > *lastItr))
break;
constexpr auto const & right(const_AssnsIter< L, R, D, Dir > const &a, const_AssnsIter< L, R, D, Dir > const &b)
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Interface for a class providing readout channel mapping to geometry.
T get(std::string const &key) const
constexpr auto const & left(const_AssnsIter< L, R, D, Dir > const &a, const_AssnsIter< L, R, D, Dir > const &b)
virtual std::vector< WireID > ChannelToWire(raw::ChannelID_t channel) const =0
This provides an interface for tools which are tasked with finding candidate hits on input waveforms...
art framework interface to geometry description