silx-kit / pyFAI

Fast Azimuthal Integration in Python
Other
104 stars 94 forks source link

tiled detectors: 3D data shape #1807

Open vespos opened 1 year ago

vespos commented 1 year ago

Hello,

It has become more common for detector data to be saved as 3D arrays, with one of the dimension being the tile index. Assuming that the (x,y,z) position of each pixel is known, is there a way to integrate 3D arrays (or flattened array), without having to make the corresponding 2D image? This would be particularly relevant to speed up cases where many images are integrated, as with FEL datasets for example. This would roughly gain a factor 2 in speed (as checked with homemade azav code) by going from a two-histogram steps to a single-histogram one.

I have tried to find a way to do that cleanly with pyFAI, but can't circumvent the need to go through the 2D image, without having to rewrite core parts of pyFAI (which I have not done).

Thanks in advance for the attention given to that question.

-- Vincent Esposito (SLAC, LCLS)

kif commented 1 year ago

Hi Vincent,

You are right, pyFAI was designed to take images as 2D arrays, and I don't know how much work it would require to change this. At many places in the upper layer of pyFAI, a 3D array is likely to be considered as a stack of images. At lower levels, this is all 1D buffer of memory. Some facilities already asked for the support of 1D detector (synchrotron Soleil with mythen detectors), and we worked around in reshaping 1d array in 2d. Easy.

I already played with images coming from the EuXFEL which had 16 tiles stacked in 3D, and the "same" trick did the job: concatenate the two first dimensions. So a detector of shape (8, 512, 1024) appears as a (4096, 1024). The major drawback is that the image, once displayed as a 2D array does not look nice. As long as the number of tiles remains reasonable (<8), it is possible to calibrate the detector position with existing tools (pyFAI-calib2). There is a reconstructed view of the detector and of the experimental setup, but it does not allow to perform the calibration on this representation (only the 2D representation of the array).

In pyFAI, detectors can be represented as a regular array of pixels (easy) or as an assembly of pixel, each pixel having [3, 4, 6] corners which position in 3D-space is known (for now the triangle and hexagonal pixel support is still partial). http://www.silx.org/doc/pyFAI/latest/usage/tutorial/Detector/Distortion/Distortion.html#WOS-XPad-detector Each pixel coordinate is known in 3D by the detector-object which can be serialized in HDF5. Converting the tile position into this array of pixel position is not very difficult ... and can be re-use as long as the detector's internal geometry does not change (likely to be an issue on XFEL). There are several tutorials on the topic which explain how to save the geometry of the detector and restore it in: http://www.silx.org/doc/pyFAI/latest/usage/tutorial/Detector/detector.html

Cheers,

Jerome

vespos commented 1 year ago

Thank you for your detailed response. I will explore the option of reshaping our 3D datasets into 2D. With accurate pixel corners' position, this should work fine.

-- Vincent Esposito

vespos commented 1 year ago

Hi Jerome,

I have a question that I could not resolve looking at the documentation: Is there a preferred order for the corner index in the pyFAI.detectors.Detector.set_pixel_corners (http://www.silx.org/doc/pyFAI/latest/api/detectors/index.html#pyFAI.detectors.Detector.set_pixel_corners) method? Or is it proofed against any order, even inconsistent across a same dataset?

Cheers, Vincent

kif commented 1 year ago

Hi Vincent,

If one assumes the images arrive in 2D arrays (and pyFAI does), those arrays are indexed in slow-speed dimension (aka dim0 or y) and fast-speed dimension (dim1 or x). The speed is relative to the stride index into the memory buffer (see https://github.com/numpy/numpy/blob/main/doc/source/reference/arrays.ndarray.rst#internal-memory-layout-of-an-ndarray) the lower the stride value, the faster. The C-memory layout convention is assumed everywhere in pyFAI and most (all) upper level function are protected against alternative memory-layout (it just slows down things by making a memory copy in the appropriate layout).

The order of this 4D array describing the position of every pixel corner is then pretty strict. :

You can have a look at the function going from a contiguous, flat detector to this 4D array: https://github.com/silx-kit/pyFAI/blob/main/pyFAI/detectors/_common.py#L649

vespos commented 1 year ago

Thanks! I guess my question was specifically if I should follow this order for the corner definition:

                    self._pixel_corners[:,:, 0, 1] = p1[:-1,:-1]
                    self._pixel_corners[:,:, 0, 2] = p2[:-1,:-1]
                    self._pixel_corners[:,:, 1, 1] = p1[1:,:-1]
                    self._pixel_corners[:,:, 1, 2] = p2[1:,:-1]
                    self._pixel_corners[:,:, 2, 1] = p1[1:, 1:]
                    self._pixel_corners[:,:, 2, 2] = p2[1:, 1:]
                    self._pixel_corners[:,:, 3, 1] = p1[:-1, 1:]
                    self._pixel_corners[:,:, 3, 2] = p2[:-1, 1:]

i.e. the corner index 0 to 3 have to be in the order: top left, bottom left, bottom right, top right; or if I can define them in any order I want?

Asking because our tiles are in random orientation, and not perfectly square. It would thus be much simpler if I don't have to care about the order of the corner definition

kif commented 1 year ago

Your tiles are rotated, probably by steps on 90°, and they are probably not inverted. If you consider the splitting scheme "no" and "bbox", the order has no importance at all (this is why hexagonal pixel are supported). It only matters for the "full" pixel splitting scheme where polygons are cut in pieces.

This is one the most tedious part of pyFAI. The corners of index [0123] should build a convex polygon, and never a crossing one. Since pixel splitting is performed by calculating area, if your turn in the wrong direction, area will change sign. The final normalization should see the same sign inversion and those should cancel out. Nevertheless if you got thousands of warning about inconstant pixel size when initializing the integrator, try to change the orientation of the polygon ABCD -> ADCB.

From: pyFAI.detector_factory("Pilatus100k").get_pixel_corners()/172e-6 If one considers the origin at the bottom of the detector, pixels are apparently rotating in clockwise direction, let A be the first corner, the closest from the origin, the second corner, B is on the same x and at larger y, C is in diagonal and the last corner D is on the same y as A but on a different x.

Apparently I wrote something wrong in my previous comment: the order in 3D space appears to be [z,y,x] and not [x,y,z]. There is a detector visualizer in the pyFAI-calib2 written in OpenGL where each pixel is 2 triangle... useful to spot errors.