added some more tests

CMakeLists.txt

@@ -471,6 +471,7 @@ if(AARE_TESTS)
             ${CMAKE_CURRENT_SOURCE_DIR}/src/ClusterFinderMT.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/Pedestal.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/JungfrauDataFile.test.cpp
+            ${CMAKE_CURRENT_SOURCE_DIR}/src/MythenFileReader.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/NumpyFile.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/NumpyHelpers.test.cpp
             ${CMAKE_CURRENT_SOURCE_DIR}/src/RawFile.test.cpp

include/aare/AngleCalibration.hpp

@@ -1,6 +1,5 @@
 #pragma once
 #include <algorithm>
-#include <bitset>
 #include <cmath>
 #include <cstdint>
 
@@ -162,13 +161,13 @@ class FlatField {
         : mythen_detector(mythen_detector_) {
 
         flat_field = NDArray<uint32_t, 1>(
-            std::array<ssize_t, 1>{mythen_detector->num_strips()});
+            std::array<ssize_t, 1>{mythen_detector->num_strips()}, 0);
     }
 
     void read_flatfield_from_file(const std::string &filename) {
 
         std::string word;
-        uint32_t module_number{};
+        uint32_t strip_number{};
 
         try {
             std::ifstream file(filename, std::ios_base::in);
@@ -181,11 +180,11 @@ class FlatField {
 
             while (file_buffer >> word) {
 
-                module_number = std::stoi(word);
+                strip_number = std::stoi(word);
 
                 file_buffer >> word;
-                if (!mythen_detector->get_bad_channels()[module_number])
-                    flat_field[module_number] = std::stod(word);
+                if (!mythen_detector->get_bad_channels()[strip_number])
+                    flat_field[strip_number] = std::stod(word);
             }
 
             file.close();
@@ -194,6 +193,8 @@ class FlatField {
         }
     }
 
+    NDView<uint32_t, 1> get_flatfield() { return flat_field.view(); }
+
     double calculate_mean(double tolerance = 0.001) {
         auto [sum, count] = std::accumulate(
             flat_field.begin(), flat_field.end(),
@@ -212,9 +213,20 @@ class FlatField {
 
         NDArray<double, 1> inverse_normalized_flatfield(flat_field.shape());
 
+        /*
         std::transform(flat_field.begin(), flat_field.end(),
                        inverse_normalized_flatfield.begin(),
-                       [&mean](const auto element) { return mean / element; });
+                       [&mean](const auto element) {
+                           return element == 0 ? 0.0 : mean / element;
+                       });
+        */
+
+        for (ssize_t i = 0; i < flat_field.size(); ++i) {
+            inverse_normalized_flatfield[i] =
+                (flat_field[i] == 0 ? 0.0 : mean / flat_field[i]);
+            if (inverse_normalized_flatfield[i] < tolerance)
+                mythen_detector->get_bad_channels()[i] = true;
+        }
 
         return inverse_normalized_flatfield; // TODO: better to have a copy in
                                              // this context but unneccessary
@@ -228,14 +240,20 @@ class FlatField {
 
         NDArray<double, 1> normalized_flatfield(flat_field.shape());
 
+        /*
         std::transform(flat_field.begin(), flat_field.end(),
                        normalized_flatfield.begin(),
                        [&mean](const auto element) { return element / mean; });
+        */
 
+        for (ssize_t i = 0; i < flat_field.size(); ++i) {
+            normalized_flatfield[i] = (flat_field[i] == flat_field[i] / mean);
+            if (normalized_flatfield[i] < tolerance)
+                mythen_detector->get_bad_channels()[i] = true;
+        }
         return normalized_flatfield;
     }
 
-    // TODO: update is bad channels
   private:
     NDArray<uint32_t, 1> flat_field; // TODO: should be 2d
     std::shared_ptr<MythenDetectorSpecifications> mythen_detector;
@@ -277,7 +295,8 @@ class AngleCalibration {
     /** reads the historical Detector Group (DG) parameters from file **/
     void read_initial_calibration_from_file(const std::string &filename);
 
-    /** converts DG parameters to easy EE parameters e.g.geometric parameters */
+    /** converts DG parameters to easy EE parameters e.g.geometric
+     * parameters */
     parameters convert_to_EE_parameters();
 
     std::tuple<double, double, double>
@@ -287,16 +306,16 @@ class AngleCalibration {
      * parameters */
     parameters convert_to_BC_parameters();
 
-    /** calculates diffraction angle from EE module parameters (used in Beer's
-     * Law)
+    /** calculates diffraction angle from EE module parameters (used in
+     * Beer's Law)
      * @param strip_index local strip index of module
      */
     double diffraction_angle_from_EE_parameters(
         const double normal_distance, const double module_center_distance,
         const double angle, const size_t strip_index);
 
-    /** calculates diffraction angle from EE module parameters (used in Beer's
-     * Law)
+    /** calculates diffraction angle from EE module parameters (used in
+     * Beer's Law)
      * @param strip_index local strip index of module
      */
     double diffraction_angle_from_DG_parameters(const double center,
@@ -324,8 +343,8 @@ class AngleCalibration {
                                               const size_t end_frame_index);
 
     /**
-     * redistributes photon counts with of histogram using one bin per strip to
-     * histogram with fixed size angle bins
+     * redistributes photon counts with of histogram using one bin per strip
+     * to histogram with fixed size angle bins
      * @param frame MythenFrame storing data from image
      * @param bin_counts accumulate new photon counts
      * @param new_statistical_weights accumulate new statistical weights
@@ -345,15 +364,15 @@ class AngleCalibration {
     // TODO: check if interpretation and units are correct
     // historical DG parameters
     // TODO change to NDArray
-    std::vector<double> centers; // orthogonal projection of sample onto
-                                 // detector (given in strip number) [mm]
-                                 // D/pitch
+    std::vector<double> centers;     // orthogonal projection of sample onto
+                                     // detector (given in strip number) [mm]
+                                     // D/pitch
+    std::vector<double> conversions; // pitch/(normal distance from sample
+                                     // to detector (R)) [mm]
+                                     // //used for easy conversion
     std::vector<double>
-        conversions; // pitch/(normal distance from sample to detector (R)) [mm]
-                     // //used for easy conversion
-    std::vector<double>
-        offsets; // position of strip zero relative to sample [degrees] phi -
-                 // 180/pi*D/R TODO: expected an arcsin(D/R)?
+        offsets; // position of strip zero relative to sample [degrees] phi
+                 // - 180/pi*D/R TODO: expected an arcsin(D/R)?
 
     std::shared_ptr<MythenDetectorSpecifications> mythen_detector;
 

include/aare/MythenFileReader.hpp

@@ -3,6 +3,7 @@
  * @short minimal file reader to read mythen files
  ***********************************************/
 
+#include <bitset>
 #include <filesystem>
 #include <string>
 
@@ -15,7 +16,7 @@ struct MythenFrame {
     NDArray<uint32_t, 1> photon_counts;
     double detector_angle{};
     // double reference_intensity{}; not needed
-    std::array<bool, 3> channel_mask{};
+    std::array<uint8_t, 3> channel_mask{};
 };
 
 /** minimal version for a mythen file reader */
@@ -30,8 +31,8 @@ class MythenFileReader : public HDF5FileReader {
         // TODO not a good design fixed number of digits in file name for frame
         // number -> pad with zeros
         //  not even sure if files have the same name
-        std::string current_file_name = m_base_path.string() + file_prefix +
-                                        std::to_string(frame_index) + ".h5";
+        std::string current_file_name =
+            m_base_path / (file_prefix + std::to_string(frame_index) + ".h5");
 
         open_file(current_file_name);
 
@@ -63,9 +64,9 @@ class MythenFileReader : public HDF5FileReader {
         // significant
         //  bit is ask Anna again
 
-        std::array<bool, 3> channel_mask{binary_channel_numbers[2],
-                                         binary_channel_numbers[1],
-                                         binary_channel_numbers[0]};
+        std::array<uint8_t, 3> channel_mask{binary_channel_numbers[0],
+                                            binary_channel_numbers[1],
+                                            binary_channel_numbers[2]};
 
         close_file();
 

src/AngleCalibration.cpp

@@ -283,8 +283,11 @@ void AngleCalibration::write_to_file(const std::string &filename) {
 
     output_file << std::fixed << std::setprecision(6);
     for (ssize_t i = 0; i < num_bins; ++i) {
-        output_file << i * histogram_bin_width << " " << new_photon_counts[i]
-                    << " " << new_photon_count_errors[i] << std::endl;
+        if (new_photon_counts[i] == 0)
+            continue;
+        output_file << mythen_detector->min_angle() + i * histogram_bin_width
+                    << " " << new_photon_counts[i] << " "
+                    << new_photon_count_errors[i] << std::endl;
     }
     output_file.close();
 }

src/AngleCalibration.test.cpp

@@ -7,21 +7,23 @@
 
 #include <filesystem>
 
+#include "test_config.hpp"
+
 #include <catch2/catch_all.hpp>
 #include <catch2/catch_test_macros.hpp>
 #include <catch2/matchers/catch_matchers_floating_point.hpp>
 
 using namespace aare;
 
-class AngleCalibrationTestClass : public aare::AngleCalibration {
+class AngleCalibrationTestClass : public AngleCalibration {
 
   public:
     AngleCalibrationTestClass(
         std::shared_ptr<MythenDetectorSpecifications> mythen_detector_,
         std::shared_ptr<FlatField> flat_field_,
         std::shared_ptr<MythenFileReader> mythen_file_reader_)
-        : aare::AngleCalibration(mythen_detector_, flat_field_,
-                                 mythen_file_reader_) {}
+        : AngleCalibration(mythen_detector_, flat_field_, mythen_file_reader_) {
+    }
     ~AngleCalibrationTestClass() = default;
 
     std::vector<double> get_centers() { return centers; }
@@ -32,7 +34,11 @@ class AngleCalibrationTestClass : public aare::AngleCalibration {
 };
 
 TEST_CASE("read initial angle calibration file",
-          "[.anglecalibration][.fileread]") {
+          "[.anglecalibration] [.files]") {
+
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
 
     std::shared_ptr<MythenDetectorSpecifications> mythen_detector_ptr =
         std::make_shared<MythenDetectorSpecifications>();
@@ -41,17 +47,13 @@ TEST_CASE("read initial angle calibration file",
         std::make_shared<FlatField>(mythen_detector_ptr);
 
     std::shared_ptr<MythenFileReader> mythen_file_reader_ptr =
-        std::make_shared<MythenFileReader>(
-            std::filesystem::path{
-                "/home/mazzola/Documents/mythen3tools/beamline/TDATA"},
-            "ang1up_22keV_LaB60p3mm_48M_a_0");
+        std::make_shared<MythenFileReader>(fpath,
+                                           "ang1up_22keV_LaB60p3mm_48M_a_0");
 
     AngleCalibrationTestClass anglecalibration(
         mythen_detector_ptr, flat_field_ptr, mythen_file_reader_ptr);
 
-    std::string filename =
-        "/home/mazzol_a/Documents/mythen3tools/beamline/"
-        "TDATA/Angcal_2E_Feb2023_P29.off"; // TODO change path upload data
+    std::string filename = fpath / "Angcal_2E_Feb2023_P29.off";
 
     REQUIRE(std::filesystem::exists(filename));
 
@@ -71,3 +73,127 @@ TEST_CASE("read initial angle calibration file",
     CHECK(offsets[47] == Catch::Approx(5.8053312));
     CHECK(conversions[27] == Catch::Approx(-0.6581179125e-4));
 }
+
+TEST_CASE("read bad channels",
+          "[.anglecalibration][.mythenspecifications][.files]") {
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
+
+    MythenDetectorSpecifications mythen_detector;
+
+    std::string bad_channels_filename = fpath / "bc2023_003_RING.chans";
+
+    REQUIRE(std::filesystem::exists(bad_channels_filename));
+
+    mythen_detector.read_bad_channels_from_file(bad_channels_filename);
+
+    CHECK(mythen_detector.get_bad_channels().size() == 61440);
+
+    CHECK(mythen_detector.get_bad_channels()[61437] == true);
+    CHECK(std::all_of(mythen_detector.get_bad_channels().begin() + 30720,
+                      mythen_detector.get_bad_channels().begin() + 61439,
+                      [](const bool element) { return element; }));
+}
+
+TEST_CASE("read unconnected modules",
+          "[.anglecalibration][.mythenspecifications][.files]") {
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
+
+    MythenDetectorSpecifications mythen_detector;
+
+    std::string unconnected_modules_filename = fpath / "ModOut.txt";
+
+    REQUIRE(std::filesystem::exists(unconnected_modules_filename));
+
+    mythen_detector.read_unconnected_modules_from_file(
+        unconnected_modules_filename);
+
+    CHECK(mythen_detector.get_connected_modules().size() == 48);
+
+    CHECK(std::all_of(mythen_detector.get_connected_modules().begin(),
+                      mythen_detector.get_connected_modules().end(),
+                      [](const bool element) { return element; }));
+}
+
+TEST_CASE("read flatfield", "[.anglecalibration][.flatfield][.files]") {
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
+
+    std::shared_ptr<MythenDetectorSpecifications> mythen_detector_ptr =
+        std::make_shared<MythenDetectorSpecifications>();
+
+    FlatField flatfield(mythen_detector_ptr);
+
+    std::string flatfield_filename =
+        fpath /
+        "Flatfield_E22p0keV_T11000eV_up_48M_a_LONG_Feb2023_open_WS_SUMC.raw";
+
+    REQUIRE(std::filesystem::exists(flatfield_filename));
+
+    flatfield.read_flatfield_from_file(flatfield_filename);
+
+    auto flatfield_data = flatfield.get_flatfield();
+
+    CHECK(flatfield_data.size() == 61440);
+
+    CHECK(flatfield_data[0] == 0);
+    CHECK(flatfield_data[21] == 4234186);
+}
+
+TEST_CASE("calculate new fixed angle width bins histogram",
+          "[.anglecalibration] [.files]") {
+
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
+
+    std::shared_ptr<MythenDetectorSpecifications> mythen_detector_ptr =
+        std::make_shared<MythenDetectorSpecifications>();
+
+    std::string bad_channels_filename = fpath / "bc2023_003_RING.chans";
+
+    REQUIRE(std::filesystem::exists(bad_channels_filename));
+
+    mythen_detector_ptr->read_bad_channels_from_file(bad_channels_filename);
+
+    std::string unconnected_modules_filename = fpath / "ModOut.txt";
+
+    REQUIRE(std::filesystem::exists(unconnected_modules_filename));
+
+    mythen_detector_ptr->read_unconnected_modules_from_file(
+        unconnected_modules_filename);
+
+    std::shared_ptr<FlatField> flat_field_ptr =
+        std::make_shared<FlatField>(mythen_detector_ptr);
+
+    std::string flatfield_filename =
+        fpath /
+        "Flatfield_E22p0keV_T11000eV_up_48M_a_LONG_Feb2023_open_WS_SUMC.raw";
+
+    REQUIRE(std::filesystem::exists(flatfield_filename));
+
+    flat_field_ptr->read_flatfield_from_file(flatfield_filename);
+
+    std::shared_ptr<MythenFileReader> mythen_file_reader_ptr =
+        std::make_shared<MythenFileReader>(fpath,
+                                           "ang1up_22keV_LaB60p3mm_48M_a_0");
+
+    AngleCalibration anglecalibration(mythen_detector_ptr, flat_field_ptr,
+                                      mythen_file_reader_ptr);
+
+    std::string initial_angles_filename = fpath / "Angcal_2E_Feb2023_P29.off";
+
+    REQUIRE(std::filesystem::exists(initial_angles_filename));
+
+    anglecalibration.read_initial_calibration_from_file(
+        initial_angles_filename);
+
+    anglecalibration.calculate_fixed_bin_angle_width_histogram(320, 340);
+
+    anglecalibration.write_to_file(
+        "cpp_new_photon_counts.xye"); // TODO adjust output path
+}

src/Hdf5FileReader.test.cpp

@@ -5,6 +5,8 @@
 
 #include <filesystem>
 
+#include "test_config.hpp"
+
 #include <H5Cpp.h>
 
 #include "aare/Hdf5FileReader.hpp"
@@ -16,11 +18,11 @@
 
 using namespace aare;
 
-TEST_CASE("read hdf5 file", "[.hdf5file]") {
+TEST_CASE("read hdf5 file", "[.hdf5file][.files]") {
 
     // TODO generalize datasetpath
-    std::string filename = "/home/mazzol_a/Documents/mythen3tools/beamline/"
-                           "TDATA/ang1up_22keV_LaB60p3mm_48M_a_0320.h5";
+    std::string filename = test_data_path() / "AngleCalibration_Test_Data" /
+                           "ang1up_22keV_LaB60p3mm_48M_a_0320.h5";
 
     REQUIRE(std::filesystem::exists(filename));
 

src/MythenFileReader.test.cpp

@@ -0,0 +1,33 @@
+/************************************************
+ * @file MythenFileReader.test.cpp
+ * @short test case for angle calibration class
+ ***********************************************/
+
+#include "aare/MythenFileReader.hpp"
+
+#include <filesystem>
+
+#include "test_config.hpp"
+
+#include <catch2/catch_all.hpp>
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_floating_point.hpp>
+
+using namespace aare;
+
+TEST_CASE("test mythenfile_reader", "[.mythenfilereader][.files]") {
+
+    auto fpath = test_data_path() / "AngleCalibration_Test_Data";
+
+    REQUIRE(std::filesystem::exists(fpath));
+
+        MythenFileReader file_reader(fpath, "ang1up_22keV_LaB60p3mm_48M_a_0");
+
+    auto frame = file_reader.read_frame(320);
+
+    CHECK(frame.detector_angle == 0.99955);
+
+    CHECK(frame.channel_mask == std::array<uint8_t, 3>{0, 0, 1});
+
+    CHECK(frame.photon_counts.size() == 61440);
+}