Analysis Refactoring - Part 2

src/analyser/CMakeLists.txt

@@ -1,7 +1,10 @@
 add_library(
   analyser
+  dataExporter.h
   dataNormaliser1D.cpp
   dataNormaliser1D.h
+  dataNormaliser2D.cpp
+  dataNormaliser2D.h
   dataNormaliser3D.cpp
   dataNormaliser3D.h
   dataNormaliserBase.h

src/analyser/dataExporter.h

@@ -0,0 +1,31 @@
+// SPDX-License-Identifier: GPL-3.0-or-later
+// Copyright (c) 2024 Team Dissolve and contributors
+
+#pragma once
+
+#include "base/processPool.h"
+
+template <typename DataND, typename DataNDExportFileFormat> class DataExporter
+{
+    public:
+    // Try to export the specified data, if a valid filename has been provided
+    static bool exportData(const DataND &targetData, DataNDExportFileFormat &fileAndFormat, const ProcessPool &procPool)
+    {
+        if (fileAndFormat.hasFilename())
+        {
+            if (procPool.isMaster())
+            {
+                if (fileAndFormat.exportData(targetData))
+                    procPool.decideTrue();
+                else
+                {
+                    procPool.decideFalse();
+                    return false;
+                }
+            }
+            else if (!procPool.decision())
+                return false;
+        }
+        return true;
+    }
+};

src/analyser/dataNormaliser1D.cpp

@@ -2,8 +2,6 @@
 // Copyright (c) 2024 Team Dissolve and contributors
 
 #include "analyser/dataNormaliser1D.h"
-#include "expression/expression.h"
-#include "expression/variable.h"
 #include "math/data1D.h"
 #include "math/integrator.h"
 
@@ -14,8 +12,10 @@ void DataNormaliser1D::normalise(NormalisationFunction1D normalisationFunction)
     const auto &xs = targetData_.xAxis();
     auto &values = targetData_.values();
 
+    const auto xDelta = xs.size() > 1 ? xs[1] - xs[0] : 1.0;
+
     for (auto i = 0; i < xs.size(); ++i)
-        values.at(i) = normalisationFunction(xs[i], values.at(i));
+        values.at(i) = normalisationFunction(xs[i], xDelta, values.at(i));
 }
 
 void DataNormaliser1D::normaliseByGrid() { Messenger::warn("Grid normalisation not implemented for 1D data."); }
@@ -57,4 +57,4 @@ void DataNormaliser1D::normaliseTo(double value, bool absolute)
     auto sum = absolute ? Integrator::absSum(targetData_) : Integrator::sum(targetData_);
     targetData_ /= sum;
     targetData_ *= value;
-}
+}

src/analyser/dataNormaliser1D.h

@@ -5,9 +5,8 @@
 
 #include "analyser/dataNormaliserBase.h"
 #include "math/data1D.h"
-#include <string_view>
 
