Merge pull request #19 from picca/master

attenuation of the images

1 files changed, 43 insertions, 8 deletions

diff --git a/binoculars/backends/sixs.py b/binoculars/backends/sixs.py
index 218cba2..4b213c1 100644
--- a/binoculars/backends/sixs.py
+++ b/binoculars/backends/sixs.py
@@ -147,6 +147,9 @@ class QparQperProjection(QxQyQzProjection):
 # Common methodes #
 ###################
 
+WRONG_ATTENUATION = -100
+
+
 def get_nxclass(hfile, nxclass, path="/"):
     """
     :param hfile: the hdf5 file.
@@ -372,6 +375,11 @@ class SIXS(backend.InputBase):
             except ValueError:
                 self.config.detrot = None
 
+        # attenuation_coefficient (Optional)
+        attenuation_coefficient = config.pop('attenuation_coefficient', None)
+        self.config.attenuation_coefficient =\
+            float(attenuation_coefficient) if attenuation_coefficient else None
+
     def get_destination_options(self, command):
         if not command:
             return False
@@ -399,6 +407,7 @@ class FlyScanUHV(SIXS):
         "omega": "UHV_OMEGA",
         "delta": "UHV_DELTA",
         "gamma": "UHV_GAMMA",
+        "attenuation": "attenuation",
     }
 
     def get_pointcount(self, scanno):
@@ -406,15 +415,24 @@ class FlyScanUHV(SIXS):
         with tables.open_file(self.get_filename(scanno), 'r') as scan:
             return get_nxclass(scan, "NXdata").xpad_image.shape[0]
 
+    def get_attenuation(self, index, h5_nodes, offset):
+        attenuation = None
+        if self.config.attenuation_coefficient is not None:
+            try:
+                attenuation = h5_nodes['attenuation'][index + offset]
+            except IndexError:
+                attenuation = WRONG_ATTENUATION
+        return attenuation
+
     def get_values(self, index, h5_nodes):
         image = h5_nodes['image'][index]
         mu = h5_nodes['mu'][index]
         omega = h5_nodes['omega'][index]
         delta = h5_nodes['delta'][index]
         gamma = h5_nodes['gamma'][index]
+        attenuation = self.get_attenuation(index, h5_nodes, 2)
 
-        return (image,
-                (mu, omega, delta, gamma))
+        return (image, attenuation, (mu, omega, delta, gamma))
 
     def process_image(self, index, dataframe, pixels):
         util.status(str(index))
@@ -428,10 +446,25 @@ class FlyScanUHV(SIXS):
         # extract the data from the h5 nodes
 
         h5_nodes = dataframe.h5_nodes
-        intensity, values = self.get_values(index, h5_nodes)
-
-        weights = numpy.ones_like(intensity)
-        weights *= ~mask
+        intensity, attenuation, values = self.get_values(index, h5_nodes)
+
+        # BEWARE in order to avoid precision problem we convert the
+        # uint16 -> float32. (the size of the mantis is on 23 bits)
+        # enought to contain the uint16. If one day we use uint32, it
+        # should be necessary to convert into float64.
+        intensity = intensity.astype('float32')
+
+        weights = None
+        if self.config.attenuation_coefficient is not None:
+            if attenuation != WRONG_ATTENUATION:
+                intensity *= self.config.attenuation_coefficient ** attenuation
+                weights = numpy.ones_like(intensity)
+                weights *= ~mask
+            else:
+                weights = numpy.zeros_like(intensity)
+        else:
+            weights = numpy.ones_like(intensity)
+            weights *= ~mask
 
         k = 2 * math.pi / dataframe.source.wavelength
 
@@ -477,6 +510,7 @@ class FlyScanUHV2(FlyScanUHV):
         "omega": "omega",
         "delta": "delta",
         "gamma": "gamma",
+        "attenuation": "attenuation",
     }
 
 
@@ -487,6 +521,7 @@ class SBSMedH(FlyScanUHV):
         "mu": "data_18",
         "gamma": "data_20",
         "delta": "data_19",
+        "attenuation": "data_xx",
         }
 
     def get_pointcount(self, scanno):
@@ -500,9 +535,9 @@ class SBSMedH(FlyScanUHV):
         mu = h5_nodes['mu'][index]
         gamma = h5_nodes['gamma'][index]
         delta = h5_nodes['delta'][index]
+        attenuation = self.get_attenuation(index, h5_nodes, 2)
 
-        return (image,
-                (pitch, mu, gamma, delta))
+        return (image, attenuation, (pitch, mu, gamma, delta))
 
 
 def load_matrix(filename):