Re-write FormatXTCJungfrau to use homogenous coordinates (LCLS Jungfr…

…au only)

Old code used psana's get_pixel_coords to read an array of pixel coordinates in lab space, then np.mean to find the centers of the objects.  This version uses the vectors in the psana calib file directly.

This change also fixes a bug where the parent frame was applied twice at the asic level.  Not noticable in our 1M experiments, but the more complicated 4M showed the bug. This new version uses PSCalib.SegGeometryStore to look up the ASIC layout and gaps in the ASIC frame of reference instead of using the GeoAccess layout, which is in the module frame of reference.

Also rename quad to module as the 1M and 4M at LCLS do not have quadrants in their calib files.

format/FormatXTCJungfrau.py

@@ -68,17 +68,17 @@ def get_raw_data(self, index):
         data = det.calib(evt)
         data = data.astype(np.float64)
         self._raw_data = []
-        for quad_count, quad in enumerate(d.hierarchy()):
-            for asic_count, asic in enumerate(quad):
+        for module_count, module in enumerate(d.hierarchy()):
+            for asic_count, asic in enumerate(module):
                 fdim, sdim = asic.get_image_size()
-                sensor_id = asic_count // 4  # There are 2X4 asics per quadrant
+                sensor_id = asic_count // 4  # There are 2X4 asics per module
                 asic_in_sensor_id = asic_count % 4  # this number will be 0,1,2 or 3
-                asic_data = data[quad_count][
+                asic_data = data[module_count][
                     sensor_id * sdim : (sensor_id + 1) * sdim,
                     asic_in_sensor_id * fdim : (asic_in_sensor_id + 1) * fdim,
-                ]  # 8 sensors per quad
+                ]  # 8 sensors per module
                 self._raw_data.append(flex.double(np.array(asic_data)))
-        assert len(d) == len(self._raw_data)
+        assert len(d) == len(self._raw_data), (len(d), len(self._raw_data))
         return tuple(self._raw_data)
 
     def get_num_images(self):
@@ -107,6 +107,9 @@ def get_scan(self, index=None):
         return None
 
     def _detector(self, index=None):
+        from xfel.cftbx.detector.cspad_cbf_tbx import basis_from_geo
+        from PSCalib.SegGeometryStore import sgs
+
         run = self.get_run_from_index(index)
         if run.run() in self._cached_detector:
             return self._cached_detector[run.run()]
@@ -124,69 +127,54 @@ def _detector(self, index=None):
         pg0 = d.hierarchy()
         # first deal with D0
         det_num = 0
-        D0 = geom.get_top_geo().get_list_of_children()[0]
-        xx, yy, zz = D0.get_pixel_coords()
-        xx = xx / 1000.0  # to mm
-        yy = yy / 1000.0  # to mm
-        zz = zz / 1000.0  # to mm
-        oriD0 = col((np.mean(xx), np.mean(yy), -np.mean(zz)))
-        fp = col((xx[0][0][1], yy[0][0][1], zz[0][0][1]))
-        sp = col((xx[0][1][0], yy[0][1][0], zz[0][1][0]))
-        op = col((xx[0][0][0], yy[0][0][0], zz[0][0][0]))
-        origin = oriD0
-        fast = (fp - op).normalize()
-        slow = (sp - op).normalize()
+        root = geom.get_top_geo()
+        root_basis = basis_from_geo(root)
+        while len(root.get_list_of_children()) == 1:
+            sub = root.get_list_of_children()[0]
+            sub_basis = basis_from_geo(sub)
+            root = sub
+            root_basis = root_basis * sub_basis
+        t = root_basis.translation
+        root_basis.translation = col((t[0], t[1], -t[2]))
+
+        origin = col((root_basis * col((0, 0, 0, 1)))[0:3])
+        fast = col((root_basis * col((1, 0, 0, 1)))[0:3]) - origin
+        slow = col((root_basis * col((0, 1, 0, 1)))[0:3]) - origin
         pg0.set_local_frame(fast.elems, slow.elems, origin.elems)
         pg0.set_name("D%d" % (det_num))
 
