final changes to bin_pixels code and mex code

CMakeLists.txt

@@ -1,4 +1,4 @@
-cmake_minimum_required(VERSION 3.15)
+cmake_minimum_required(VERSION 3.29)
 
 set(Horace_ROOT ${CMAKE_CURRENT_LIST_DIR})
 file(READ "${Horace_ROOT}/VERSION" _version)
@@ -66,7 +66,7 @@ include(PACE_FindMatlab)
 #
 message(" MATLAB CURRENT VERSION: ${MATLAB_CURRENT_VERSION}")
 set(BUILD_HDF_MEX_PLUGIN TRUE)
-if ((${MATLAB_CURRENT_VERSION} VERSION_LESS "7.11") OR (${MATLAB_CURRENT_VERSION} VERSION_GREATER "9.4"))
+if ((${MATLAB_CURRENT_VERSION} VERSION_LESS_EQUAL "7.11") OR (${MATLAB_CURRENT_VERSION} VERSION_GREATER "9.4"))
     set(BUILD_HDF_MEX_PLUGIN FALSE)
 endif()
 
@@ -75,7 +75,7 @@ include(PACE_AddMex)
 include(PACE_CodeAnalysis)
 include(PACE_Docs)
 if (${BUILD_HDF_MEX_PLUGIN})
-	include(horace_FindHDF5)
+    include(horace_FindHDF5)
 endif()
 
 include(herbert_FindMPI)

_LowLevelCode/cpp/accumulate_cut_c/accumulate_cut_c.h

