/** * @file SOnlineNeutronID.cpp * @author Dan Saunders, on behalf of the SoLid collaboration. * @date 17 Feb 2016 */ #include "SOnlineNeutronID.h" //============================================================================== //! Constructor setting up default values. SOnlineNeutronID::SOnlineNeutronID(SDetector * dtr, SClipboard * cb) : ISAlgorithm(dtr, cb, "SOnlineNeutronID"), m_totSampleWindow(40), // As in SM1. m_totThreshold(50), m_nPeaksHtimeWindow(1000) { m_options.push_back(new SOptionInt("TotSampleWindow", &m_totSampleWindow)); m_options.push_back(new SOptionInt("TotThreshold", &m_totThreshold)); m_options.push_back(new SOptionDouble("nPeaksHtimeWindow", &m_nPeaksHtimeWindow)); } //============================================================================== void SOnlineNeutronID::initialize() {} //============================================================================== //! Online neutron ID. /*! Current implimentation is very simple. Loop over all events, and draw histos * of usual nID parameters (tot, nPeaks, I/A). The parameters are summed across * all channels involved in the event (no selections placed of the channels yet). * Since nPeaks, tot, and IonA are not yet found for peaks, this needs to be done * manually. Things will likely move around after the starterkit. */ void SOnlineNeutronID::execute() { for (auto event : (*cb()->events())) { findIonA(event); findTot(event); findNpeaks(event); } } //============================================================================== void SOnlineNeutronID::finalize() {} //============================================================================== //! For each event, assign a value of IonA. Take as the sum from all channels. void SOnlineNeutronID::findIonA(SEvent * event) { double IonA = 0; for (auto p : (*event->peaks())) IonA += p->integral()/p->amplitude(); event->setNID_IonA(IonA); } //============================================================================== //! For each event, assign a value of tot. /*! Find number of samples above some threshold for a preset window (in units of * sampled). Need to loop over the waveforms in the vicinity of the event for * m_totSampleWindow samples, and see if over threshold. Sum over channels. * Use the iterator provided by each peak. Need to be careful not to access * something that doesn't exist. */ void SOnlineNeutronID::findTot(SEvent * event) { // Has to be updated for NZS data. unsigned int tot = 0; for (auto p : (*event->peaks())) { uint thres = uint(m_totThreshold) + uint(p->channel()->baseline()); // Loop over the samples at the peak. auto upperBound = std::min(p->waveform()->samples()->end(), p->waveformPos()+m_totSampleWindow); for (auto iSamp = p->waveformPos(); iSamp != upperBound; iSamp++) { if ((*iSamp) > thres) tot++; } } event->setNID_tot(tot); } //============================================================================== //! For each event, assign a value of nPeaks. /*! Need to look at the channels forming the event, and look for other peaks. * Since Saffron2 has already done peak finding at a low threshold, we'll * count those peaks from the forming channels in the htime window. */ void SOnlineNeutronID::findNpeaks(SEvent * event) { unsigned int nPeaks = 0; for (auto c : (*event->channels())) { // Loop over peaks on this channel. for (auto p : (*c->peaks())) { double delt = event->htime() - p->htime(); if (delt >= 0 && delt < m_nPeaksHtimeWindow) nPeaks++; } } event->setNID_nPeaks(nPeaks); } //==============================================================================