-        # Now deal with Qx
-        for quad_num in range(2):
+        # Now deal with modules
+        for module_num in range(len(root.get_list_of_children())):
             pg1 = pg0.add_group()
-            Qx = D0.get_list_of_children()[quad_num]
-            xx, yy, zz = Qx.get_pixel_coords()
-            xx = xx / 1000.0  # to mm
-            yy = yy / 1000.0  # to mm
-            zz = zz / 1000.0  # to mm
-            oriQx = col((np.mean(xx), np.mean(yy), np.mean(zz)))
-            fp = col((xx[0][1], yy[0][1], zz[0][1]))
-            sp = col((xx[1][0], yy[1][0], zz[1][0]))
-            op = col((xx[0][0], yy[0][0], zz[0][0]))
-            origin = oriQx
-            fast = (fp - op).normalize()
-            slow = (sp - op).normalize()
+            module = root.get_list_of_children()[module_num]
+            module_basis = basis_from_geo(module)
+            origin = col((module_basis * col((0, 0, 0, 1)))[0:3])
+            fast = col((module_basis * col((1, 0, 0, 1)))[0:3]) - origin
+            slow = col((module_basis * col((0, 1, 0, 1)))[0:3]) - origin
             pg1.set_local_frame(fast.elems, slow.elems, origin.elems)
-            pg1.set_name("D%dQ%d" % (det_num, quad_num))
+            pg1.set_name("D%dM%d" % (det_num, module_num))
+
+            # Read the known layout of the Jungfrau 2x4 module
+            sg = sgs.Create(segname=module.oname)
+            xx, yy = sg.get_seg_xy_maps_um()
+            xx = xx / 1000
+            yy = yy / 1000
 
-            # Now deal with Az
+            # Now deal with ASICs
             for asic_num in range(8):
-                val = "ARRAY_D0Q%dA%d" % (quad_num, asic_num)
+                val = "ARRAY_D0M%dA%d" % (module_num, asic_num)
                 p = pg1.add_panel()
                 dim_slow = xx.shape[0]
                 dim_fast = xx.shape[1]
-                sensor_id = asic_num // 4  # There are 2X4 asics per quadrant
+                sensor_id = asic_num // 4  # There are 2X4 asics per module
                 asic_in_sensor_id = asic_num % 4  # this number will be 0,1,2 or 3
                 id_slow = sensor_id * (dim_slow // 2)
                 id_fast = asic_in_sensor_id * (dim_fast // 4)
-                oriAy = col(
-                    (xx[id_slow][id_fast], yy[id_slow][id_fast], zz[id_slow][id_fast])
-                )
-                fp = col(
-                    (
-                        xx[id_slow][id_fast + 1],
-                        yy[id_slow][id_fast + 1],
-                        zz[id_slow][id_fast + 1],
-                    )
-                )
-                sp = col(
-                    (
-                        xx[id_slow + 1][id_fast],
-                        yy[id_slow + 1][id_fast],
-                        zz[id_slow + 1][id_fast],
-                    )
-                )
-                origin = oriAy - oriQx
-                fast = (fp - oriAy).normalize()
-                slow = (sp - oriAy).normalize()
+                origin = col((xx[id_slow][id_fast], yy[id_slow][id_fast], 0,))
+                fp = col((xx[id_slow][id_fast + 1], yy[id_slow][id_fast + 1], 0))
+                sp = col((xx[id_slow + 1][id_fast], yy[id_slow + 1][id_fast], 0))
+                fast = (fp - origin).normalize()
+                slow = (sp - origin).normalize()
                 p.set_local_frame(fast.elems, slow.elems, origin.elems)
                 p.set_pixel_size((pixel_size, pixel_size))
                 p.set_image_size((dim_fast // 4, dim_slow // 2))