@@ -36,7 +36,7 @@ mwSize accumulate_cut(double* s, double* e, double* npix,
 
     T Inf(0);
     double ebin_inv = (1 / ebin);
-    bool ignore_something, ignote_all;
+    bool ignore_something, ignore_all;
 
     //if we want to ignore nan and inf in the data
     bool ignore_nan(false);
@@ -50,7 +50,7 @@ mwSize accumulate_cut(double* s, double* e, double* npix,
         ignore_inf = true;
     }
     ignore_something = ignore_nan || ignore_inf;
-    ignote_all = ignore_nan && ignore_inf;
+    ignore_all = ignore_nan && ignore_inf;
     if (ignore_inf)
     {
         Inf = static_cast<T>(mxGetInf());
@@ -143,7 +143,7 @@ mwSize accumulate_cut(double* s, double* e, double* npix,
            pStor)                                                                            \
         firstprivate(data_size, distribution_size,                                           \
                      trans_elo, ebin_inv, Inf, PIXEL_data_width,                             \
-                     ignote_all, ignore_nan, ignore_inf, ignore_something, transform_energy, \
+                     ignore_all, ignore_nan, ignore_inf, ignore_something, transform_energy, \
                      nDimX, nDimY, nDimZ, nDimE,                                             \
                      s, e, npix,pixel_data)                                                  \
         reduction(+: nPixel_retained)
@@ -158,7 +158,7 @@ mwSize accumulate_cut(double* s, double* e, double* npix,
             ok[i] = false;
             if (ignore_something)
             {
-                if (ignote_all)
+                if (ignore_all)
                 {
                     if (pixel_data[j0 + 7] == Inf || isNaN(pixel_data[j0 + 7]) ||
                         pixel_data[j0 + 8] == Inf || isNaN(pixel_data[j0 + 8]))

_LowLevelCode/cpp/bin_pixels_c/BinningArg.cpp

@@ -264,7 +264,9 @@ void BinningArg::set_all_pix(mxArray const* const pField)
             auto data_type = mxGetClassID(a_cell_ptr);
             if (data_type != mxDOUBLE_CLASS) {
                 std::stringstream buf;
-                buf << "Binning the dataype N " << data_type << " have not been implemented yet";
+                buf << " Setting input pixel data to bin\n";
+                buf << " Binning for input datatype N: " << data_type << " have not been implemented yet\n";
+                buf << " You can bin datatype N " << mxDOUBLE_CLASS << " only";
                 mexErrMsgIdAndTxt("HORACE:bin_pixels_c:not_implemented",
                     buf.str().c_str());
             }
@@ -389,8 +391,9 @@ void BinningArg::return_pix_range(mxArray* pFieldName, mxArray* pFieldValue, int
         pix_range = this->pix_data_range_ptr;
     }
     mxSetCell(pFieldValue, fld_idx, pix_range);
+    this->pix_data_range_ptr = nullptr;
 };
-// return pixel obtained after binning and may be sorting.Sets up empty matrix if algorithm have not been using pixels
+// return pixel obtained after binning and may be sorting. Sets up empty matrix if algorithm have not been using pixels
 void BinningArg::return_pix_ok_data(mxArray* pFieldName, mxArray* pFieldValue, int fld_idx, const std::string& field_name)
 {
     mxSetCell(pFieldName, fld_idx, mxCreateString(field_name.c_str()));
@@ -401,7 +404,36 @@ void BinningArg::return_pix_ok_data(mxArray* pFieldName, mxArray* pFieldValue, i
         pix_ok = mxCreateDoubleMatrix(0, 0, mxREAL);
     }
     mxSetCell(pFieldValue, fld_idx, pix_ok);
+    this->pix_ok_ptr = nullptr;
 };
+// return array of pixels indices, which specify position of pixels within image cell
+void BinningArg::return_pix_img_idx(mxArray* pFieldName, mxArray* pFieldValue, int fld_idx, const std::string& field_name)
+{
+    mxSetCell(pFieldName, fld_idx, mxCreateString(field_name.c_str()));
+    mxArray* pix_img_idx(nullptr);
+    if (this->pix_img_idx_ptr) {
+        pix_img_idx = this->pix_img_idx_ptr;
+    } else {
+        pix_img_idx = mxCreateNumericMatrix(0, 0, mxINT64_CLASS, mxREAL);
+    }
+    mxSetCell(pFieldValue, fld_idx, pix_img_idx);
+    this->pix_img_idx_ptr = nullptr;
+};
+// return array of pixels indices, which specify position of pixels within image cell
+void BinningArg::return_is_pix_selected(mxArray* pFieldName, mxArray* pFieldValue, int fld_idx, const std::string& field_name)
+{
+    mxSetCell(pFieldName, fld_idx, mxCreateString(field_name.c_str()));
+    mxArray* is_pix_selected(nullptr);
+    if (this->is_pix_selected_ptr) {
+        is_pix_selected = this->is_pix_selected_ptr;
+    } else {
+        is_pix_selected = mxCreateLogicalMatrix(0, 0);
+    }
+    mxSetCell(pFieldValue, fld_idx, is_pix_selected);
+    this->is_pix_selected_ptr = nullptr;
+};
+
+
 //===================================================================================
 // calculate steps used in binning over non-unit directions and numbers of these dimensions
 void BinningArg::calc_step_sizes_pax_and_strides()
@@ -458,11 +490,24 @@ void BinningArg::register_output_methods()
     this->Mode5ParList = this->Mode4ParList;
     this->Mode5ParList["unique_runid"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_unique_runid(p1, p2, idx, name); };
 
+    this->Mode6ParList = this->Mode5ParList;
+    this->Mode6ParList["pix_img_idx"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_pix_img_idx(p1, p2, idx, name); };
+
+    this->Mode7ParList = this->Mode6ParList;
+    this->Mode7ParList["is_pix_selected"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_is_pix_selected(p1, p2, idx, name); };
+
+    this->Mode8ParList = this->Mode0ParList;
+    this->Mode8ParList["is_pix_selected"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_is_pix_selected(p1, p2, idx, name); };
+
+
     this->out_handlers[opModes::npix_only] = &Mode0ParList;
     this->out_handlers[opModes::sig_err] = &Mode0ParList;
     this->out_handlers[opModes::sigerr_cell] = &Mode0ParList;
     this->out_handlers[opModes::sort_pix] = &Mode4ParList;
     this->out_handlers[opModes::sort_and_uid] = &Mode5ParList;
+    this->out_handlers[opModes::nosort] = &Mode6ParList;
+    this->out_handlers[opModes::nosort_sel] = &Mode7ParList;
+    this->out_handlers[opModes::siger_selected] = &Mode8ParList;
 };
 /**  Parse input binning arguments and set new BinningArg from MATLAB input arguments
  *    structure.
@@ -694,6 +739,8 @@ void BinningArg::return_test_inputs(mxArray* plhs[], int nlhs)
     this->pix_ok_ptr = mxDuplicateArray(this->all_pix_ptr);
     this->OutParList["pix_ok_data"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_pix_ok_data(p1, p2, idx, name); };
     this->OutParList["unique_runid"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_unique_runid(p1, p2, idx, name); };
+    this->OutParList["pix_img_idx"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_pix_img_idx(p1, p2, idx, name); };
+    this->OutParList["is_pix_selected"] = [this](mxArray* p1, mxArray* p2, int idx, const std::string& name) { this->return_is_pix_selected(p1, p2, idx, name); };
 
     /* ********************************************************************************
      * retrieve binning parameters form BinningArg class and copy them into output array
@@ -813,28 +860,29 @@ void BinningArg::check_and_init_accumulators(mxArray* plhs[], mxArray const* prh
         plhs[out_arg::Signal] = this->signal_ptr;
         plhs[out_arg::Error] = this->error_ptr;
     }
-    if (this->binMode >= opModes::sort_pix) {
+    // pixels modes
+    if (this->binMode >= opModes::sort_pix && this->binMode<opModes::siger_selected) {
         if (this->n_data_points > this->pix_ok_bin_idx.size()) {
             this->pix_ok_bin_idx.resize(this->n_data_points);
         }
         // fill all positions of the pix_ok vector with certainly invalid value. Index can not be negative
         // this will indicate invalid elements
         std::fill(this->pix_ok_bin_idx.begin(), this->pix_ok_bin_idx.end(), -1);
-        if (this->npix_bin_start.size() != this->distr_size) {
-            this->npix_bin_start.resize(this->distr_size);
-            this->npix1.resize(this->distr_size);
-        }
-        std::fill(this->npix1.begin(), this->npix1.end(), 0); // nullify accumulators for npix1
-        // ranges calculated per each pixels block, i.e. calculations per call to bin_pixels_c
+        // allocate space for pixel data range, which is always calculated for
+        // any pixel mode
         this->pix_data_range_ptr = mxCreateDoubleMatrix(2, pix_flds::PIX_WIDTH, mxREAL);
         if (init_new_accumulators) {
             // just in case, clear unique run-id-s set (should be empty anyway at this stage)
             this->unique_runID.clear();
         }
-    } else {
-        // clear ranges left from previous call to binning functions
-        // as they are not needed for this call but may confuse output routines
-        this->pix_data_range_ptr = nullptr;
+        if (this->binMode < opModes::nosort) {
+            // accumulators needed for sorting pixels according to bins
+            if (this->npix_bin_start.size() != this->distr_size) {
+                this->npix_bin_start.resize(this->distr_size);
+                this->npix1.resize(this->distr_size);
+            }
+            std::fill(this->npix1.begin(), this->npix1.end(), 0); // nullify accumulators for npix1
+        }
     }
 }
 
@@ -905,6 +953,8 @@ BinningArg::BinningArg()
     , n_pix_retained(0)
     , pix_data_range_ptr(nullptr)
     , pix_ok_ptr(nullptr)
+    , pix_img_idx_ptr(nullptr)
+    , is_pix_selected_ptr(nullptr)
 {
     /* initialize input Matlab parameters map with methods which  associate
      * Matlab field names with the methods, which set appropriate property value  */

_LowLevelCode/cpp/bin_pixels_c/BinningArg.h

@@ -54,7 +54,7 @@ enum opModes {
     sort_and_uid = 5, // in additional to binning and sorting, return unique pixels id
     nosort = 6, // do binning but do not sort pixels but return array which defines pixels position
     //              within the image grid
-    nosort_sel = 7, // like 6, but return ?logical? array which specifies what pixels have been selected
+    nosort_sel = 7, // like 6, but return logical array which specifies what pixels have been selected
     //                   and what were rejected by binning operations
     siger_selected = 8, // the same as sig_err but return logical array of selected piels instead of pix_ok array
     test_inputs = 9, // do not do calculations but just return parsed inputs for
@@ -111,6 +111,9 @@ class BinningArg {
     mxArray* pix_data_range_ptr;
     mxArray* pix_ok_ptr; // pointer to array of all pixels retained after binning
     std::unordered_set<uint32_t> unique_runID; // set containing unique run_id-s of the
+    mxArray* pix_img_idx_ptr; // pointer to array of pixel indices within the image cell 
+    mxArray* is_pix_selected_ptr; // pointer to logical array containing true where pixel 
+    //                               was selected for binning and false otherwise
     // processed pixels
     //********************************************************************************
     // helper values
@@ -119,7 +122,7 @@ class BinningArg {
     std::vector<size_t> stride; // vector, which describes binning steps reflecting multidimensional array strides
     std::vector<size_t> bin_cell_idx_range; // vector containing allowed maximal indixes of the binning (with nbins_all_dims>1) cells in binning directions
     // auxiliary array containing pixel indices over bins
-    std::vector<long> pix_ok_bin_idx;
+    std::vector<mxInt64> pix_ok_bin_idx;
     // auxiliary array defining ranges of the bins to sort pixels over
     std::vector<size_t> npix_bin_start;
     std::vector<size_t> npix1; // pixel distribution over bins calculated in single call to bin_pixels routine;
@@ -149,13 +152,12 @@ class BinningArg {
     void set_nbins_all_dims(mxArray const* const pField); //
     void set_unique_runid(mxArray const* const pField); // holder for the information about unique run_id-s present in the data. Set procedure is non-standard
     void set_force_double(mxArray const* const pField); // boolean parameters which would request output transformed pixels always been double regardless of input pixels
-    void set_return_selected(mxArray const* const pField);
     void set_test_input_mode(mxArray const* const pField); // intialize testing mode (or not)
     void set_all_pix(mxArray const* const pField); //  pointer to all pixels to sort or use as binning arguments
     void set_alignment_matrix(mxArray const* const pField); // matrix which have to be applied to raw pixels to bring them into Crystal Cartesian coordinate system
     void set_check_pix_selection(mxArray const* const pField); // if true, check if detector_id are negative which may suggest that pixels have been alreary used in previous binning operation
 
-    // register with parameters map all methods which
+    // register with parameters map all methods which return variable results to MATLAB
     void register_output_methods();
     // setters for binning results returned to MATLAB in output structure
     void return_npix_retained(mxArray* p1, mxArray* p2, int idx, const std::string& name);
@@ -165,7 +167,11 @@ class BinningArg {
     void return_pix_ok_data(mxArray* p1, mxArray* p2, int idx, const std::string& name);
     // setter to return unique run_id,.calculated in the call
     void return_unique_runid(mxArray* p1, mxArray* p2, int idx, const std::string& name);
-
+    // setter to return pixels indices within the image cell
+    void return_pix_img_idx(mxArray* p1, mxArray* p2, int fld_idx, const std::string& name);
+    // setter to return array containing true for selected pixels and false if they have not been selected
+    void return_is_pix_selected(mxArray* p1, mxArray* p2, int fld_idx, const std::string& name);
+    //
 public:
     BinningArg(); // construction
     // process binning arguments input values for new binning arguments cycle
@@ -198,6 +204,9 @@ class BinningArg {
     OutHandlerMap Mode0ParList;
     OutHandlerMap Mode4ParList;
     OutHandlerMap Mode5ParList;
+    OutHandlerMap Mode6ParList;
+    OutHandlerMap Mode7ParList;
+    OutHandlerMap Mode8ParList;
 
     std::unordered_map<opModes, OutHandlerMap*> out_handlers;
 };