Implement unweighted average for 2D integration

src/pyFAI/azimuthalIntegrator.py

@@ -2320,8 +2320,7 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                                                      deviceid=method.target[1],
                                                                      checksum=cython_integr.lut_checksum,
                                                                      unit=unit, empty=empty,
-                                                                     mask_checksum=mask_crc
-                                                                     )
+                                                                     mask_checksum=mask_crc)
 
                         elif (method.impl_lower == "python"):
                             with ocl_py_engine.lock:
@@ -2348,7 +2347,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                                        polarization=polarization,
                                                        polarization_checksum=polarization_crc,
                                                        safe=safe,
-                                                       normalization_factor=normalization_factor)
+                                                       normalization_factor=normalization_factor,
+                                                       weighted_average=method.weighted_average,)
             if intpl is None:  # fallback if OpenCL failed or default cython
                 # The integrator has already been initialized previously
                 intpl = cython_integr.integrate_ng(data,
@@ -2360,7 +2360,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                                    flat=flat,
                                                    solidangle=solidangle,
                                                    polarization=polarization,
-                                                   normalization_factor=normalization_factor)
+                                                   normalization_factor=normalization_factor,
+                                                   weighted_average=method.weighted_average,)
 
         elif method.algo_lower == "histogram":
             if method.split_lower in ("pseudo", "full"):
@@ -2384,7 +2385,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                    empty=empty,
                                    variance=variance,
                                    error_model=error_model,
-                                   allow_pos0_neg=not radial_unit.positive)
+                                   allow_pos0_neg=not radial_unit.positive,
+                                   weighted_average=method.weighted_average,)
 
             elif method.split_lower == "bbox":
                 logger.debug("integrate2d uses BBox implementation")
@@ -2414,7 +2416,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                      variance=variance,
                                      error_model=error_model,
                                      allow_pos0_neg=not radial_unit.positive,
-                                     clip_pos1=bool(azimuth_unit.period))
+                                     clip_pos1=bool(azimuth_unit.period),
+                                     weighted_average=method.weighted_average,)
             elif method.split_lower == "no":
                 if method.impl_lower == "opencl":
                     logger.debug("integrate2d uses OpenCL histogram implementation")
@@ -2486,7 +2489,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                                  normalization_factor=normalization_factor,
                                                  radial_range=radial_range,
                                                  azimuthal_range=azimuth_range,
-                                                 error_model=error_model)
+                                                 error_model=error_model,
+                                                 weighted_average=method.weighted_average,)
 ####################
                 else:  # if method.impl_lower in ["python", "cython"]:
                     logger.debug("integrate2d uses [CP]ython histogram implementation")
@@ -2499,10 +2503,6 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                                 radial_range=radial_range,
                                                 azimuth_range=azimuth_range,
                                                 mode="normal").ravel()
-                        # if radial_range is None:
-                        #     radial_range = [radial.min(), radial.max()]
-                        # if azimuth_range is None:
-                        #     azimuth_range = [azim.min(), azim.max()]
                     histogrammer = method.class_funct_ng.function
                     intpl = histogrammer(radial=radial,
                                          azimuthal=azim,
@@ -2524,7 +2524,8 @@ def integrate2d_ng(self, data, npt_rad, npt_azim=360,
                                          radial_range=radial_range,
                                          azimuth_range=azimuth_range,
                                          allow_radial_neg=not radial_unit.positive,
-                                         clip_pos1=bool(azimuth_unit.period))
+                                         clip_pos1=bool(azimuth_unit.period),
+                                         weighted_average=method.weighted_average,)
 
         I = intpl.intensity
         bins_azim = intpl.azimuthal

src/pyFAI/engines/CSC_engine.py

@@ -1,5 +1,5 @@
 #
-#    Copyright (C) 2017-2023 European Synchrotron Radiation Facility, Grenoble, France
+#    Copyright (C) 2017-2024 European Synchrotron Radiation Facility, Grenoble, France
 #
 #  Permission is hereby granted, free of charge, to any person obtaining a copy
 #  of this software and associated documentation files (the "Software"), to deal
