diff --git a/inc/TRestRawToDetectorSignalProcess.h b/inc/TRestRawToDetectorSignalProcess.h index 4351c18..3530c7c 100644 --- a/inc/TRestRawToDetectorSignalProcess.h +++ b/inc/TRestRawToDetectorSignalProcess.h @@ -73,13 +73,20 @@ class TRestRawToDetectorSignalProcess : public TRestEventProcess { /// A parameter to determine if baseline correction has been applied by a previous process Bool_t fBaseLineCorrection = false; + // Perform signal smoothing to the data + Bool_t fSignalSmoothing = false; + + // In case signal smoothing is active provide averaging points + Int_t fAveragingPoints = 5; + public: RESTValue GetInputEvent() const override { return fInputSignalEvent; } RESTValue GetOutputEvent() const override { return fOutputSignalEvent; } TRestEvent* ProcessEvent(TRestEvent* inputEvent) override; - void ZeroSuppresion(TRestRawSignal* rawSignal, TRestDetectorSignal& sgnl); + void ProcessSignal(TRestRawSignal* rawSignal, TRestDetectorSignal& sgnl); + void ProcessSignalSmoothed(TRestRawSignal* rawSignal, TRestDetectorSignal& sgnl); /// It prints out the process parameters stored in the metadata structure void PrintMetadata() override { @@ -90,6 +97,7 @@ class TRestRawToDetectorSignalProcess : public TRestEventProcess { RESTMetadata << "Gain : " << fGain << RESTendl; if (fZeroSuppression) { + RESTMetadata << "ZeroSuppression is enabled " << RESTendl; RESTMetadata << "Base line range definition : ( " << fBaseLineRange.X() << " , " << fBaseLineRange.Y() << " ) " << RESTendl; RESTMetadata << "Integral range : ( " << fIntegralRange.X() << " , " << fIntegralRange.Y() @@ -99,8 +107,11 @@ class TRestRawToDetectorSignalProcess : public TRestEventProcess { RESTMetadata << "Number of points over threshold : " << fNPointsOverThreshold << RESTendl; } - if (fBaseLineCorrection) - RESTMetadata << "BaseLine correction is enabled for TRestRawSignalAnalysisProcess" << RESTendl; + if (fBaseLineCorrection) RESTMetadata << "BaseLine correction is enabled" << RESTendl; + if (fSignalSmoothing) { + RESTMetadata << "Signal smoothing is enabled" << RESTendl; + RESTMetadata << "Averaging points " << fAveragingPoints << RESTendl; + } EndPrintProcess(); } @@ -117,6 +128,6 @@ class TRestRawToDetectorSignalProcess : public TRestEventProcess { // Destructor ~TRestRawToDetectorSignalProcess(); - ClassDefOverride(TRestRawToDetectorSignalProcess, 2); + ClassDefOverride(TRestRawToDetectorSignalProcess, 3); }; #endif diff --git a/src/TRestRawToDetectorSignalProcess.cxx b/src/TRestRawToDetectorSignalProcess.cxx index ff96a35..cc900d8 100644 --- a/src/TRestRawToDetectorSignalProcess.cxx +++ b/src/TRestRawToDetectorSignalProcess.cxx @@ -55,6 +55,10 @@ /// * **signalThreshold**: The number of sigmas a set of consecutive points /// identified over threshold must be over the baseline fluctuations to be /// finally considered a physical signal. +/// * **signalSmoothing**: Perform smoothing of the signal before saving as +/// a TRestDetectorSignalEvent +/// * **averagingPoints**: In case smoothing is performed set the number of +/// points to be averaged to perform the smoothing /// /// List of observables: /// @@ -69,7 +73,8 @@ /// // A raw signal with 200ns binning will be translated to a /// // TRestDetectorSignalEvent. The new signal will start at time=20us, and its /// // amplitude will be reduced a factor 50. If zeroSuppression is true it will -/// // perform +/// // perform the substraction of the points out of the range of interest. A +/// // smoothing of the signal will be performed using 3 averaging points. /// /// /// @@ -81,6 +86,8 @@ /// /// /// +/// +/// /// /// /// @@ -110,6 +117,8 @@ /// Javier Galan /// 2022-January: Adding ZeroSuppression method /// JuanAn Garcia +/// 2023-April: Added new method to process smoothed data +/// JuanAn Garcia /// /// \class TRestRawToDetectorSignalProcess /// \author Javier Gracia @@ -119,6 +128,8 @@ ///
/// #include "TRestRawToDetectorSignalProcess.h" + +#include "TRestPulseShapeAnalysis.h" using namespace std; ClassImp(TRestRawToDetectorSignalProcess); @@ -163,18 +174,17 @@ TRestEvent* TRestRawToDetectorSignalProcess::ProcessEvent(TRestEvent* inputEvent TRestRawSignal* rawSignal = fInputSignalEvent->GetSignal(n); signal.SetID(rawSignal->GetID()); - if (fZeroSuppression) { - ZeroSuppresion(rawSignal, signal); + if (fSignalSmoothing) { + ProcessSignalSmoothed(rawSignal, signal); } else { - for (int p = 0; p < rawSignal->GetNumberOfPoints(); p++) - if (rawSignal->GetData(p) > fThreshold) - signal.NewPoint(fTriggerStarts + fSampling * p, fGain * rawSignal->GetData(p)); + ProcessSignal(rawSignal, signal); } - if (signal.GetNumberOfPoints() > 0) + if (signal.GetNumberOfPoints() > 0) { fOutputSignalEvent->AddSignal(signal); - else + } else { rejectedSignal++; + } } SetObservableValue("NSignalsRejected", rejectedSignal); @@ -184,13 +194,47 @@ TRestEvent* TRestRawToDetectorSignalProcess::ProcessEvent(TRestEvent* inputEvent return fOutputSignalEvent; } -void TRestRawToDetectorSignalProcess::ZeroSuppresion(TRestRawSignal* rawSignal, TRestDetectorSignal& sgnl) { - rawSignal->InitializePointsOverThreshold(TVector2(fPointThreshold, fSignalThreshold), - fNPointsOverThreshold, 512); +/////////////////////////////////////////////// +/// \brief This function performs the smoothing +/// of the TRestRawSignalEvent before transforming +/// on a TRestDetectorSignalEvent +/// +void TRestRawToDetectorSignalProcess::ProcessSignalSmoothed(TRestRawSignal* rawSignal, + TRestDetectorSignal& sgnl) { + auto smoothed = TRestPulseShapeAnalysis::GetSignalSmoothed(rawSignal->GetData(), fAveragingPoints); + const double baselineSigma = rawSignal->GetBaseLineSigma(); + if (fZeroSuppression) { + auto pOver = TRestPulseShapeAnalysis::GetPointsOverThreshold( + smoothed, fIntegralRange, TVector2(fPointThreshold, fSignalThreshold), fNPointsOverThreshold, 512, + baselineSigma); + for (const auto& [j, data] : pOver) { + sgnl.NewPoint(fTriggerStarts + fSampling * j, fGain * data); + } + } else { + for (size_t p = 0; p < smoothed.size(); p++) { + const double data = smoothed[p]; + if (data > fThreshold) sgnl.NewPoint(fTriggerStarts + fSampling * p, fGain * data); + } + } +} - std::vector pOver = rawSignal->GetPointsOverThreshold(); - for (unsigned int n = 0; n < pOver.size(); n++) { - int j = pOver[n]; - sgnl.NewPoint(fTriggerStarts + fSampling * j, fGain * rawSignal->GetData(j)); +/////////////////////////////////////////////// +/// \brief This function transform a TRestRawSignalEvent +/// into a TRestDetectorSignalEvent without +/// any processing, despite Zero suppression +/// +void TRestRawToDetectorSignalProcess::ProcessSignal(TRestRawSignal* rawSignal, TRestDetectorSignal& sgnl) { + if (fZeroSuppression) { + rawSignal->InitializePointsOverThreshold(TVector2(fPointThreshold, fSignalThreshold), + fNPointsOverThreshold, 512); + + auto pOver = rawSignal->GetPointsOverThreshold(); + for (const auto& [j, data] : pOver) { + sgnl.NewPoint(fTriggerStarts + fSampling * j, fGain * data); + } + } else { + for (int p = 0; p < rawSignal->GetNumberOfPoints(); p++) + if (rawSignal->GetData(p) > fThreshold) + sgnl.NewPoint(fTriggerStarts + fSampling * p, fGain * rawSignal->GetData(p)); } }