From 9b6bfd0ccf30a10c595023263434a5c60c4a85b9 Mon Sep 17 00:00:00 2001 From: juanan Date: Wed, 15 Feb 2023 18:44:52 +0100 Subject: [PATCH 01/13] Moving to generic methods inside TRestRawSignal --- inc/TRestRawSignal.h | 8 +- src/TRestRawSignal.cxx | 201 +++-------------------------------------- 2 files changed, 17 insertions(+), 192 deletions(-) diff --git a/inc/TRestRawSignal.h b/inc/TRestRawSignal.h index 1359d292..a193f85e 100644 --- a/inc/TRestRawSignal.h +++ b/inc/TRestRawSignal.h @@ -36,10 +36,6 @@ class TRestRawSignal : public TObject { private: void CalculateThresholdIntegral(); - void CalculateBaseLineSigmaSD(Int_t startBin, Int_t endBin); - - void CalculateBaseLineSigmaIQR(Int_t startBin, Int_t endBin); - std::vector GetSignalSmoothed_ExcludeOutliers(Int_t averagingPoints); protected: @@ -111,6 +107,8 @@ class TRestRawSignal : public TObject { /// Returns the range defined by user inline TVector2 GetRange() const { return fRange; } + inline std::vector GetSignalData() const { return fSignalData; } + /// Returns false if the baseline and its baseline fluctuation was not initialized. inline Bool_t isBaseLineInitialized() { if (fBaseLineSigma == 0 && fBaseLine == 0) return false; @@ -187,8 +185,6 @@ class TRestRawSignal : public TObject { void GetDifferentialSignal(TRestRawSignal* diffSignal, Int_t smearPoints); - void GetSignalSmoothed(TRestRawSignal* smoothedSignal, Int_t averagingPoints); - std::vector GetSignalSmoothed(Int_t averagingPoints, std::string option = ""); void GetWhiteNoiseSignal(TRestRawSignal* noiseSignal, Double_t noiseLevel = 1.); diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 0eb276f7..2be20f2c 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -59,6 +59,8 @@ /// #include "TRestRawSignal.h" +#include "TRestSignalAnalysis.h" + #include #include #include @@ -221,57 +223,8 @@ void TRestRawSignal::IncreaseBinBy(Int_t bin, Double_t data) { /// void TRestRawSignal::InitializePointsOverThreshold(const TVector2& thrPar, Int_t nPointsOver, Int_t nPointsFlat) { - if (fRange.X() < 0) fRange.SetX(0); - if (fRange.Y() <= 0) fRange.SetY(GetNumberOfPoints()); - - fPointsOverThreshold.clear(); - - double pointTh = thrPar.X(); - double signalTh = thrPar.Y(); - - double threshold = pointTh * fBaseLineSigma; - - for (int i = fRange.X(); i < fRange.Y(); i++) { - // Filling a pulse with consecutive points that are over threshold - if (this->GetData(i) > threshold) { - int pos = i; - std::vector pulse; - pulse.push_back(this->GetData(i)); - i++; - - // If the pulse ends in a flat end above the threshold, the parameter - // nPointsFlat will serve to artificially end the pulse. - // If nPointsFlat is big enough, this parameter will not affect the - // decision to cut this anomalous behaviour. And all points over threshold - // will be added to the pulse vector. - int flatN = 0; - while (i < fRange.Y() && this->GetData(i) > threshold) { - if (TMath::Abs(this->GetData(i) - this->GetData(i - 1)) > threshold) { - flatN = 0; - } else { - flatN++; - } - - if (flatN < nPointsFlat) { - pulse.push_back(this->GetData(i)); - i++; - } else { - break; - } - } - - if (pulse.size() >= (unsigned int)nPointsOver) { - // auto stdev = TMath::StdDev(begin(pulse), end(pulse)); - // calculate stdev - double mean = std::accumulate(pulse.begin(), pulse.end(), 0.0) / pulse.size(); - double sq_sum = std::inner_product(pulse.begin(), pulse.end(), pulse.begin(), 0.0); - double stdev = std::sqrt(sq_sum / pulse.size() - mean * mean); - - if (stdev > signalTh * fBaseLineSigma) - for (int j = pos; j < i; j++) fPointsOverThreshold.push_back(j); - } - } - } + fPointsOverThreshold = TRestSignalAnalysis::GetPointsOverThreshold( + fSignalData, fRange, thrPar, nPointsOver, nPointsFlat, fBaseLine, fBaseLineSigma); CalculateThresholdIntegral(); } @@ -599,32 +552,6 @@ void TRestRawSignal::GetWhiteNoiseSignal(TRestRawSignal* noiseSignal, Double_t n delete fRandom; } -/////////////////////////////////////////////// -/// \brief It smoothes the existing signal and places it at the signal pointer -/// given by argument. -/// -/// \param averagingPoints It defines the number of neighbour consecutive -/// points used to average the signal -/// -void TRestRawSignal::GetSignalSmoothed(TRestRawSignal* smoothedSignal, Int_t averagingPoints) { - smoothedSignal->Initialize(); - - averagingPoints = (averagingPoints / 2) * 2 + 1; // make it odd >= averagingPoints - - Double_t sumAvg = GetIntegralInRange(0, averagingPoints) / averagingPoints; - - for (int i = 0; i <= averagingPoints / 2; i++) smoothedSignal->AddPoint((Short_t)sumAvg); - - for (int i = averagingPoints / 2 + 1; i < GetNumberOfPoints() - averagingPoints / 2; i++) { - sumAvg -= this->GetRawData(i - (averagingPoints / 2 + 1)) / averagingPoints; - sumAvg += this->GetRawData(i + averagingPoints / 2) / averagingPoints; - smoothedSignal->AddPoint((Short_t)sumAvg); - } - - for (int i = GetNumberOfPoints() - averagingPoints / 2; i < GetNumberOfPoints(); i++) - smoothedSignal->AddPoint(sumAvg); -} - /////////////////////////////////////////////// /// \brief It smoothes the existing signal and returns it in a vector of Float_t values /// @@ -635,31 +562,15 @@ void TRestRawSignal::GetSignalSmoothed(TRestRawSignal* smoothedSignal, Int_t ave /// baseline will be ignored to improve the smoothing result /// std::vector TRestRawSignal::GetSignalSmoothed(Int_t averagingPoints, std::string option) { - std::vector result; - if (option == "") { - result.resize(GetNumberOfPoints()); - - averagingPoints = (averagingPoints / 2) * 2 + 1; // make it odd >= averagingPoints - - Float_t sumAvg = (Float_t)GetIntegralInRange(0, averagingPoints) / averagingPoints; - - for (int i = 0; i <= averagingPoints / 2; i++) result[i] = sumAvg; - - for (int i = averagingPoints / 2 + 1; i < GetNumberOfPoints() - averagingPoints / 2; i++) { - sumAvg -= this->GetRawData(i - (averagingPoints / 2 + 1)) / averagingPoints; - sumAvg += this->GetRawData(i + averagingPoints / 2) / averagingPoints; - result[i] = sumAvg; - } - - for (int i = GetNumberOfPoints() - averagingPoints / 2; i < GetNumberOfPoints(); i++) - result[i] = sumAvg; + return TRestSignalAnalysis::GetSignalSmoothed(fSignalData, averagingPoints); } else if (ToUpper(option) == "EXCLUDE OUTLIERS") { - result = GetSignalSmoothed_ExcludeOutliers(averagingPoints); + return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints); } else { - cout << "TRestRawSignal::GetSignalSmoothed. Error! No such option!" << endl; + std::cout << "TRestRawSignal::GetSignalSmoothed. Error! No such option!" << std::endl; + std::vector result; + return result; } - return result; } /////////////////////////////////////////////// @@ -672,52 +583,8 @@ std::vector TRestRawSignal::GetSignalSmoothed(Int_t averagingPoints, st /// points used to average the signal /// std::vector TRestRawSignal::GetSignalSmoothed_ExcludeOutliers(Int_t averagingPoints) { - std::vector result(GetNumberOfPoints()); - - if (fBaseLine == 0) CalculateBaseLine(5, GetNumberOfPoints() - 5, "ROBUST"); - - averagingPoints = (averagingPoints / 2) * 2 + 1; // make it odd >= averagingPoints - - Float_t sumAvg = (Float_t)GetIntegralInRange(0, averagingPoints) / averagingPoints; - - // Points at the beginning, where we can calculate a moving average - for (int i = 0; i <= averagingPoints / 2; i++) result[i] = sumAvg; - - // Points in the middle - float_t amplitude; - for (int i = averagingPoints / 2 + 1; i < GetNumberOfPoints() - averagingPoints / 2; i++) { - amplitude = this->GetRawData(i - (averagingPoints / 2 + 1)); - sumAvg -= (std::abs(amplitude - fBaseLine) > 3 * fBaseLineSigma) ? fBaseLine / averagingPoints - : amplitude / averagingPoints; - amplitude = this->GetRawData(i + averagingPoints / 2); - sumAvg += (std::abs(amplitude - fBaseLine) > 3 * fBaseLineSigma) ? fBaseLine / averagingPoints - : amplitude / averagingPoints; - result[i] = sumAvg; - } - - // Points at the end, where we can calculate a moving average - for (int i = GetNumberOfPoints() - averagingPoints / 2; i < GetNumberOfPoints(); i++) result[i] = sumAvg; - return result; -} - -/////////////////////////////////////////////// -/// \brief It applies the moving average filter (GetSignalSmoothed) to the signal, which is then subtracted -/// from the raw data, resulting in a corrected baseline. The returned signal is placed at the signal pointer -/// given by the argument. -/// -/// \param smoothedSignal The pointer to the TRestRawSignal which will contain the corrected signal -/// -/// \param averagingPoints It defines the number of neighbour consecutive -/// points used to average the signal -/// -void TRestRawSignal::GetBaseLineCorrected(TRestRawSignal* smoothedSignal, Int_t averagingPoints) { - smoothedSignal->Initialize(); - - std::vector averagedSignal = GetSignalSmoothed(averagingPoints, "EXCLUDE OUTLIERS"); - - for (int i = 0; i < GetNumberOfPoints(); i++) { - smoothedSignal->AddPoint(GetRawData(i) - averagedSignal[i]); - } + return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, fBaseLine, + fBaseLineSigma); } /////////////////////////////////////////////// @@ -772,49 +639,11 @@ void TRestRawSignal::CalculateBaseLineMedian(Int_t startBin, Int_t endBin) { /// void TRestRawSignal::CalculateBaseLine(Int_t startBin, Int_t endBin, const std::string& option) { if (ToUpper(option) == "ROBUST") { - CalculateBaseLineMedian(startBin, endBin); - CalculateBaseLineSigmaIQR(startBin, endBin); + TRestSignalAnalysis::CalculateBaselineAndSigmaSD(fSignalData, startBin, endBin, fBaseLine, + fBaseLineSigma); } else { - CalculateBaseLineMean(startBin, endBin); - CalculateBaseLineSigmaSD(startBin, endBin); - } -} - -/////////////////////////////////////////////// -/// \brief This method is called by CalculateBaseLine to -/// determine the value of the baseline -/// fluctuation as its standard deviation in the baseline range provided. -/// -void TRestRawSignal::CalculateBaseLineSigmaSD(Int_t startBin, Int_t endBin) { - if (endBin - startBin <= 0) { - fBaseLineSigma = 0; - } else { - Double_t baseLineSigma = 0; - for (int i = startBin; i < endBin; i++) - baseLineSigma += (fBaseLine - fSignalData[i]) * (fBaseLine - fSignalData[i]); - fBaseLineSigma = TMath::Sqrt(baseLineSigma / (endBin - startBin)); - } -} - -/////////////////////////////////////////////// -/// \brief This method is called by CalculateBaseLine with the "ROBUST"-option to -/// determine the value of the baseline -/// fluctuation as its interquartile range (IQR) in the baseline range provided. The IQR is more robust -/// towards outliers than the standard deviation. -/// -void TRestRawSignal::CalculateBaseLineSigmaIQR(Int_t startBin, Int_t endBin) { - if (endBin - startBin <= 0) { - fBaseLineSigma = 0; - } else { - vector::const_iterator first = fSignalData.begin() + startBin; - vector::const_iterator last = fSignalData.begin() + endBin; - vector v(first, last); - std::sort(v.begin(), v.end()); - Short_t Q1 = v[(int)(endBin - startBin) * 0.25]; - Short_t Q3 = v[(int)(endBin - startBin) * 0.75]; - Double_t IQR = Q3 - Q1; - fBaseLineSigma = - IQR / 1.349; // IQR/1.349 equals the standard deviation in case of normally distributed data + TRestSignalAnalysis::CalculateBaselineAndSigmaIQR(fSignalData, startBin, endBin, fBaseLine, + fBaseLineSigma); } } From 9d6054e235ad237ab3700d02e1e40e0a91a6dd93 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 15 Feb 2023 17:48:38 +0000 Subject: [PATCH 02/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/TRestRawSignal.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 2be20f2c..0029559b 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -59,8 +59,6 @@ /// #include "TRestRawSignal.h" -#include "TRestSignalAnalysis.h" - #include #include #include @@ -68,6 +66,8 @@ #include +#include "TRestSignalAnalysis.h" + using namespace std; ClassImp(TRestRawSignal); From 5282c8aa7d20b5e4341ad668499f97254203dab6 Mon Sep 17 00:00:00 2001 From: juanan Date: Wed, 15 Feb 2023 22:52:50 +0100 Subject: [PATCH 03/13] Addressing minor minor bugs in the new implementation --- src/TRestRawSignal.cxx | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 0029559b..d40e7871 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -565,7 +565,8 @@ std::vector TRestRawSignal::GetSignalSmoothed(Int_t averagingPoints, st if (option == "") { return TRestSignalAnalysis::GetSignalSmoothed(fSignalData, averagingPoints); } else if (ToUpper(option) == "EXCLUDE OUTLIERS") { - return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints); + return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, fBaseLine, + fBaseLineSigma); } else { std::cout << "TRestRawSignal::GetSignalSmoothed. Error! No such option!" << std::endl; std::vector result; @@ -639,11 +640,11 @@ void TRestRawSignal::CalculateBaseLineMedian(Int_t startBin, Int_t endBin) { /// void TRestRawSignal::CalculateBaseLine(Int_t startBin, Int_t endBin, const std::string& option) { if (ToUpper(option) == "ROBUST") { - TRestSignalAnalysis::CalculateBaselineAndSigmaSD(fSignalData, startBin, endBin, fBaseLine, - fBaseLineSigma); - } else { TRestSignalAnalysis::CalculateBaselineAndSigmaIQR(fSignalData, startBin, endBin, fBaseLine, fBaseLineSigma); + } else { + TRestSignalAnalysis::CalculateBaselineAndSigmaSD(fSignalData, startBin, endBin, fBaseLine, + fBaseLineSigma); } } From 9475fa4869dbd5d3ad34c1766e6e685871d3e6d4 Mon Sep 17 00:00:00 2001 From: juanan Date: Wed, 15 Feb 2023 23:30:44 +0100 Subject: [PATCH 04/13] Reverting removal of GetBaseLineCorredted erased by mistake --- src/TRestRawSignal.cxx | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index d40e7871..ec28d9c3 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -574,6 +574,26 @@ std::vector TRestRawSignal::GetSignalSmoothed(Int_t averagingPoints, st } } +/////////////////////////////////////////////// +/// \brief It applies the moving average filter (GetSignalSmoothed) to the signal, which is then subtracted +/// from the raw data, resulting in a corrected baseline. The returned signal is placed at the signal pointer +/// given by the argument. +/// +/// \param smoothedSignal The pointer to the TRestRawSignal which will contain the corrected signal +/// +/// \param averagingPoints It defines the number of neighbour consecutive +/// points used to average the signal +/// +void TRestRawSignal::GetBaseLineCorrected(TRestRawSignal* smoothedSignal, Int_t averagingPoints) { + smoothedSignal->Initialize(); + + auto averagedSignal = GetSignalSmoothed(averagingPoints, "EXCLUDE OUTLIERS"); + + for (int i = 0; i < GetNumberOfPoints(); i++) { + smoothedSignal->AddPoint(GetRawData(i) - averagedSignal[i]); + } +} + /////////////////////////////////////////////// /// \brief It smooths the existing signal and returns it in a vector of Float_t values. This method excludes /// points which are far off from the BaseLine IQR (e.g. signals). In case the baseline parameters were not From d914d972b9b8aa110bc08cbcfe2a93ede5807f08 Mon Sep 17 00:00:00 2001 From: juanan Date: Thu, 16 Feb 2023 18:31:00 +0100 Subject: [PATCH 05/13] Moving pointsOverThreshold to a vector of pairs instead of a vector of Ints --- inc/TRestRawSignal.h | 8 +++--- src/TRestRawSignal.cxx | 39 ++++++++++++++++++----------- src/TRestRawSignalEvent.cxx | 26 ++++++------------- src/TRestRawSignalViewerProcess.cxx | 30 +++++++--------------- 4 files changed, 47 insertions(+), 56 deletions(-) diff --git a/inc/TRestRawSignal.h b/inc/TRestRawSignal.h index a193f85e..b6863660 100644 --- a/inc/TRestRawSignal.h +++ b/inc/TRestRawSignal.h @@ -52,7 +52,7 @@ class TRestRawSignal : public TObject { TGraph* fGraph; //! /// A std::vector containing the index of points that are identified over threshold. - std::vector fPointsOverThreshold; //! + std::vector > fPointsOverThreshold; //! /// It stores the integral value obtained from the points identified over threshold. Double_t fThresholdIntegral = -1; //! @@ -76,13 +76,15 @@ class TRestRawSignal : public TObject { inline Int_t GetSignalID() const { return fSignalID; } /// Returns the value of signal ID - inline Int_t GetID() const { return fSignalID; } + inline Int_t GetID() const { return GetSignalID(); } /// Returns the actual number of points, or size of the signal inline Int_t GetNumberOfPoints() const { return fSignalData.size(); } /// Returns a std::vector containing the indexes of data points over threshold - inline std::vector GetPointsOverThreshold() const { return fPointsOverThreshold; } + inline std::vector > GetPointsOverThreshold() const { + return fPointsOverThreshold; + } /// Returns the maximum value found in the data points. It includes baseline correction inline Double_t GetMaxValue() { return GetMaxPeakValue(); } diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index ec28d9c3..9ffe55b0 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -240,9 +240,9 @@ void TRestRawSignal::CalculateThresholdIntegral() { fThresholdIntegral = 0; - for (unsigned int n = 0; n < fPointsOverThreshold.size(); n++) { - if (fPointsOverThreshold[n] >= fRange.X() && fPointsOverThreshold[n] < fRange.Y()) { - fThresholdIntegral += GetData(fPointsOverThreshold[n]); + for (const auto& [bin, value] : fPointsOverThreshold) { + if (bin >= fRange.X() && bin < fRange.Y()) { + fThresholdIntegral += value; } } } @@ -303,12 +303,11 @@ Double_t TRestRawSignal::GetSlopeIntegral() { << endl; Double_t sum = 0; - for (unsigned int n = 0; n < fPointsOverThreshold.size(); n++) { - if (n + 1 >= fPointsOverThreshold.size()) return sum; - - sum += GetData(fPointsOverThreshold[n]); - - if (GetData(fPointsOverThreshold[n + 1]) - GetData(fPointsOverThreshold[n]) < 0) break; + Double_t pVal = 0; + for (const auto& [index, val] : fPointsOverThreshold) { + if (val - pVal < 0) break; + sum += val; + pVal = val; } return sum; } @@ -330,13 +329,18 @@ Double_t TRestRawSignal::GetRiseSlope() { return 0; } - Int_t maxBin = GetMaxPeakBin() - 1; + auto max = std::max_element(std::begin(fPointsOverThreshold), std::end(fPointsOverThreshold), + [](const auto& p1, const auto& p2) { return p1.second < p2.second; }); + + Int_t maxBin = max->first; + + Int_t startBin = fPointsOverThreshold.front().first; - Double_t hP = GetData(maxBin); + Double_t hP = max->second; - Double_t lP = GetData(fPointsOverThreshold[0]); + Double_t lP = fPointsOverThreshold.front().second; - return (hP - lP) / (maxBin - fPointsOverThreshold[0]); + return (hP - lP) / (maxBin - startBin - 1); } /////////////////////////////////////////////// @@ -354,7 +358,14 @@ Int_t TRestRawSignal::GetRiseTime() { return 0; } - return GetMaxPeakBin() - fPointsOverThreshold[0]; + auto max = std::max_element(std::begin(fPointsOverThreshold), std::end(fPointsOverThreshold), + [](const auto& p1, const auto& p2) { return p1.second < p2.second; }); + + Int_t maxBin = max->first; + + Int_t startBin = fPointsOverThreshold.front().first; + + return maxBin - startBin; } /////////////////////////////////////////////// diff --git a/src/TRestRawSignalEvent.cxx b/src/TRestRawSignalEvent.cxx index f49cf150..b6b8a229 100644 --- a/src/TRestRawSignalEvent.cxx +++ b/src/TRestRawSignalEvent.cxx @@ -757,30 +757,20 @@ TPad* TRestRawSignalEvent::DrawSignal(Int_t signalID, TString option) { gr2->Draw("CP"); - vector pOver = signal->GetPointsOverThreshold(); + auto pOver = signal->GetPointsOverThreshold(); - TGraph* gr3[5]; - Int_t nGraphs = 0; - gr3[nGraphs] = new TGraph(); - gr3[nGraphs]->SetLineWidth(2); - gr3[nGraphs]->SetLineColor(3); + TGraph* gr3; + gr3 = new TGraph(); + gr3->SetLineWidth(2); + gr3->SetLineColor(3); Int_t point = 0; Int_t nPoints = pOver.size(); - for (int n = 0; n < nPoints; n++) { - gr3[nGraphs]->SetPoint(point, pOver[n], signal->GetData(pOver[n])); + for (const auto& [index, data] : pOver) { + gr3->SetPoint(point, index, data); point++; - if (n + 1 < nPoints && pOver[n + 1] - pOver[n] > 1) { - gr3[nGraphs]->Draw("CP"); - nGraphs++; - if (nGraphs > 4) cout << "Ngraphs : " << nGraphs << endl; - point = 0; - gr3[nGraphs] = new TGraph(); - gr3[nGraphs]->SetLineWidth(2); - gr3[nGraphs]->SetLineColor(3); // Green - } } - if (nPoints > 0) gr3[nGraphs]->Draw("CP"); + if (nPoints > 0) gr3->Draw("CP"); return fPad; } diff --git a/src/TRestRawSignalViewerProcess.cxx b/src/TRestRawSignalViewerProcess.cxx index 5705e341..372c6093 100644 --- a/src/TRestRawSignalViewerProcess.cxx +++ b/src/TRestRawSignalViewerProcess.cxx @@ -256,32 +256,20 @@ TPad* TRestRawSignalViewerProcess::DrawSignal(Int_t signal) { gr2->Draw("CP"); - vector pOver = sgnl->GetPointsOverThreshold(); - - TGraph* gr3[5]; - Int_t nGraphs = 0; - gr3[nGraphs] = new TGraph(); - fDrawingObjects.push_back((TObject*)gr3[nGraphs]); - gr3[nGraphs]->SetLineWidth(2); - gr3[nGraphs]->SetLineColor(3); + auto pOver = sgnl->GetPointsOverThreshold(); + + TGraph* gr3; + gr3 = new TGraph(); + gr3->SetLineWidth(2); + gr3->SetLineColor(3); Int_t point = 0; Int_t nPoints = pOver.size(); - for (int n = 0; n < nPoints; n++) { - gr3[nGraphs]->SetPoint(point, pOver[n], sgnl->GetData(pOver[n])); + for (const auto& [index, data] : pOver) { + gr3->SetPoint(point, index, data); point++; - if (n + 1 < nPoints && pOver[n + 1] - pOver[n] > 1) { - gr3[nGraphs]->Draw("CP"); - nGraphs++; - if (nGraphs > 4) cout << "Ngraphs : " << nGraphs << endl; - point = 0; - gr3[nGraphs] = new TGraph(); - fDrawingObjects.push_back((TObject*)gr3[nGraphs]); - gr3[nGraphs]->SetLineWidth(2); - gr3[nGraphs]->SetLineColor(3); - } } - if (nPoints > 0) gr3[nGraphs]->Draw("CP"); + if (nPoints > 0) gr3->Draw("CP"); return pad; } From 384530bb89d5bc0cee103e3f67e7f60f090d96d9 Mon Sep 17 00:00:00 2001 From: juanan Date: Fri, 31 Mar 2023 10:23:59 +0200 Subject: [PATCH 06/13] Implementing further generic methods inside TRestRawSignal --- inc/TRestRawSignal.h | 6 +- src/TRestRawSignal.cxx | 209 ++++++----------------------------------- 2 files changed, 33 insertions(+), 182 deletions(-) diff --git a/inc/TRestRawSignal.h b/inc/TRestRawSignal.h index b6863660..132452fd 100644 --- a/inc/TRestRawSignal.h +++ b/inc/TRestRawSignal.h @@ -117,6 +117,8 @@ class TRestRawSignal : public TObject { return true; } + std::vector GetData() const; + Double_t GetData(Int_t n) const; Double_t GetRawData(Int_t n) const; @@ -191,10 +193,6 @@ class TRestRawSignal : public TObject { void GetWhiteNoiseSignal(TRestRawSignal* noiseSignal, Double_t noiseLevel = 1.); - void CalculateBaseLineMean(Int_t startBin, Int_t endBin); - - void CalculateBaseLineMedian(Int_t startBin, Int_t endBin); - void CalculateBaseLine(Int_t startBin, Int_t endBin, const std::string& option = ""); void GetBaseLineCorrected(TRestRawSignal* smoothedSignal, Int_t averagingPoints); diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 9ffe55b0..36e389ae 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -65,6 +65,7 @@ #include #include +#include #include "TRestSignalAnalysis.h" @@ -170,6 +171,18 @@ Short_t TRestRawSignal::operator[](Int_t n) { /// Double_t TRestRawSignal::GetData(Int_t n) const { return (Double_t)fSignalData[n] - fBaseLine; } +/////////////////////////////////////////////// +/// \brief It returns a vector of the data with the +/// baseline substracted +/// +std::vector TRestRawSignal::GetData() const { + std::vector signalData(GetNumberOfPoints()); + + for (int i = 0; i < GetNumberOfPoints(); i++) signalData[i] = GetData(i); + + return signalData; +} + /////////////////////////////////////////////// /// \brief It returns the original data value of point *n* without baseline /// correction. @@ -256,9 +269,7 @@ Double_t TRestRawSignal::GetIntegral() { if (fRange.X() < 0) fRange.SetX(0); if (fRange.Y() <= 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); - Double_t sum = 0; - for (int i = fRange.X(); i < fRange.Y(); i++) sum += GetData(i); - return sum; + return TRestSignalAnalysis::GetIntegral(GetData(), fRange.X(), fRange.Y()); } /////////////////////////////////////////////// @@ -266,12 +277,7 @@ Double_t TRestRawSignal::GetIntegral() { /// by (startBin,endBin). /// Double_t TRestRawSignal::GetIntegralInRange(Int_t startBin, Int_t endBin) { - if (startBin < 0) startBin = 0; - if (endBin <= 0 || endBin > GetNumberOfPoints()) endBin = GetNumberOfPoints(); - - Double_t sum = 0; - for (int i = startBin; i < endBin; i++) sum += GetRawData(i); - return sum; + return TRestSignalAnalysis::GetIntegral(fSignalData, startBin, endBin); } /////////////////////////////////////////////// @@ -297,19 +303,13 @@ Double_t TRestRawSignal::GetThresholdIntegral() { /// signal. InitializePointsOverThreshold should have been called first. /// Double_t TRestRawSignal::GetSlopeIntegral() { - if (fThresholdIntegral == -1) + if (fThresholdIntegral == -1) { cout << "REST Warning. TRestRawSignal::GetSlopeIntegral. " "InitializePointsOverThreshold should be called first." << endl; - - Double_t sum = 0; - Double_t pVal = 0; - for (const auto& [index, val] : fPointsOverThreshold) { - if (val - pVal < 0) break; - sum += val; - pVal = val; } - return sum; + + return TRestSignalAnalysis::GetSlopeIntegral(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -324,23 +324,7 @@ Double_t TRestRawSignal::GetRiseSlope() { "InitializePointsOverThreshold should be called first." << endl; - if (fPointsOverThreshold.size() < 2) { - // cout << "REST Warning. TRestRawSignal::GetRiseSlope. Less than 2 points!." << endl; - return 0; - } - - auto max = std::max_element(std::begin(fPointsOverThreshold), std::end(fPointsOverThreshold), - [](const auto& p1, const auto& p2) { return p1.second < p2.second; }); - - Int_t maxBin = max->first; - - Int_t startBin = fPointsOverThreshold.front().first; - - Double_t hP = max->second; - - Double_t lP = fPointsOverThreshold.front().second; - - return (hP - lP) / (maxBin - startBin - 1); + return TRestSignalAnalysis::GetRiseSlope(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -353,19 +337,7 @@ Int_t TRestRawSignal::GetRiseTime() { "InitializePointsOverThreshold should be called first." << endl; - if (fPointsOverThreshold.size() < 2) { - // cout << "REST Warning. TRestRawSignal::GetRiseTime. Less than 2 points!." << endl; - return 0; - } - - auto max = std::max_element(std::begin(fPointsOverThreshold), std::end(fPointsOverThreshold), - [](const auto& p1, const auto& p2) { return p1.second < p2.second; }); - - Int_t maxBin = max->first; - - Int_t startBin = fPointsOverThreshold.front().first; - - return maxBin - startBin; + return TRestSignalAnalysis::GetRiseTime(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -373,33 +345,8 @@ Int_t TRestRawSignal::GetRiseTime() { /// and its neightbour points. /// Double_t TRestRawSignal::GetTripleMaxIntegral() { - if (fRange.X() < 0) fRange.SetX(0); - if (fRange.Y() <= 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); - - if (fThresholdIntegral == -1) { - cout << "REST Warning. TRestRawSignal::GetTripleMaxIntegral. " - "InitializePointsOverThreshold should be called first." - << endl; - return 0; - } - - if (fPointsOverThreshold.size() < 2) { - // cout << "REST Warning. TRestRawSignal::GetTripleMaxIntegral. Points over - // " - // "threshold = " - // << fPointsOverThreshold.size() << endl; - return 0; - } - - Int_t cBin = GetMaxPeakBin(); - - if (cBin + 1 >= GetNumberOfPoints()) return 0; - - Double_t a1 = GetData(cBin); - Double_t a2 = GetData(cBin - 1); - Double_t a3 = GetData(cBin + 1); - - return a1 + a2 + a3; + auto gr = GetGraph(); + return TRestSignalAnalysis::GetTripleMaxIntegral(gr); } /////////////////////////////////////////////// @@ -409,36 +356,14 @@ Double_t TRestRawSignal::GetTripleMaxIntegral() { Double_t TRestRawSignal::GetAverageInRange(Int_t startBin, Int_t endBin) { if (startBin < 0) startBin = 0; if (endBin <= 0 || endBin > GetNumberOfPoints()) endBin = GetNumberOfPoints(); - - Double_t sum = 0; - for (int i = startBin; i <= endBin; i++) sum += this->GetData(i); - - return sum / (endBin - startBin + 1); + return TRestSignalAnalysis::GetAverage(GetData(), startBin, endBin); } /////////////////////////////////////////////// /// \brief It returns the temporal width of the peak with maximum amplitude /// inside the signal /// -Int_t TRestRawSignal::GetMaxPeakWidth() { - Int_t mIndex = this->GetMaxPeakBin(); - Double_t maxValue = this->GetData(mIndex); - - Double_t value = maxValue; - Int_t rightIndex = mIndex; - while (value > maxValue / 2) { - value = this->GetData(rightIndex); - rightIndex++; - } - Int_t leftIndex = mIndex; - value = maxValue; - while (value > maxValue / 2) { - value = this->GetData(leftIndex); - leftIndex--; - } - - return rightIndex - leftIndex; -} +Int_t TRestRawSignal::GetMaxPeakWidth() { return TRestSignalAnalysis::GetMaxPeakWidth(GetData()); } /////////////////////////////////////////////// /// \brief It returns the amplitude of the signal maximum, baseline will be @@ -450,23 +375,7 @@ Double_t TRestRawSignal::GetMaxPeakValue() { return GetData(GetMaxPeakBin()); } /////////////////////////////////////////////// /// \brief It returns the bin at which the maximum peak amplitude happens /// -Int_t TRestRawSignal::GetMaxPeakBin() { - Double_t max = numeric_limits::min(); - - Int_t index = 0; - - if (fRange.Y() == 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); - if (fRange.X() < 0) fRange.SetX(0); - - for (int i = fRange.X(); i < fRange.Y(); i++) { - if (this->GetData(i) > max) { - max = GetData(i); - index = i; - } - } - - return index; -} +Int_t TRestRawSignal::GetMaxPeakBin() { return TRestSignalAnalysis::GetMaxBin(fSignalData); } /////////////////////////////////////////////// /// \brief It returns the amplitude of the signal minimum, baseline will be @@ -478,22 +387,7 @@ Double_t TRestRawSignal::GetMinPeakValue() { return GetData(GetMinPeakBin()); } /////////////////////////////////////////////// /// \brief It returns the bin at which the minimum peak amplitude happens /// -Int_t TRestRawSignal::GetMinPeakBin() { - Double_t min = numeric_limits::max(); - Int_t index = 0; - - if (fRange.Y() == 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); - if (fRange.X() < 0) fRange.SetX(0); - - for (int i = fRange.X(); i < fRange.Y(); i++) { - if (this->GetData(i) < min) { - min = GetData(i); - index = i; - } - } - - return index; -} +Int_t TRestRawSignal::GetMinPeakBin() { return TRestSignalAnalysis::GetMinBin(fSignalData); } /////////////////////////////////////////////// /// \brief It returns whether the signal has ADC saturation @@ -552,15 +446,14 @@ void TRestRawSignal::GetDifferentialSignal(TRestRawSignal* diffSignal, Int_t sme /// as its standard deviation. /// void TRestRawSignal::GetWhiteNoiseSignal(TRestRawSignal* noiseSignal, Double_t noiseLevel) { - double* dd = new double(); - uintptr_t seed = (uintptr_t)dd + (uintptr_t)this; - delete dd; - TRandom3* fRandom = new TRandom3(seed); + std::mt19937::result_type seed = std::random_device()(); + std::mt19937 gen(seed); + + TRandom3 fRandom(seed); for (int i = 0; i < GetNumberOfPoints(); i++) { - noiseSignal->AddPoint(this->GetData(i) + (Short_t)fRandom->Gaus(0, noiseLevel)); + noiseSignal->AddPoint(this->GetData(i) + (Short_t)fRandom.Gaus(0, noiseLevel)); } - delete fRandom; } /////////////////////////////////////////////// @@ -619,46 +512,6 @@ std::vector TRestRawSignal::GetSignalSmoothed_ExcludeOutliers(Int_t ave fBaseLineSigma); } -/////////////////////////////////////////////// -/// \brief This method is called by CalculateBaseLine and is used to determine the value of the baseline as -/// average (arithmetic mean) of the data points found -/// in the range defined between startBin and endBin. -/// -void TRestRawSignal::CalculateBaseLineMean(Int_t startBin, Int_t endBin) { - if (endBin - startBin <= 0) { - fBaseLine = 0.; - } else if (endBin > (int)fSignalData.size()) { - cout << "TRestRawSignal::CalculateBaseLine. Error! Baseline range exceeds the rawdata depth!!" - << endl; - endBin = fSignalData.size(); - } else { - Double_t baseLine = 0; - for (int i = startBin; i < endBin; i++) baseLine += fSignalData[i]; - fBaseLine = baseLine / (endBin - startBin); - } -} - -/////////////////////////////////////////////// -/// \brief This method is called by CalculateBaseLine with the "ROBUST"-option and is used to determine the -/// value of the baseline as the median of the data points found in the range defined between startBin and -/// endBin. -/// -void TRestRawSignal::CalculateBaseLineMedian(Int_t startBin, Int_t endBin) { - if (endBin - startBin <= 0) { - fBaseLine = 0.; - } else if (endBin > (int)fSignalData.size()) { - cout << "TRestRawSignal::CalculateBaseLine. Error! Baseline range exceeds the rawdata depth!!" - << endl; - endBin = fSignalData.size(); - } else { - vector::const_iterator first = fSignalData.begin() + startBin; - vector::const_iterator last = fSignalData.begin() + endBin; - vector v(first, last); - const Short_t* signalInRange = &v[0]; - fBaseLine = TMath::Median(endBin - startBin, signalInRange); - } -} - /////////////////////////////////////////////// /// \brief This method calculates the average and fluctuation of the baseline in the /// specified range and writes the values to fBaseLine and fBaseLineSigma respectively. From 76b8d0dc97708541c3550d2d3c6be93fce7cfa99 Mon Sep 17 00:00:00 2001 From: juanan Date: Fri, 31 Mar 2023 11:41:02 +0200 Subject: [PATCH 07/13] Addressing pipeline failure --- src/TRestRawSignal.cxx | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 36e389ae..ba3d816a 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -375,7 +375,12 @@ Double_t TRestRawSignal::GetMaxPeakValue() { return GetData(GetMaxPeakBin()); } /////////////////////////////////////////////// /// \brief It returns the bin at which the maximum peak amplitude happens /// -Int_t TRestRawSignal::GetMaxPeakBin() { return TRestSignalAnalysis::GetMaxBin(fSignalData); } +Int_t TRestRawSignal::GetMaxPeakBin() { + if (fRange.Y() <= 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); + if (fRange.X() < 0) fRange.SetX(0); + + return TRestSignalAnalysis::GetMaxBin(fSignalData, fRange.X(), fRange.Y()); +} /////////////////////////////////////////////// /// \brief It returns the amplitude of the signal minimum, baseline will be From f805aac95cfc57c888ac8691577db51049b1f694 Mon Sep 17 00:00:00 2001 From: juanan Date: Fri, 31 Mar 2023 23:35:06 +0200 Subject: [PATCH 08/13] Updating points over threshold calculation after changes in TRestRawSignalAnalysis --- src/TRestRawSignal.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index ba3d816a..c7b5ed33 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -236,8 +236,8 @@ void TRestRawSignal::IncreaseBinBy(Int_t bin, Double_t data) { /// void TRestRawSignal::InitializePointsOverThreshold(const TVector2& thrPar, Int_t nPointsOver, Int_t nPointsFlat) { - fPointsOverThreshold = TRestSignalAnalysis::GetPointsOverThreshold( - fSignalData, fRange, thrPar, nPointsOver, nPointsFlat, fBaseLine, fBaseLineSigma); + fPointsOverThreshold = TRestSignalAnalysis::GetPointsOverThreshold(GetData(), fRange, thrPar, nPointsOver, + nPointsFlat, fBaseLineSigma); CalculateThresholdIntegral(); } From 9fc02a325096c3c8e8abac9e9722d945538bc905 Mon Sep 17 00:00:00 2001 From: juanan Date: Sat, 1 Apr 2023 19:15:54 +0200 Subject: [PATCH 09/13] Minor update --- src/TRestRawSignalEvent.cxx | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/TRestRawSignalEvent.cxx b/src/TRestRawSignalEvent.cxx index b6b8a229..57d444d8 100644 --- a/src/TRestRawSignalEvent.cxx +++ b/src/TRestRawSignalEvent.cxx @@ -764,13 +764,12 @@ TPad* TRestRawSignalEvent::DrawSignal(Int_t signalID, TString option) { gr3->SetLineWidth(2); gr3->SetLineColor(3); Int_t point = 0; - Int_t nPoints = pOver.size(); for (const auto& [index, data] : pOver) { gr3->SetPoint(point, index, data); point++; } - if (nPoints > 0) gr3->Draw("CP"); + if (!pOver.empty()) gr3->Draw("CP"); return fPad; } From 6027b916e9491f9b57ddf2b0714d239477149560 Mon Sep 17 00:00:00 2001 From: juanan Date: Mon, 3 Apr 2023 20:29:52 +0200 Subject: [PATCH 10/13] Renaming TRestSignalAnalysis to TRestPulseShapeAnalysis --- src/TRestRawSignal.cxx | 44 +++++++++++++++++++++--------------------- 1 file changed, 22 insertions(+), 22 deletions(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index c7b5ed33..c252d272 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -67,7 +67,7 @@ #include #include -#include "TRestSignalAnalysis.h" +#include "TRestPulseShapeAnalysis.h" using namespace std; @@ -236,8 +236,8 @@ void TRestRawSignal::IncreaseBinBy(Int_t bin, Double_t data) { /// void TRestRawSignal::InitializePointsOverThreshold(const TVector2& thrPar, Int_t nPointsOver, Int_t nPointsFlat) { - fPointsOverThreshold = TRestSignalAnalysis::GetPointsOverThreshold(GetData(), fRange, thrPar, nPointsOver, - nPointsFlat, fBaseLineSigma); + fPointsOverThreshold = TRestPulseShapeAnalysis::GetPointsOverThreshold( + GetData(), fRange, thrPar, nPointsOver, nPointsFlat, fBaseLineSigma); CalculateThresholdIntegral(); } @@ -269,7 +269,7 @@ Double_t TRestRawSignal::GetIntegral() { if (fRange.X() < 0) fRange.SetX(0); if (fRange.Y() <= 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); - return TRestSignalAnalysis::GetIntegral(GetData(), fRange.X(), fRange.Y()); + return TRestPulseShapeAnalysis::GetIntegral(GetData(), fRange.X(), fRange.Y()); } /////////////////////////////////////////////// @@ -277,7 +277,7 @@ Double_t TRestRawSignal::GetIntegral() { /// by (startBin,endBin). /// Double_t TRestRawSignal::GetIntegralInRange(Int_t startBin, Int_t endBin) { - return TRestSignalAnalysis::GetIntegral(fSignalData, startBin, endBin); + return TRestPulseShapeAnalysis::GetIntegral(fSignalData, startBin, endBin); } /////////////////////////////////////////////// @@ -309,7 +309,7 @@ Double_t TRestRawSignal::GetSlopeIntegral() { << endl; } - return TRestSignalAnalysis::GetSlopeIntegral(fPointsOverThreshold); + return TRestPulseShapeAnalysis::GetSlopeIntegral(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -324,7 +324,7 @@ Double_t TRestRawSignal::GetRiseSlope() { "InitializePointsOverThreshold should be called first." << endl; - return TRestSignalAnalysis::GetRiseSlope(fPointsOverThreshold); + return TRestPulseShapeAnalysis::GetRiseSlope(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -337,7 +337,7 @@ Int_t TRestRawSignal::GetRiseTime() { "InitializePointsOverThreshold should be called first." << endl; - return TRestSignalAnalysis::GetRiseTime(fPointsOverThreshold); + return TRestPulseShapeAnalysis::GetRiseTime(fPointsOverThreshold); } /////////////////////////////////////////////// @@ -346,7 +346,7 @@ Int_t TRestRawSignal::GetRiseTime() { /// Double_t TRestRawSignal::GetTripleMaxIntegral() { auto gr = GetGraph(); - return TRestSignalAnalysis::GetTripleMaxIntegral(gr); + return TRestPulseShapeAnalysis::GetTripleMaxIntegral(gr); } /////////////////////////////////////////////// @@ -356,14 +356,14 @@ Double_t TRestRawSignal::GetTripleMaxIntegral() { Double_t TRestRawSignal::GetAverageInRange(Int_t startBin, Int_t endBin) { if (startBin < 0) startBin = 0; if (endBin <= 0 || endBin > GetNumberOfPoints()) endBin = GetNumberOfPoints(); - return TRestSignalAnalysis::GetAverage(GetData(), startBin, endBin); + return TRestPulseShapeAnalysis::GetAverage(GetData(), startBin, endBin); } /////////////////////////////////////////////// /// \brief It returns the temporal width of the peak with maximum amplitude /// inside the signal /// -Int_t TRestRawSignal::GetMaxPeakWidth() { return TRestSignalAnalysis::GetMaxPeakWidth(GetData()); } +Int_t TRestRawSignal::GetMaxPeakWidth() { return TRestPulseShapeAnalysis::GetMaxPeakWidth(GetData()); } /////////////////////////////////////////////// /// \brief It returns the amplitude of the signal maximum, baseline will be @@ -379,7 +379,7 @@ Int_t TRestRawSignal::GetMaxPeakBin() { if (fRange.Y() <= 0 || fRange.Y() > GetNumberOfPoints()) fRange.SetY(GetNumberOfPoints()); if (fRange.X() < 0) fRange.SetX(0); - return TRestSignalAnalysis::GetMaxBin(fSignalData, fRange.X(), fRange.Y()); + return TRestPulseShapeAnalysis::GetMaxBin(fSignalData, fRange.X(), fRange.Y()); } /////////////////////////////////////////////// @@ -392,7 +392,7 @@ Double_t TRestRawSignal::GetMinPeakValue() { return GetData(GetMinPeakBin()); } /////////////////////////////////////////////// /// \brief It returns the bin at which the minimum peak amplitude happens /// -Int_t TRestRawSignal::GetMinPeakBin() { return TRestSignalAnalysis::GetMinBin(fSignalData); } +Int_t TRestRawSignal::GetMinPeakBin() { return TRestPulseShapeAnalysis::GetMinBin(fSignalData); } /////////////////////////////////////////////// /// \brief It returns whether the signal has ADC saturation @@ -472,10 +472,10 @@ void TRestRawSignal::GetWhiteNoiseSignal(TRestRawSignal* noiseSignal, Double_t n /// std::vector TRestRawSignal::GetSignalSmoothed(Int_t averagingPoints, std::string option) { if (option == "") { - return TRestSignalAnalysis::GetSignalSmoothed(fSignalData, averagingPoints); + return TRestPulseShapeAnalysis::GetSignalSmoothed(fSignalData, averagingPoints); } else if (ToUpper(option) == "EXCLUDE OUTLIERS") { - return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, fBaseLine, - fBaseLineSigma); + return TRestPulseShapeAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, + fBaseLine, fBaseLineSigma); } else { std::cout << "TRestRawSignal::GetSignalSmoothed. Error! No such option!" << std::endl; std::vector result; @@ -513,8 +513,8 @@ void TRestRawSignal::GetBaseLineCorrected(TRestRawSignal* smoothedSignal, Int_t /// points used to average the signal /// std::vector TRestRawSignal::GetSignalSmoothed_ExcludeOutliers(Int_t averagingPoints) { - return TRestSignalAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, fBaseLine, - fBaseLineSigma); + return TRestPulseShapeAnalysis::GetSignalSmoothed_ExcludeOutliers(fSignalData, averagingPoints, fBaseLine, + fBaseLineSigma); } /////////////////////////////////////////////// @@ -529,11 +529,11 @@ std::vector TRestRawSignal::GetSignalSmoothed_ExcludeOutliers(Int_t ave /// void TRestRawSignal::CalculateBaseLine(Int_t startBin, Int_t endBin, const std::string& option) { if (ToUpper(option) == "ROBUST") { - TRestSignalAnalysis::CalculateBaselineAndSigmaIQR(fSignalData, startBin, endBin, fBaseLine, - fBaseLineSigma); + TRestPulseShapeAnalysis::CalculateBaselineAndSigmaIQR(fSignalData, startBin, endBin, fBaseLine, + fBaseLineSigma); } else { - TRestSignalAnalysis::CalculateBaselineAndSigmaSD(fSignalData, startBin, endBin, fBaseLine, - fBaseLineSigma); + TRestPulseShapeAnalysis::CalculateBaselineAndSigmaSD(fSignalData, startBin, endBin, fBaseLine, + fBaseLineSigma); } } From 1f4521c3db99731fdb0b3b0070c01a3526d36050 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Juan=20Antonio=20Garc=C3=ADa?= <80903717+juanangp@users.noreply.github.com> Date: Tue, 23 May 2023 12:26:35 +0200 Subject: [PATCH 11/13] Update src/TRestRawSignal.cxx Co-authored-by: Luis Antonio Obis Aparicio <35803280+lobis@users.noreply.github.com> --- src/TRestRawSignal.cxx | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index c252d272..0ae69dbd 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -178,7 +178,9 @@ Double_t TRestRawSignal::GetData(Int_t n) const { return (Double_t)fSignalData[n std::vector TRestRawSignal::GetData() const { std::vector signalData(GetNumberOfPoints()); - for (int i = 0; i < GetNumberOfPoints(); i++) signalData[i] = GetData(i); + for (int i = 0; i < GetNumberOfPoints(); i++) { + signalData[i] = GetData(i); + } return signalData; } From 13eb9eff8ad0701eee92bc7acd5fb8a6345393f7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Juan=20Antonio=20Garc=C3=ADa?= <80903717+juanangp@users.noreply.github.com> Date: Tue, 23 May 2023 12:26:43 +0200 Subject: [PATCH 12/13] Update src/TRestRawSignal.cxx Co-authored-by: Luis Antonio Obis Aparicio <35803280+lobis@users.noreply.github.com> --- src/TRestRawSignal.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/TRestRawSignal.cxx b/src/TRestRawSignal.cxx index 0ae69dbd..ce6824b1 100644 --- a/src/TRestRawSignal.cxx +++ b/src/TRestRawSignal.cxx @@ -347,8 +347,8 @@ Int_t TRestRawSignal::GetRiseTime() { /// and its neightbour points. /// Double_t TRestRawSignal::GetTripleMaxIntegral() { - auto gr = GetGraph(); - return TRestPulseShapeAnalysis::GetTripleMaxIntegral(gr); + auto graph = GetGraph(); + return TRestPulseShapeAnalysis::GetTripleMaxIntegral(graph); } /////////////////////////////////////////////// From 5f9ead23f86c113e761a809dabebace8af9232d8 Mon Sep 17 00:00:00 2001 From: juanan Date: Thu, 25 May 2023 13:20:18 +0200 Subject: [PATCH 13/13] =?UTF-8?q?Moving=20GetPointsOverThreshold=20from=20?= =?UTF-8?q?=C2=B4std::pair=C2=B4=20to=20=C2=B4std:?= =?UTF-8?q?:pair=C2=B4?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- inc/TRestRawSignal.h | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/inc/TRestRawSignal.h b/inc/TRestRawSignal.h index 132452fd..342f44d3 100644 --- a/inc/TRestRawSignal.h +++ b/inc/TRestRawSignal.h @@ -52,7 +52,7 @@ class TRestRawSignal : public TObject { TGraph* fGraph; //! /// A std::vector containing the index of points that are identified over threshold. - std::vector > fPointsOverThreshold; //! + std::vector > fPointsOverThreshold; //! /// It stores the integral value obtained from the points identified over threshold. Double_t fThresholdIntegral = -1; //! @@ -82,9 +82,7 @@ class TRestRawSignal : public TObject { inline Int_t GetNumberOfPoints() const { return fSignalData.size(); } /// Returns a std::vector containing the indexes of data points over threshold - inline std::vector > GetPointsOverThreshold() const { - return fPointsOverThreshold; - } + inline auto GetPointsOverThreshold() const { return fPointsOverThreshold; } /// Returns the maximum value found in the data points. It includes baseline correction inline Double_t GetMaxValue() { return GetMaxPeakValue(); }