@@ -26,7 +26,7 @@
 __contact__ = "Jerome.Kieffer@ESRF.eu"
 __license__ = "MIT"
 __copyright__ = "European Synchrotron Radiation Facility, Grenoble, France"
-__date__ = "04/10/2023"
+__date__ = "26/04/2024"
 __status__ = "development"
 
 import logging
@@ -87,7 +87,7 @@ def set_matrix(self, data, indices, indptr):
         self.indptr = indptr
         self.lut_size = len(indices)
         if self.size > len(indptr) - 1:
-            new_indptr = numpy.empty(self.size+1, indptr.dtype)
+            new_indptr = numpy.empty(self.size + 1, indptr.dtype)
             new_indptr[:] = indptr[-1]
             new_indptr[:len(indptr)] = indptr
             indptr = new_indptr
@@ -147,7 +147,7 @@ def integrate(self,
                        variance=variance,
                        dtype=numpy.float32,
                        error_model=error_model,
-                       apply_normalization = not weighted_average,
+                       apply_normalization=not weighted_average,
                        out=self.preprocessed)
         prep.shape = numpy.prod(shape), 4
         flat_sig, flat_var, flat_nrm, flat_cnt = prep.T  # should create views!
@@ -450,6 +450,7 @@ def set_matrix(self, data, indices, indptr):
         :param indptr: the index of the start of line"""
 
         CSCIntegrator.set_matrix(self, data, indices, indptr)
+
         assert self.size == len(indptr) - 1
 
     def integrate(self,
@@ -464,6 +465,7 @@ def integrate(self,
                   polarization=None,
                   absorption=None,
                   normalization_factor=1.0,
+                  weighted_average=True,
                   **kwargs):
         """Actually perform the 2D integration
 
@@ -478,14 +480,15 @@ def integrate(self,
         :param polarization: :solidangle normalization array
         :param absorption: :absorption normalization array
         :param normalization_factor: scale all normalization with this scalar
+        :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average
         :return: Integrate2dtpl namedtuple: "radial azimuthal intensity error signal variance normalization count"
 
         """
         error_model = ErrorModel.parse(error_model)
         do_variance = variance is not None or  error_model.do_variance
         trans = CSCIntegrator.integrate(self, signal, variance, error_model, dummy, delta_dummy,
                                         dark, flat, solidangle, polarization,
-                                        absorption, normalization_factor)
+                                        absorption, normalization_factor, weighted_average=weighted_average)
         trans.shape = self.bins + (-1,)
 
         signal = trans[..., 0]

src/pyFAI/engines/CSR_engine.py

@@ -26,7 +26,7 @@
 __contact__ = "Jerome.Kieffer@ESRF.eu"
 __license__ = "MIT"
 __copyright__ = "European Synchrotron Radiation Facility, Grenoble, France"
-__date__ = "24/04/2024"
+__date__ = "25/04/2024"
 __status__ = "development"
 
 import logging
@@ -138,7 +138,7 @@ def integrate(self,
                        variance=variance,
                        dtype=numpy.float32,
                        error_model=error_model,
-                       apply_normalization = not weighted_average,
+                       apply_normalization=not weighted_average,
                        out=self.preprocessed)
         prep.shape = numpy.prod(shape), 4
         flat_sig, flat_var, flat_nrm, flat_cnt = prep.T  # should create views!
@@ -431,8 +431,6 @@ def __init__(self,
             self.checksum = checksum
         CSRIntegrator.__init__(self, image_size, lut, empty)
 
-
-
     def set_matrix(self, data, indices, indptr):
         """Actually set the CSR sparse matrix content
 
@@ -456,6 +454,7 @@ def integrate(self,
                   polarization=None,
                   absorption=None,
                   normalization_factor=1.0,
+                  weighted_average=True,
                   **kwargs):
         """Actually perform the 2D integration
 
@@ -470,14 +469,15 @@ def integrate(self,
         :param polarization: :solidangle normalization array
         :param absorption: :absorption normalization array
         :param normalization_factor: scale all normalization with this scalar
+        :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average
         :return: Integrate2dtpl namedtuple: "radial azimuthal intensity error signal variance normalization count"
 
         """
         error_model = ErrorModel.parse(error_model)
         do_variance = variance is not None or  error_model.do_variance
         trans = CSRIntegrator.integrate(self, signal, variance, error_model, dummy, delta_dummy,
                                         dark, flat, solidangle, polarization,
-                                        absorption, normalization_factor)
+                                        absorption, normalization_factor, weighted_average=weighted_average)
         trans.shape = self.bins + (-1,)
 
         signal = trans[..., 0]