-using NormalisationFunction1D = std::function<double(const double &, const double &)>;
+using NormalisationFunction1D = std::function<double(const double &, const double &, const double &)>;
 class DataNormaliser1D : public DataNormaliserBase<Data1D, NormalisationFunction1D>
 {
     public:

src/analyser/dataNormaliser2D.cpp

@@ -0,0 +1,72 @@
+// SPDX-License-Identifier: GPL-3.0-or-later
+// Copyright (c) 2024 Team Dissolve and contributors
+
+#include "analyser/dataNormaliser2D.h"
+#include "math/data2D.h"
+#include "math/integrator.h"
+
+DataNormaliser2D::DataNormaliser2D(Data2D &targetData) : DataNormaliserBase<Data2D, NormalisationFunction2D>(targetData) {}
+
+void DataNormaliser2D::normalise(NormalisationFunction2D normalisationFunction)
+{
+    const auto &xs = targetData_.xAxis();
+    const auto &ys = targetData_.yAxis();
+    auto &values = targetData_.values();
+
+    const auto xDelta = xs.size() > 1 ? xs[1] - xs[0] : 1.0;
+    const auto yDelta = ys.size() > 1 ? ys[1] - ys[0] : 1.0;
+
+    for (auto i = 0; i < xs.size(); ++i)
+    {
+        for (auto j = 0; j < ys.size(); ++j)
+        {
+            values[{i, j}] = normalisationFunction(xs[i], xDelta, ys[j], yDelta, values[{i, j}]);
+        }
+    }
+}
+
+void DataNormaliser2D::normaliseByGrid() { Messenger::warn("Grid normalisation not implemented for 2D data."); }
+
+void DataNormaliser2D::normaliseBySphericalShell()
+{
+    // We expect x values to be centre-bin values, and regularly spaced
+    const auto &xAxis = targetData_.xAxis();
+    const auto &yAxis = targetData_.yAxis();
+    auto &values = targetData_.values();
+
+    if (xAxis.size() < 2)
+        return;
+
+    // Derive first left-bin boundary from the delta between points 0 and 1
+    auto leftBin = xAxis[0] - (xAxis[1] - xAxis[0]) * 0.5;
+    auto r1Cubed = pow(leftBin, 3);
+    for (auto n = 0; n < xAxis.size(); ++n)
+    {
+        auto r2Cubed = 0.0, rightBin = 0.0;
+        for (auto m = 0; m < yAxis.size(); ++m)
+        {
+            // Get new right-bin from existing left bin boundary and current bin centre
+            auto rightBin = leftBin + 2 * (xAxis[n] - leftBin);
+            auto r2Cubed = pow(rightBin, 3);
+
+            // Calculate divisor for normalisation
+            auto divisor = (4.0 / 3.0) * PI * (r2Cubed - r1Cubed);
+
+            // Peform normalisation step
+            values[{n, m}] /= divisor;
+            if (targetData_.valuesHaveErrors())
+                targetData_.error(n, m) /= divisor;
+        }
+
+        // Overwrite old values
+        r1Cubed = r2Cubed;
+        leftBin = rightBin;
+    }
+}
+
+void DataNormaliser2D::normaliseTo(double value, bool absolute)
+{
+    auto sum = absolute ? Integrator::absSum(targetData_) : Integrator::sum(targetData_);
+    targetData_ /= sum;
+    targetData_ *= value;
+}

src/analyser/dataNormaliser2D.h

@@ -0,0 +1,24 @@
+// SPDX-License-Identifier: GPL-3.0-or-later
+// Copyright (c) 2024 Team Dissolve and contributors
+
+#pragma once
+
+#include "analyser/dataNormaliserBase.h"
+#include "math/data2D.h"
+
+using NormalisationFunction2D =
+    std::function<double(const double &, const double &, const double &, const double &, const double &)>;
+class DataNormaliser2D : public DataNormaliserBase<Data2D, NormalisationFunction2D>
+{
+    public:
+    DataNormaliser2D(Data2D &targetData);
+
+    /*
+     * Normalisation Functions
+     */
+    public:
+    void normalise(NormalisationFunction2D normalisationFunction) override;
+    void normaliseByGrid() override;
+    void normaliseBySphericalShell() override;
+    void normaliseTo(double value = 1.0, bool absolute = true) override;
+};

src/analyser/dataNormaliser3D.h

@@ -5,9 +5,9 @@
 
 #include "analyser/dataNormaliserBase.h"
 #include "math/data3D.h"
-#include <string_view>
 
-using NormalisationFunction3D = std::function<double(const double &, const double &, const double &, const double &)>;
+using NormalisationFunction3D = std::function<double(const double &, const double &, const double &, const double &,
+                                                     const double &, const double &, const double &)>;
 class DataNormaliser3D : public DataNormaliserBase<Data3D, NormalisationFunction3D>
 {
     public:

src/modules/angle/angle.cpp

@@ -8,140 +8,19 @@
 #include "keywords/speciesSiteVector.h"
 #include "keywords/vec3Double.h"
 #include "module/context.h"
-#include "procedure/nodes/calculateAngle.h"
-#include "procedure/nodes/calculateDistance.h"
-#include "procedure/nodes/collect1D.h"
-#include "procedure/nodes/collect2D.h"
-#include "procedure/nodes/collect3D.h"
-#include "procedure/nodes/operateExpression.h"
-#include "procedure/nodes/operateNormalise.h"
-#include "procedure/nodes/operateNumberDensityNormalise.h"
-#include "procedure/nodes/operateSitePopulationNormalise.h"
-#include "procedure/nodes/operateSphericalShellNormalise.h"
-#include "procedure/nodes/process1D.h"
-#include "procedure/nodes/process2D.h"
-#include "procedure/nodes/process3D.h"
-#include "procedure/nodes/select.h"
 
-AngleModule::AngleModule() : Module(ModuleTypes::Angle), analyser_(ProcedureNode::AnalysisContext)
+AngleModule::AngleModule() : Module(ModuleTypes::Angle)
 {
-    // Select: Site 'A'
-    selectA_ = analyser_.createRootNode<SelectProcedureNode>("A");
-    auto &forEachA = selectA_->branch()->get();
-
-    // -- Select: Site 'B'
-    selectB_ = forEachA.create<SelectProcedureNode>("B");
-    auto &forEachB = selectB_->branch()->get();
-
-    // -- -- Calculate: 'rAB'
-    auto calcAB = forEachB.create<CalculateDistanceProcedureNode>("rAB", selectA_, selectB_);
-
-    // -- -- Collect1D:  'RDF(AB)'
-    collectAB_ =
-        forEachB.create<Collect1DProcedureNode>({}, calcAB, ProcedureNode::AnalysisContext, rangeAB_.x, rangeAB_.y, rangeAB_.z);
-
-    // -- -- Select: Site 'C'
-    selectC_ = forEachB.create<SelectProcedureNode>("C");
-    auto &forEachC = selectC_->branch()->get();
-
-    // -- -- -- Calculate: 'rBC'
-    auto calcBC = forEachC.create<CalculateDistanceProcedureNode>({}, selectB_, selectC_);
-
-    // -- -- -- Calculate: 'aABC'
-    calculateAngle_ = forEachC.create<CalculateAngleProcedureNode>({}, selectA_, selectB_, selectC_);
-    calculateAngle_->keywords().set("Symmetric", symmetric_);
-
-    // -- -- -- Collect1D:  'RDF(BC)'
-    collectBC_ =
-        forEachC.create<Collect1DProcedureNode>({}, calcBC, ProcedureNode::AnalysisContext, rangeBC_.x, rangeBC_.y, rangeBC_.z);
-
-    // -- -- -- Collect1D:  'ANGLE(ABC)'
-    collectABC_ = forEachC.create<Collect1DProcedureNode>({}, calculateAngle_, ProcedureNode::AnalysisContext, angleRange_.x,
-                                                          angleRange_.y, angleRange_.z);
-
-    // -- -- -- Collect2D:  'DAngle (A-B)-C'
-    collectDAngleAB_ =
-        forEachC.create<Collect2DProcedureNode>({}, calcAB, calculateAngle_, ProcedureNode::AnalysisContext, rangeAB_.x,
-                                                rangeAB_.y, rangeAB_.z, angleRange_.x, angleRange_.y, angleRange_.z);
-
-    // -- -- -- Collect2D:  'DAngle A-(B-C)'
-    collectDAngleBC_ =
-        forEachC.create<Collect2DProcedureNode>({}, calcBC, calculateAngle_, ProcedureNode::AnalysisContext, rangeBC_.x,
-                                                rangeBC_.y, rangeBC_.z, angleRange_.x, angleRange_.y, angleRange_.z);
-
-    // -- -- -- Collect3D:  'rAB vs rBC vs aABC'
-    collectDDA_ = forEachC.create<Collect3DProcedureNode>({}, calcAB, calcBC, calculateAngle_, ProcedureNode::AnalysisContext,
-                                                          rangeAB_, rangeBC_, angleRange_);
-
-    // Process1D: 'RDF(AB)'
-    processAB_ = analyser_.createRootNode<Process1DProcedureNode>("RDF(AB)", collectAB_);
-    processAB_->keywords().set("LabelValue", std::string("g\\sub{AB}(r)"));
-    processAB_->keywords().set("LabelX", std::string("r, \\symbol{Angstrom}"));
-    auto &rdfABNormalisation = processAB_->branch()->get();
-    rdfABNormalisation.create<OperateSitePopulationNormaliseProcedureNode>({},
-                                                                           ConstNodeVector<SelectProcedureNode>({selectA_}));
-    rdfABNormalisation.create<OperateNumberDensityNormaliseProcedureNode>({}, ConstNodeVector<SelectProcedureNode>({selectB_}));
-    rdfABNormalisation.create<OperateSphericalShellNormaliseProcedureNode>({});
-
-    // Process1D: 'RDF(BC)'
-    processBC_ = analyser_.createRootNode<Process1DProcedureNode>("RDF(BC)", collectBC_);
-    processBC_->keywords().set("LabelValue", std::string("g\\sub{BC}(r)"));
-    processBC_->keywords().set("LabelX", std::string("r, \\symbol{Angstrom}"));
-    auto &rdfBCNormalisation = processBC_->branch()->get();
-    rdfBCNormalisation.create<OperateSitePopulationNormaliseProcedureNode>(
-        {}, ConstNodeVector<SelectProcedureNode>({selectB_, selectA_}));
-    rdfBCNormalisation.create<OperateNumberDensityNormaliseProcedureNode>({}, ConstNodeVector<SelectProcedureNode>({selectC_}));
-    rdfBCNormalisation.create<OperateSphericalShellNormaliseProcedureNode>({});
-
-    // Process1D: 'ANGLE(ABC)'
-    processAngle_ = analyser_.createRootNode<Process1DProcedureNode>("Angle(ABC)", collectABC_);
-    processAngle_->keywords().set("LabelValue", std::string("Normalised Frequency"));
-    processAngle_->keywords().set("LabelX", std::string("\\symbol{theta}, \\symbol{degrees}"));
-    auto &angleNormalisation = processAngle_->branch()->get();
-    angleNormalisation.create<OperateExpressionProcedureNode>({}, "value/sin(toRad(x))");
-    angleNormalisation.create<OperateNormaliseProcedureNode>({}, 1.0);
-
-    // Process2D: 'DAngle (A-B)-C'
-    processDAngleAB_ = analyser_.createRootNode<Process2DProcedureNode>("DAngle((A-B)-C)", collectDAngleAB_);
-    processDAngleAB_->keywords().set("LabelValue", std::string("g\\sub{AB}(r)"));
-    processDAngleAB_->keywords().set("LabelX", std::string("r, \\symbol{Angstrom}"));
-    processDAngleAB_->keywords().set("LabelY", std::string("\\symbol{theta}, \\symbol{degrees}"));
-    auto &dAngleABNormalisation = processDAngleAB_->branch()->get();
-    dAngleABNormalisationExpression_ = dAngleABNormalisation.create<OperateExpressionProcedureNode>(
-        {}, fmt::format("{} * value/sin(toRad(y))/sin(toRad(yDelta))", symmetric_ ? 1.0 : 2.0));
-    dAngleABNormalisation.create<OperateSitePopulationNormaliseProcedureNode>(
-        {}, ConstNodeVector<SelectProcedureNode>({selectA_, selectC_}));
-    dAngleABNormalisation.create<OperateNumberDensityNormaliseProcedureNode>({},
-                                                                             ConstNodeVector<SelectProcedureNode>({selectB_}));
-    dAngleABNormalisation.create<OperateSphericalShellNormaliseProcedureNode>({});
-
-    // Process2D: 'DAngle A-(B-C)'
-    processDAngleBC_ = analyser_.createRootNode<Process2DProcedureNode>("DAngle(A-(B-C))", collectDAngleBC_);
-    processDAngleBC_->keywords().set("LabelValue", std::string("g\\sub{BC}(r)"));
-    processDAngleBC_->keywords().set("LabelX", std::string("r, \\symbol{Angstrom}"));
-    processDAngleBC_->keywords().set("LabelY", std::string("\\symbol{theta}, \\symbol{degrees}"));
-    auto &dAngleBCNormalisation = processDAngleBC_->branch()->get();
-    dAngleBCNormalisationExpression_ = dAngleBCNormalisation.create<OperateExpressionProcedureNode>(
-        {}, fmt::format("{} * value/sin(toRad(y))/sin(toRad(yDelta))", symmetric_ ? 1.0 : 2.0));
-    dAngleBCNormalisation.create<OperateSitePopulationNormaliseProcedureNode>(
-        {}, ConstNodeVector<SelectProcedureNode>({selectB_, selectA_}));
-    dAngleBCNormalisation.create<OperateNumberDensityNormaliseProcedureNode>({},
-                                                                             ConstNodeVector<SelectProcedureNode>({selectC_}));
-    dAngleBCNormalisation.create<OperateSphericalShellNormaliseProcedureNode>({});
-
     /*
      * Keywords
      */
 
     keywords_.addTarget<ConfigurationKeyword>("Configuration", "Set target configuration for the module", targetConfiguration_);
 
     keywords_.setOrganisation("Options", "Sites", "Specify sites defining the angle interaction A-B-C.");
-    keywords_.add<SpeciesSiteVectorKeyword>("SiteA", "Specify site(s) which represent 'A' in the interaction A-B-C",
-                                            selectA_->speciesSites(), selectA_->axesRequired());
-    keywords_.add<SpeciesSiteVectorKeyword>("SiteB", "Specify site(s) which represent 'B' in the interaction A-B-C",
-                                            selectB_->speciesSites(), selectB_->axesRequired());
-    keywords_.add<SpeciesSiteVectorKeyword>("SiteC", "Specify site(s) which represent 'C' in the interaction A-B-C",
-                                            selectC_->speciesSites(), selectC_->axesRequired());
+    keywords_.add<SpeciesSiteVectorKeyword>("SiteA", "Specify site(s) which represent 'A' in the interaction A-B-C", a_);
+    keywords_.add<SpeciesSiteVectorKeyword>("SiteB", "Specify site(s) which represent 'B' in the interaction A-B-C", b_);
+    keywords_.add<SpeciesSiteVectorKeyword>("SiteC", "Specify site(s) which represent 'C' in the interaction A-B-C", c_);
 
     keywords_.setOrganisation("Options", "Ranges", "Ranges over which to bin quantities from the calculation.");
     keywords_.add<Vec3DoubleKeyword>("RangeAB", "Range (min, max, binwidth) of A-B distance binning", rangeAB_,
@@ -167,16 +46,16 @@ AngleModule::AngleModule() : Module(ModuleTypes::Angle), analyser_(ProcedureNode
 
     keywords_.setOrganisation("Export");
     keywords_.add<FileAndFormatKeyword>("ExportAB", "File format and file name under which to save calculated A-B RDF data",
-                                        processAB_->exportFileAndFormat(), "EndExportAB");
+                                        exportFileAndFormatAB_, "EndExportAB");
     keywords_.add<FileAndFormatKeyword>("ExportBC", "File format and file name under which to save calculated B-C RDF data",
-                                        processBC_->exportFileAndFormat(), "EndExportBC");
+                                        exportFileAndFormatBC_, "EndExportBC");
     keywords_.add<FileAndFormatKeyword>("ExportAngle",
                                         "File format and file name under which to save calculated A-B-C angle histogram",
-                                        processAngle_->exportFileAndFormat(), "EndExportAngle");
+                                        exportFileAndFormatAngle_, "EndExportAngle");
     keywords_.add<FileAndFormatKeyword>("ExportDAngleAB",
                                         "File format and file name under which to save calculated (A-B)-C distance-angle map",
-                                        processDAngleAB_->exportFileAndFormat(), "EndExportDAngleAB");
+                                        exportFileAndFormatDAngleAB_, "EndExportDAngleAB");
     keywords_.add<FileAndFormatKeyword>("ExportDAngleBC",
                                         "File format and file name under which to save calculated A-(B-C) distance-angle map",
-                                        processDAngleBC_->exportFileAndFormat(), "EndExportDAngleBC");
+                                        exportFileAndFormatDAngleBC_, "EndExportDAngleBC");
 }