From 104ec86af4c58ac3074a55abdf95b976dbb70c3d Mon Sep 17 00:00:00 2001
From: Filip Leonarski <filip.leonarski@psi.ch>
Date: Tue, 8 Aug 2023 11:13:14 +0200
Subject: [PATCH] RadialIntegration: Minor improvements

---
 image_analysis/RadialIntegration.cpp | 31 +++++++++++++++++++++-------
 image_analysis/RadialIntegration.h   |  6 +++---
 2 files changed, 26 insertions(+), 11 deletions(-)

diff --git a/image_analysis/RadialIntegration.cpp b/image_analysis/RadialIntegration.cpp
index 73895449..4ec71cdb 100644
--- a/image_analysis/RadialIntegration.cpp
+++ b/image_analysis/RadialIntegration.cpp
@@ -5,25 +5,38 @@
 #include "../common/JFJochException.h"
 
 RadialIntegration::RadialIntegration(const std::vector<uint16_t>& in_mapping, uint32_t in_nbins, uint32_t in_pixel_split) :
-        pixel_to_bin(in_mapping), nbins(in_nbins), sum(in_nbins, 0), count(in_nbins, 0), coeff(in_mapping.size(), 1.0f),
+        pixel_to_bin(in_mapping), nbins(in_nbins), sum(in_nbins, 0), count(in_nbins, 0),
         pixel_split(in_pixel_split)
 {
-    if (pixel_split == 1)
-        pixel_spit_weight = 1.0f;
-    else if (pixel_split == 4) {
-        pixel_spit_weight = 0.25f;
+
+    if (pixel_split == 4) {
         if (pixel_to_bin.size() % 4 != 0)
             throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
                                   "With pixel split of 4 input array must be of size multiple of 4");
-    } else
+    } else if (pixel_split != 1)
         throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
                               "Only pixel split of 1 and 4 allowed at the moment for radial integration");
+
+    coeff = (float *) std::aligned_alloc(64, in_mapping.size() * sizeof(float));
+    for (int i = 0; i < in_mapping.size(); i++)
+        coeff[i] = 1.0f;
 }
 
 RadialIntegration::RadialIntegration(const RadialIntegrationMapping &mapping, uint32_t in_pixel_split) :
         RadialIntegration(mapping.GetPixelToBinMapping(), mapping.GetBinNumber(), in_pixel_split)
 {}
 
+RadialIntegration::~RadialIntegration() {
+    std::free(coeff);
+}
+
+void RadialIntegration::LoadRadialIntegrationCorr(const std::vector<float> &v) {
+    if (v.size() != pixel_to_bin.size() * pixel_split)
+        throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
+                              "Mismatch in size of pixel-to-bin mapping and correction");
+    memcpy(coeff, v.data(), pixel_to_bin.size() * pixel_split * sizeof(float));
+}
+
 void RadialIntegration::Clear() {
     for (auto &i : sum)
         i = 0;
@@ -37,13 +50,15 @@ void RadialIntegration::Process(const int16_t *__restrict data, size_t npixel) {
         throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
                               "Mismatch in size of pixel-to-bin mapping and image");
 
+    auto coeff_aligned = std::assume_aligned<64>(coeff);
+
     if (pixel_split == 1) {
         for (int i = 0; i < npixel; i++) {
             auto value = data[i];
             if ((value > INT16_MIN + 4) && (value < INT16_MAX - 4)) {
                 auto bin = pixel_to_bin[i];
                 if (bin < nbins) {
-                    sum[bin] += coeff[i] * value;
+                    sum[bin] += coeff_aligned[i] * value;
                     count[bin] += 1.0f;
                 }
             }
@@ -55,7 +70,7 @@ void RadialIntegration::Process(const int16_t *__restrict data, size_t npixel) {
                 for (int p = 0; p < 4; p++) {
                     auto bin = pixel_to_bin[i * pixel_split + p];
                     if (bin < nbins) {
-                        sum[bin] += coeff[i] * value * 0.25f;
+                        sum[bin] += coeff_aligned[i * pixel_split + p] * value * 0.25f;
                         count[bin] += 0.25f;
                     }
                 }
diff --git a/image_analysis/RadialIntegration.h b/image_analysis/RadialIntegration.h
index d77ddb88..47e4bd41 100644
--- a/image_analysis/RadialIntegration.h
+++ b/image_analysis/RadialIntegration.h
@@ -13,20 +13,20 @@
 
 class RadialIntegration {
     const std::vector<uint16_t> pixel_to_bin;
-    const std::vector<float> coeff;
+    float *coeff;
     const uint32_t nbins;
     const uint32_t pixel_split;
-    float pixel_spit_weight;
-    float one_over_pixel_spit_weight;
     std::vector<float> sum;
     std::vector<float> count;
 public:
     RadialIntegration(const RadialIntegrationMapping& mapping, uint32_t pixel_split = 1);
     RadialIntegration(const std::vector<uint16_t>& mapping, uint32_t nbins, uint32_t pixel_split = 1);
+    ~RadialIntegration();
     void Clear();
     void Process(const int16_t *data, size_t npixel);
     void ProcessOneImage(const int16_t *data, size_t npixel); // Process + Clear
     void GetResult(std::vector<float> &result) const;
+    void LoadRadialIntegrationCorr(const std::vector<float> &v);
     [[nodiscard]] float GetRangeValue(uint32_t min_bin, uint32_t max_bin);
     [[nodiscard]] const std::vector<float>& GetSum() const;
     [[nodiscard]] const std::vector<float>& GetCount() const;