src/pyFAI/engines/histogram_engine.py

@@ -1,5 +1,5 @@
 #
-#    Copyright (C) 2019-2022 European Synchrotron Radiation Facility, Grenoble, France
+#    Copyright (C) 2019-2024 European Synchrotron Radiation Facility, Grenoble, France
 #
 #  Permission is hereby granted, free of charge, to any person obtaining a copy
 #  of this software and associated documentation files (the "Software"), to deal
@@ -22,11 +22,11 @@
 """simple histogram rebinning engine implemented in pure python (with the help of numpy !)
 """
 
-__author__ = "Jerome Kieffer"
+__author__ = "Jérôme Kieffer"
 __contact__ = "Jerome.Kieffer@ESRF.eu"
 __license__ = "MIT"
 __copyright__ = "European Synchrotron Radiation Facility, Grenoble, France"
-__date__ = "24/04/2024"
+__date__ = "26/04/2024"
 __status__ = "development"
 
 import logging
@@ -112,7 +112,7 @@ def histogram1d_engine(radial, npt,
                    dark_variance=dark_variance,
                    error_model=error_model,
                    empty=0,
-                   apply_normalization = not weighted_average,
+                   apply_normalization=not weighted_average,
                    )
     radial = radial.ravel()
     prep.shape = -1, 4
