Merge branch 'master' of github.com:rest-for-physics/framework

* 'master' of github.com:rest-for-physics/framework: (45 commits) Update README.md Add new presentation to README Update README.md Update README.md Adding dimensions to the NuSTAR mask example [pre-commit.ci] auto fixes from pre-commit.com hooks Adding example and new image TRestRadialStrippedMask. Adding an internal region Setting the spacers of the new NuSTAR like optics as an example [pre-commit.ci] auto fixes from pre-commit.com hooks TRestPatternMask. Solving drawing issue TRestPatternMask::fOffset adding comment Testing pipeline Fixing issues when drawing the mask with an offset and/or rotation TRestRadialStrippedMask minor fixes Adding new radial stripped mask TRestSensitivity::ExportCurve is now renormalized by factor and power [pre-commit.ci] auto fixes from pre-commit.com hooks TRestComponentDataSet::fDFRange added to allow reducing the component dataset statistics TRestDataSet::RegenerateTree method added and used when invoking ApplyRange ...
rest-for-physics · Jun 21, 2024 · 187791e · 187791e
2 parents 0131a06 + fc06755
commit 187791e
Show file tree

Hide file tree

Showing 22 changed files with 629 additions and 114 deletions.
diff --git a/README.md b/README.md
@@ -1,6 +1,6 @@
 # The REST Framework
 [![DOI](https://zenodo.org/badge/324291710.svg)](http://doi.org/10.5281/zenodo.4528985)
-[![pipeline status](https://gitlab.cern.ch/rest-for-physics/framework/badges/master/pipeline.svg)](https://gitlab.cern.ch/rest-for-physics/framework/-/commits/master)
+[![pipeline status](https://github.com/rest-for-physics/framework/actions/workflows/validation.yml/badge.svg)](https://github.com/rest-for-physics/framework)
 [![website](https://img.shields.io/badge/user-guide-E8B6FF.svg)](https://rest-for-physics.github.io)
 [![api](https://img.shields.io/badge/user-API-FFCA78.svg)](https://rest-for-physics.github.io/framework/)
 [![forum](https://img.shields.io/badge/user-forum-AAFF90.svg)](https://rest-forum.unizar.es/)
@@ -77,6 +77,8 @@ A major change at 2.3 will prevent from backwards compatibility, since class nam
 - PandaX-III: Searching for neutrinoless double beta decay with high pressure 136Xe gas time projection chambers. [X. Chen et al., Science China Physics, Mechanics & Astronomy 60, 061011 (2017)](https://doi.org/10.1007/s11433-017-9028-0) [arXiv:1610.08883](https://arxiv.org/abs/1610.08883).
 
 ## Presentations
+- Background model studies for IAXO using the REST-for-Physics framework, Luis Obis, [2024-Apr, VIEnna Workshop on Simulations 2024  (VIEWS24), Vienna](https://indico.cern.ch/event/1275551/contributions/5858816/).
+- REST-for-Physics framework for Geant4 simulations and data analysis, Álvaro Ezquerro, [2024-Apr, VIEnna Workshop on Simulations 2024  (VIEWS24), Vienna](https://indico.cern.ch/event/1275551/contributions/5858814/).
 - REST-for-Physics, Luis Obis, [2022-May, ROOT Users Workshop, FermiLab](https://indico.fnal.gov/event/23628/contributions/240755/).
 - REST v2.0 : A data analysis and simulation framework for micro-patterned readout detectors., Javier Galan, [2016-Dec, 8th Symposium on Large TPCs for low-energy rare event detection, Paris](https://indico.cern.ch/event/473362/contributions/2334838/).
 

diff --git a/data/distributions/CosmicsCry.root b/data/distributions/CosmicsCry.root
diff --git a/doc/doxygen/images/radialstrippedmask.png b/doc/doxygen/images/radialstrippedmask.png
diff --git a/examples/masks.rml b/examples/masks.rml
@@ -14,6 +14,17 @@
         <parameter name="stripsGap" value="12mm"/>
         <parameter name="stripsThickness" value="1.2mm"/>
     </TRestStrippedMask>
+	<TRestRadialStrippedMask name="radialStrips" verboseLevel="warning">
+		<parameter name="maskRadius" value="19.4cm" />
+		<parameter name="offset" value="(0,0)cm" />
+		<parameter name="rotationAngle" value="30degrees" />
+
+		<parameter name="initialRadius" value="5.4cm" />
+		<parameter name="internalRegionRadius" value="0cm" />
+
+		<parameter name="stripsAngle" value="60degrees" />
+		<parameter name="stripsThickness" value="0.7cm" />
+	</TRestRadialStrippedMask>
     <TRestRingsMask name="rings" verboseLevel="warning">
         <parameter name="maskRadius" value="9cm"/>
         <parameter name="offset" value="(3,3)mm"/>

diff --git a/source/framework/core/inc/TRestDataSet.h b/source/framework/core/inc/TRestDataSet.h
@@ -112,7 +112,7 @@ class TRestDataSet : public TRestMetadata {
     Bool_t fExternal = false;  //<
 
     /// The resulting RDF::RNode object after initialization
-    ROOT::RDF::RNode fDataSet = ROOT::RDataFrame(0);  //!
+    ROOT::RDF::RNode fDataFrame = ROOT::RDataFrame(0);  //!
 
     /// A pointer to the generated tree
     TChain* fTree = nullptr;  //!
@@ -122,12 +122,14 @@ class TRestDataSet : public TRestMetadata {
    protected:
     virtual std::vector<std::string> FileSelection();
 
+    void RegenerateTree(std::vector<std::string> finalList = {});
+
    public:
     /// Gives access to the RDataFrame
     ROOT::RDF::RNode GetDataFrame() const {
         if (!fExternal && fTree == nullptr)
             RESTWarning << "DataFrame has not been yet initialized" << RESTendl;
-        return fDataSet;
+        return fDataFrame;
     }
 
     void EnableMultiThreading(Bool_t enable = true) { fMT = enable; }
@@ -152,7 +154,7 @@ class TRestDataSet : public TRestMetadata {
     }
 
     /// Number of variables (or observables)
-    size_t GetNumberOfColumns() { return fDataSet.GetColumnNames().size(); }
+    size_t GetNumberOfColumns() { return fDataFrame.GetColumnNames().size(); }
 
     /// Number of variables (or observables)
     size_t GetNumberOfBranches() { return GetNumberOfColumns(); }
@@ -187,7 +189,7 @@ class TRestDataSet : public TRestMetadata {
 
     void SetTotalTimeInSeconds(Double_t seconds) { fTotalDuration = seconds; }
     void SetDataFrame(const ROOT::RDF::RNode& dS) {
-        fDataSet = dS;
+        fDataFrame = dS;
         fExternal = true;
     }
 
@@ -198,8 +200,12 @@ class TRestDataSet : public TRestMetadata {
     void Export(const std::string& filename, std::vector<std::string> excludeColumns = {});
 
     ROOT::RDF::RNode MakeCut(const TRestCut* cut);
+    ROOT::RDF::RNode ApplyRange(size_t from, size_t to);
+    ROOT::RDF::RNode Range(size_t from, size_t to);
     ROOT::RDF::RNode DefineColumn(const std::string& columnName, const std::string& formula);
 
+    size_t GetEntries();
+
     void PrintMetadata() override;
     void Initialize() override;
 
@@ -209,6 +215,6 @@ class TRestDataSet : public TRestMetadata {
     TRestDataSet(const char* cfgFileName, const std::string& name = "");
     ~TRestDataSet();
 
-    ClassDefOverride(TRestDataSet, 7);
+    ClassDefOverride(TRestDataSet, 8);
 };
 #endif
diff --git a/source/framework/core/src/TRestDataSet.cxx b/source/framework/core/src/TRestDataSet.cxx
@@ -382,30 +382,40 @@ void TRestDataSet::GenerateDataSet() {
         ROOT::DisableImplicitMT();
 
     RESTInfo << "Initializing dataset" << RESTendl;
-    fDataSet = ROOT::RDataFrame("AnalysisTree", fFileSelection);
+    fDataFrame = ROOT::RDataFrame("AnalysisTree", fFileSelection);
 
     RESTInfo << "Making cuts" << RESTendl;
-    fDataSet = MakeCut(fCut);
+    fDataFrame = MakeCut(fCut);
 
     // Adding new user columns added to the dataset
     for (const auto& [cName, cExpression] : fColumnNameExpressions) {
         RESTInfo << "Adding column to dataset: " << cName << RESTendl;
         finalList.emplace_back(cName);
-        fDataSet = DefineColumn(cName, cExpression);
+        fDataFrame = DefineColumn(cName, cExpression);
     }
 
+    RegenerateTree(finalList);
+
+    RESTInfo << " - Dataset generated!" << RESTendl;
+}
+
+///////////////////////////////////////////////
+/// \brief It regenerates the tree so that it is an exact copy of the present DataFrame
+///
+void TRestDataSet::RegenerateTree(std::vector<std::string> finalList) {
     RESTInfo << "Generating snapshot." << RESTendl;
     std::string user = getenv("USER");
     std::string fOutName = "/tmp/rest_output_" + user + ".root";
-    fDataSet.Snapshot("AnalysisTree", fOutName, finalList);
+    if (!finalList.empty())
+        fDataFrame.Snapshot("AnalysisTree", fOutName, finalList);
+    else
+        fDataFrame.Snapshot("AnalysisTree", fOutName);
 
     RESTInfo << "Re-importing analysis tree." << RESTendl;
-    fDataSet = ROOT::RDataFrame("AnalysisTree", fOutName);
+    fDataFrame = ROOT::RDataFrame("AnalysisTree", fOutName);
 
     TFile* f = TFile::Open(fOutName.c_str());
     fTree = (TChain*)f->Get("AnalysisTree");
-
-    RESTInfo << " - Dataset generated!" << RESTendl;
 }
 
 ///////////////////////////////////////////////
@@ -517,14 +527,32 @@ std::vector<std::string> TRestDataSet::FileSelection() {
     return fFileSelection;
 }
 
+///////////////////////////////////////////////
+/// \brief This method returns a RDataFrame node with the number of
+/// samples inside the dataset by selecting a range. It will not
+/// modify internally the dataset. See ApplyRange to modify internally
+/// the dataset.
+///
+ROOT::RDF::RNode TRestDataSet::Range(size_t from, size_t to) { return fDataFrame.Range(from, to); }
+
+///////////////////////////////////////////////
+/// \brief This method reduces the number of samples inside the
+/// dataset by selecting a range.
+///
+ROOT::RDF::RNode TRestDataSet::ApplyRange(size_t from, size_t to) {
+    fDataFrame = fDataFrame.Range(from, to);
+    RegenerateTree();
+    return fDataFrame;
+}
+
 ///////////////////////////////////////////////
 /// \brief This function applies a TRestCut to the dataframe
 /// and returns a dataframe with the applied cuts. Note that
 /// the cuts are not applied directly to the dataframe on
-/// TRestDataSet, to do so you should do fDataSet = MakeCut(fCut);
+/// TRestDataSet, to do so you should do fDataFrame = MakeCut(fCut);
 ///
 ROOT::RDF::RNode TRestDataSet::MakeCut(const TRestCut* cut) {
-    auto df = fDataSet;
+    auto df = fDataFrame;
 
     if (cut == nullptr) return df;
 
@@ -561,6 +589,20 @@ ROOT::RDF::RNode TRestDataSet::MakeCut(const TRestCut* cut) {
     return df;
 }
 
+///////////////////////////////////////////////
+/// \brief It returns the number of entries found inside fDataFrame
+/// and prints out a warning if the number of entries inside the
+/// tree is not the same.
+///
+size_t TRestDataSet::GetEntries() {
+    auto nEntries = fDataFrame.Count();
+    if (*nEntries == (long long unsigned int)GetTree()->GetEntries()) return *nEntries;
+    RESTWarning << "TRestDataSet::GetEntries. Number of tree entries is not the same as RDataFrame entries."
+                << RESTendl;
+    RESTWarning << "Returning RDataFrame entries" << RESTendl;
+    return *nEntries;
+}
+
 ///////////////////////////////////////////////
 /// \brief This function will add a new column to the RDataFrame using
 /// the same scheme as the usual RDF::Define method, but it will on top of
@@ -574,7 +616,7 @@ ROOT::RDF::RNode TRestDataSet::MakeCut(const TRestCut* cut) {
 /// \endcode
 ///
 ROOT::RDF::RNode TRestDataSet::DefineColumn(const std::string& columnName, const std::string& formula) {
-    auto df = fDataSet;
+    auto df = fDataFrame;
 
     std::string evalFormula = formula;
     for (auto const& [name, properties] : fQuantity)
@@ -819,7 +861,7 @@ void TRestDataSet::InitFromConfigFile() {
 void TRestDataSet::Export(const std::string& filename, std::vector<std::string> excludeColumns) {
     RESTInfo << "Exporting dataset" << RESTendl;
 
-    std::vector<std::string> columns = fDataSet.GetColumnNames();
+    std::vector<std::string> columns = fDataFrame.GetColumnNames();
     if (!excludeColumns.empty()) {
         columns.erase(std::remove_if(columns.begin(), columns.end(),
                                      [&excludeColumns](std::string elem) {
@@ -831,10 +873,10 @@ void TRestDataSet::Export(const std::string& filename, std::vector<std::string>
         RESTInfo << "Re-Generating snapshot." << RESTendl;
         std::string user = getenv("USER");
         std::string fOutName = "/tmp/rest_output_" + user + ".root";
-        fDataSet.Snapshot("AnalysisTree", fOutName, columns);
+        fDataFrame.Snapshot("AnalysisTree", fOutName, columns);
 
         RESTInfo << "Re-importing analysis tree." << RESTendl;
-        fDataSet = ROOT::RDataFrame("AnalysisTree", fOutName);
+        fDataFrame = ROOT::RDataFrame("AnalysisTree", fOutName);
 
         TFile* f = TFile::Open(fOutName.c_str());
         fTree = (TChain*)f->Get("AnalysisTree");
@@ -846,7 +888,7 @@ void TRestDataSet::Export(const std::string& filename, std::vector<std::string>
             RESTInfo << "Re-Generating snapshot." << RESTendl;
             std::string user = getenv("USER");
             std::string fOutName = "/tmp/rest_output_" + user + ".root";
-            fDataSet.Snapshot("AnalysisTree", fOutName);
+            fDataFrame.Snapshot("AnalysisTree", fOutName);
 
             TFile* f = TFile::Open(fOutName.c_str());
             fTree = (TChain*)f->Get("AnalysisTree");
@@ -910,7 +952,7 @@ void TRestDataSet::Export(const std::string& filename, std::vector<std::string>
         fprintf(f, "###\n");
         fprintf(f, "### Data starts here\n");
 
-        auto obsNames = fDataSet.GetColumnNames();
+        auto obsNames = fDataFrame.GetColumnNames();
         std::string obsListStr = "";
         for (const auto& l : obsNames) {
             if (!obsListStr.empty()) obsListStr += ":";
@@ -938,7 +980,7 @@ void TRestDataSet::Export(const std::string& filename, std::vector<std::string>
 
         return;
     } else if (TRestTools::GetFileNameExtension(filename) == "root") {
-        fDataSet.Snapshot("AnalysisTree", filename);
+        fDataFrame.Snapshot("AnalysisTree", filename);
 
         TFile* f = TFile::Open(filename.c_str(), "UPDATE");
         std::string name = this->GetName();
@@ -1038,7 +1080,7 @@ void TRestDataSet::Import(const std::string& fileName) {
     else
         ROOT::DisableImplicitMT();
 
-    fDataSet = ROOT::RDataFrame("AnalysisTree", fileName);
+    fDataFrame = ROOT::RDataFrame("AnalysisTree", fileName);
 
     fTree = (TChain*)file->Get("AnalysisTree");
 }
@@ -1104,7 +1146,7 @@ void TRestDataSet::Import(std::vector<std::string> fileNames) {
     }
 
     RESTInfo << "Opening list of files. First file: " << fileNames[0] << RESTendl;
-    fDataSet = ROOT::RDataFrame("AnalysisTree", fileNames);
+    fDataFrame = ROOT::RDataFrame("AnalysisTree", fileNames);
 
     if (fTree != nullptr) {
         delete fTree;

diff --git a/source/framework/masks/inc/TRestPatternMask.h b/source/framework/masks/inc/TRestPatternMask.h
@@ -29,7 +29,7 @@
 /// An abstract class used to encapsulate different mask pattern class definitions.
 class TRestPatternMask : public TRestMetadata {
    private:
-    /// It is used to introduce an offset on the pattern
+    /// It is used to introduce an offset on the pattern (not the mask, mask is always centered)
     TVector2 fOffset = TVector2(0, 0);  //<
 
     /// An angle (in radians) used to introduce a rotation to the pattern

diff --git a/source/framework/masks/inc/TRestRadialStrippedMask.h b/source/framework/masks/inc/TRestRadialStrippedMask.h
@@ -0,0 +1,70 @@
+/*************************************************************************
+ * This file is part of the REST software framework.                     *
+ *                                                                       *
+ * Copyright (C) 2016 GIFNA/TREX (University of Zaragoza)                *
+ * For more information see https://gifna.unizar.es/trex                 *
+ *                                                                       *
+ * REST is free software: you can redistribute it and/or modify          *
+ * it under the terms of the GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3 of the License, or     *
+ * (at your option) any later version.                                   *
+ *                                                                       *
+ * REST is distributed in the hope that it will be useful,               *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of        *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the          *
+ * GNU General Public License for more details.                          *
+ *                                                                       *
+ * You should have a copy of the GNU General Public License along with   *
+ * REST in $REST_PATH/LICENSE.                                           *
+ * If not, see https://www.gnu.org/licenses/.                            *
+ * For the list of contributors see $REST_PATH/CREDITS.                  *
+ *************************************************************************/
+
+#ifndef REST_TRestRadialStrippedMask
+#define REST_TRestRadialStrippedMask
+
+#include <TRestPatternMask.h>
+
+/// A class used to define a stripped mask pattern
+class TRestRadialStrippedMask : public TRestPatternMask {
+   private:
+    void Initialize() override;
+
+    /// The periodity of the stripped structure in radians
+    Double_t fStripsAngle = TMath::Pi() / 3;  //<
+
+    /// The width of the stripped structure in mm
+    Double_t fStripsThickness = 0.5;  //<
+
+    /// The spacers structure will be effective from this radius, in mm. Default is from 20 mm.
+    Double_t fInitialRadius = 20.;  //<
+
+    /// Radius of an internal circular region defined inside the fInitialRadius. If 0, there will be no region
+    Double_t fInternalRegionRadius = 0.;  //<
+
+    /// It defines the maximum number of cells/regions in each axis
+    Int_t fModulus = 10;
+
+   public:
+    virtual Int_t GetRegion(Double_t& x, Double_t& y) override;
+
+    /// It returns the gap/periodicity of the strips in degrees
+    Double_t GetStripsAngle() { return fStripsAngle * units("degrees"); }
+
+    /// It returns the thickness of the strips in mm
+    Double_t GetStripsThickness() { return fStripsThickness; }
+
+    /// It returns the modulus used to define a finite set of ids
+    Int_t GetModulus() { return fModulus; }
+
+    void PrintMetadata() override;
+    void PrintMaskMembers() override;
+    void PrintMask() override;
+
+    TRestRadialStrippedMask();
+    TRestRadialStrippedMask(const char* cfgFileName, std::string name = "");
+    ~TRestRadialStrippedMask();
+
+    ClassDefOverride(TRestRadialStrippedMask, 1);
+};
+#endif
diff --git a/source/framework/masks/src/TRestPatternMask.cxx b/source/framework/masks/src/TRestPatternMask.cxx
@@ -167,7 +167,7 @@ TCanvas* TRestPatternMask::DrawMonteCarlo(Int_t nSamples) {
         delete fCanvas;
         fCanvas = NULL;
     }
-    fCanvas = new TCanvas("canv", "This is the canvas title", 1400, 1200);
+    fCanvas = new TCanvas("canv", "This is the canvas title", 1200, 1200);
     fCanvas->Draw();
 
     TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
@@ -183,17 +183,19 @@ TCanvas* TRestPatternMask::DrawMonteCarlo(Int_t nSamples) {
     TRandom3* rnd = new TRandom3(0);
 
     for (int n = 0; n < nSamples; n++) {
-        Double_t x = 2.5 * (rnd->Rndm() - 0.5) * fMaskRadius;
-        Double_t y = 2.5 * (rnd->Rndm() - 0.5) * fMaskRadius;
+        Double_t xO = 2.5 * (rnd->Rndm() - 0.5) * fMaskRadius;
+        Double_t yO = 2.5 * (rnd->Rndm() - 0.5) * fMaskRadius;
+        Double_t x = xO;
+        Double_t y = yO;
 
         Int_t id = GetRegion(x, y);
 
         if (points.count(id) == 0) {
             std::vector<TVector2> a;
-            a.push_back(TVector2(x, y));
+            a.push_back(TVector2(xO, yO));
             points[id] = a;
         } else {
-            points[id].push_back(TVector2(x, y));
+            points[id].push_back(TVector2(xO, yO));
         }
     }