@@ -173,6 +173,7 @@ def histogram2d_engine(radial, azimuthal, bins,
                        variance=None,
                        dark_variance=None,
                        error_model=ErrorModel.NO,
+                       weighted_average=True,
                        radial_range=None,
                        azimuth_range=None,
                        allow_radial_neg=False,
@@ -198,6 +199,7 @@ def histogram2d_engine(radial, azimuthal, bins,
     :param variance: provide an estimation of the variance
     :param dark_variance: provide an estimation of the variance of the dark_current,
     :param error_model: set to "poisson" for assuming the detector is poissonian and variance = raw + dark
+    :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average
     :param radial_range: enforce boundaries in radial dimention, 2tuple with lower and upper bound
     :param azimuth_range: enforce boundaries in azimuthal dimention, 2tuple with lower and upper bound
     :param allow_radial_neg: clip negative radial position (can a dimention be negative ?)
@@ -233,7 +235,8 @@ def histogram2d_engine(radial, azimuthal, bins,
                    variance=variance,
                    dark_variance=dark_variance,
                    error_model=error_model,
-                   empty=0
+                   empty=0,
+                   apply_normalization=not weighted_average,
                    )
     radial = radial.ravel()
     azimuthal = azimuthal.ravel()

src/pyFAI/engines/preproc.py

@@ -31,7 +31,7 @@
 __contact__ = "Jerome.Kieffer@ESRF.eu"
 __license__ = "MIT"
 __copyright__ = "European Synchrotron Radiation Facility, Grenoble, France"
-__date__ = "23/04/2024"
+__date__ = "26/04/2024"
 __status__ = "development"
 
 import warnings
@@ -214,7 +214,7 @@ def preproc(raw,
                 numpy.divide(signal, normalization, out=signal)
                 if variance is not None:
                     numpy.divide(variance, normalization*normalization, out=variance)
-                normalization[:] = 1.0
+                normalization[...] = 1.0
 
             signal[mask] = 0.0
             normalization[mask] = 0.0

src/pyFAI/ext/LUT_common.pxi

@@ -29,7 +29,7 @@
 
 __author__ = "Jerome Kieffer"
 __contact__ = "Jerome.kieffer@esrf.fr"
-__date__ = "23/01/2024"
+__date__ = "26/04/2024"
 __status__ = "stable"
 __license__ = "MIT"
 
@@ -89,10 +89,9 @@ cdef class LutIntegrator(object):
         """Getter a copy of the LUT as an actual numpy array"""
         cdef double[:, ::1] tmp_ary = numpy.empty((self.output_size, self.lut_size), dtype=numpy.float64)
         memcpy(&tmp_ary[0, 0], &self._lut[0, 0], self._lut.nbytes)
-
-        return numpy.core.records.array(numpy.asarray(tmp_ary).view(dtype=lut_d),
-                                        shape=(self.output_size, self.lut_size), dtype=lut_d,
-                                        copy=True)
+        return numpy.recarray(buf=numpy.array(tmp_ary, copy=True).view(dtype=lut_d),
+                              shape=(self.output_size, self.lut_size),
+                              dtype=lut_d)
 
 
     def integrate_legacy(self, weights,

src/pyFAI/ext/histogram.pyx

@@ -38,7 +38,7 @@ Can be replaced by ``silx.math.histogramnd``.
 """
 
 __author__ = "Jérôme Kieffer"
-__date__ = "24/04/2024"
+__date__ = "25/04/2024"
 __license__ = "MIT"
 __copyright__ = "2011-2022, ESRF"
 __contact__ = "jerome.kieffer@esrf.fr"
@@ -475,7 +475,7 @@ def histogram1d_engine(radial, int npt,
     :param variance: provide an estimation of the variance
     :param dark_variance: provide an estimation of the variance of the dark_current,
     :param error_model: One of the several ErrorModel, only variance and Poisson are implemented.
-    :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average. WIP
+    :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average.
 
 
     NaN are always considered as invalid values
@@ -514,7 +514,7 @@ def histogram1d_engine(radial, int npt,
                    variance=variance,
                    dark_variance=dark_variance,
                    error_model=error_model,
-                   apply_normalization = not weighted_average,
+                   apply_normalization=not weighted_average,
                    ).reshape(-1, 4)
     res, position = histogram_preproc(radial.ravel(),
                                       prep,
@@ -580,6 +580,7 @@ def histogram2d_engine(radial, azimuthal,
                        variance=None,
                        dark_variance=None,
                        int error_model=ErrorModel.NO,
+                       bint weighted_average=True,
                        radial_range=None,
                        azimuth_range=None,
                        bint allow_radial_neg=False,
@@ -605,6 +606,7 @@ def histogram2d_engine(radial, azimuthal,
     :param variance: provide an estimation of the variance
     :param dark_variance: provide an estimation of the variance of the dark_current,
     :param error_model: set to "poisson" for assuming the detector is poissonian and variance = raw + dark
+    :param bool weighted_average: set to False to use an unweighted mean (similar to legacy) instead of the weighted average.
     :param radial_range: enforce boundaries in radial dimention, 2tuple with lower and upper bound
     :param azimuth_range: enforce boundaries in azimuthal dimention, 2tuple with lower and upper bound
     :param allow_radial_neg: clip negative radial position (can a dimention be negative ?)
@@ -740,7 +742,8 @@ def histogram2d_engine(radial, azimuthal,
                                              check_dummy=check_dummy,
                                              normalization_factor=normalization_factor,
                                              dark_variance=cdark_variance[idx] if do_dark_variance else 0.0,
-                                             error_model=error_model)
+                                             error_model=error_model,
+                                             apply_normalization = not weighted_average,)
 
             if not is_valid:
                 continue

src/pyFAI/ext/preproc.pyx

@@ -38,7 +38,7 @@ flat-field normalization... taking care of masked values and normalization.
 
 __author__ = "Jerome Kieffer"
 __license__ = "MIT"
-__date__ = "23/04/2024"
+__date__ = "26/04/2024"
 __copyright__ = "2011-2022, ESRF"
 __contact__ = "jerome.kieffer@esrf.fr"