diff --git a/CHANGELOG.md b/CHANGELOG.md index 537c312..475379a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -80,6 +80,121 @@ A ***critical bug*** was fixed that affects bulk CTR calculations: released versions, up to and including v1.5.0. See the CTR structure-factor documentation for details. +GUI changes: + +- The position readout of the image plot now shows ``HKL`` as a single + bracketed triplet instead of three separate ``H``, ``K`` and ``L`` fields, + and reports the momentum transfer of the selected reciprocal-space frame as + ``Q[alpha]``, ``Q[lab]``, ``Q[omega]``, ``Q[chi]``, ``Q[phi]`` or + ``Q[cryst]``. The field is relabelled when the frame selection changes. + Both triplets are shown with five decimals, and the two fields are sized so + that the full triplet is visible instead of being elided. Pixel coordinates + are shown with two decimals, which is the space this needs. + +GUI fixes: + +- The Q-plot now restores the plot orientation when it is switched off. The + x-axis origin was flipped on entering but never recorded or restored, and + the y-axis restore was tied to the maximum/sum image, so it was skipped + whenever no additional image was displayed. The conversion itself is now + guarded as well: a failure reports an error and leaves the displayed image + and both axes untouched, instead of leaving the plot half-converted with the + scan image hidden. +- The Q-plot ordinate now follows the geometry instead of always pointing + upwards. In an inverted, downward-scattering setup such as the one at ID31 — + recognised from the azimuthal reference, not hardcoded — the out-of-plane + coordinate grows towards the bottom of the detector image, and drawing it + upwards showed the reciprocal-space image mirrored with respect to the + detector image. +- The position readout is now correct while the Q-plot is switched on. The + plot axes hold momentum transfer there and not detector pixels, so the + ``X``/``Y`` labels were wrong and ``HKL``, ``del``, ``gam`` and ``Q`` were + computed as if the cursor were on a pixel. The two axis fields are now + labelled ``q_par[]`` and ``q_perp[]`` and shown with five + decimals, and the remaining fields are obtained by inverting the conversion, + so hovering a point of the Q-plot reports the same ``HKL`` and the same + detector angles as hovering the pixel it came from. The reported ``Q`` is + the reconstructed vector, so its components are exactly the coordinates + under the cursor; it is unrefracted in the Q-plot because the axes are + unrefracted as well, while ``HKL`` stays refraction corrected in both modes. + Positions that no pixel maps to — the corners of the rebinned grid, or + anything outside the Ewald sphere — are shown as ``------``. +- Cancelling or resetting the machine parameter or the crystal parameter + dialog now restores the configuration that was active when the dialog was + opened. Both dialogs apply every edit immediately, so restoring the widgets + alone left the edited values active, and the discarded configuration stayed + in use until it was overwritten or a config file was loaded. + +ESRF ID31 beamline support and reciprocal-space display: + +- Detector images can now be displayed in reciprocal space. The new ``Q-plot`` + toolbar action rebins the currently shown image (single image, maximum, or + sum) onto a regular grid of in-plane and out-of-plane momentum transfer + using pyFAI's ``FiberIntegrator``, which requires pyFAI >= 2025.1. With an + older pyFAI the action reports an error and stays disabled. A drop-down next + to the action selects the reciprocal-space frame: ``Q_alpha`` (the surface + frame returned by ``QAlpha``, the default), ``Q_lab``, ``Q_omega``, + ``Q_chi``, ``Q_phi`` and ``Q_cryst``. The frames that undo the ``omega`` + rotation are only defined for a single image and are refused for maximum and + sum images. **This feature is experimental and its conventions may still + change.** +- The reciprocal-space conversion lives in the new module + ``orgui.app.qconversion`` and is exact for arbitrary azimuthal references. + It is deliberately part of the application layer rather than of + ``orgui.datautils.xrayutils``, so that it is not mistaken for production + reciprocal-space code. pyFAI's ``sample_orientation`` flag can only express + quarter turns and its rotations are composed about fixed axes, neither of + which matches orGUI's continuous azimuth convention, so orGUI supplies + ``FiberIntegrator`` with its own unit definitions. The result agrees with the + per-pixel ``QAlpha`` calculation to numerical precision; see the geometry + documentation for when to use which. +- Added the compiled extension ``_qconversion_cpp``, which converts a full + detector image in a single pass. A 6.2 megapixel Pilatus 6M image is + converted in roughly 50 ms; a numpy implementation is used when the extension + has not been built. The pyFAI integrator is reused while the conversion is + unchanged, so stepping through images no longer resets it for every image. +- HDF5 files that are still being written can now be refreshed from the GUI. + The tree view is rebuilt and the current scan reloaded without restarting + the application. +- Added a ``BlissScan_EBS_p4`` backend for the ID31 Pilatus4 detector, an + example backend under ``examples/backend/ID31_EBS_p4_backend.py``, and the + ``ch8153`` beamtime. +- The default backend is now ``id31_default_p4`` and the fallback geometry + describes a Pilatus4 4M CdTe detector. +- Maximum and sum images are now untoggled when the image number changes, and + the maximum/sum toolbar icons no longer get out of sync with the displayed + image. +- The Q-plot now stays switched on and follows the display instead of being + turned off. It is rebuilt when the image number changes, when the maximum or + sum image is toggled, when the frame is changed, and when the machine angles, + the azimuthal reference, the energy or the orientation matrix are edited. +- The ``Q-plot`` toolbar action now has an icon, a capital Q, instead of being + the only text button in that toolbar. The text is kept as the action name + and the tooltip is unchanged. +- **BREAKING CHANGE:** HDF5 file locking is now disabled by default + (``HDF5_USE_FILE_LOCKING=False``), so that files still open for writing by + the acquisition system can be read. A manually set environment variable + still wins, and ``--hdflocking`` / ``-l`` restores the previous behavior. +- **BREAKING CHANGE:** the segmented ("interlaced") scan loader no longer + guesses a file's layout from a hardcoded list of beamtime ids. Backends can + now answer directly with the new optional ``Scan.listScans`` classmethod, + which returns the scan identifiers -- numbers such as ``[1, 2, 10]``, names + such as ``["ascan_12", "dscan_3"]``, or ``(identifier, label)`` pairs to + label the rows of the selection dialog. Backends that do not implement it + are handled by applying their own ``parse_h5_node`` to every entry of the + file root, which requires no backend change but relies on ``parse_h5_node`` + raising for entries that are not scans. Either way the loader can no longer + disagree with how the same backend opens a single scan. This fixes the + ``id31_default_p4`` backend, which the old list did not cover, and makes + custom backends work regardless of what their class is called -- the example + backend under ``examples/backend/ID31_EBS_p4_backend.py`` was itself + affected. Backends addressed by a name rather than by a number are now + supported here as well, which makes segmented scans work for the legacy + ``ch5523`` beamtime for the first time: it was listed as ID31-style, so the + loader looked for a ``"."`` name its files never had and + raised before showing the dialog. Subscans such as ``1.10`` are no longer + mistaken for the ``1.1`` fast-counter subscan. + ## [1.5.0] (2026-06-07) diff --git a/doc/source/beamline_backends.rst b/doc/source/beamline_backends.rst index 0843c96..9f0aade 100644 --- a/doc/source/beamline_backends.rst +++ b/doc/source/beamline_backends.rst @@ -110,7 +110,9 @@ Required Scan Methods ``parse_h5_node(cls, node)`` Return a dictionary with at least ``scanno`` and ``name`` for a scan node - double clicked in the Nexus tree. + double clicked in the Nexus tree. It is expected to **raise** for a node + that is not a scan: the segmented scan loader uses that to tell scans apart + from the other entries of a file, see below. ``__len__()`` Return the number of images or points in the scan. @@ -123,6 +125,63 @@ Required Scan Methods Optional but convenient. The examples implement this by returning ``get_raw_img(i)``. +How a Scan Is Identified +------------------------ + +orGUI carries two fields for a scan, and a backend may need both when they +differ: + +``scanno`` + The scan number shown in the GUI. It is an integer, because the scan number + selector is a spin box. + +``name`` + A free-form identifier. Most backends do not need it, but a backend whose + constructor takes a *name* rather than a number is opened with this one -- + ``BlissScan`` (beamtime ``ch5523``) reads its scans as ``ascan_12``, + ``dscan_3`` and so on, and uses the name as an HDF5 group key. + +Which of the two is passed to the constructor is decided per beamtime in +``orgui.backend.backends.openScan``. + +Listing Scans for the Segmented Scan Loader +------------------------------------------- + +The segmented ("interlaced") scan loader has to know every scan in a file +*before* any scan object exists, so it cannot use ``parse_h5_node`` on a single +selected node the way normal scan loading does. Backends can answer directly +with the optional ``listScans`` classmethod: + +.. code-block:: python + + from orgui.backend.scans import Scan + + class MyBackend(Scan): + + @classmethod + def listScans(cls, h5file): + return [1, 2, 10] # scan numbers + # or + return ["ascan_12", "dscan_3"] # scan names + # or, with a label for the selection dialog + return [(1, "ascan th 0 90 90 1"), (2, "dscan mu 0 1 20 1")] + +Return whatever your ``__init__`` accepts as its scan identifier. It is passed +back untouched when the scan is opened, and only converted to ``str`` for +display, so a backend addressed by name works exactly like one addressed by +number. Do not use ``float`` for BLISS style ``"."`` names -- +``1.1`` and ``1.10`` are the same float. + +**Implementing this is optional.** A backend that does not implement it gets +the fallback: orGUI applies the backend's own ``parse_h5_node`` to every entry +of the file root, keeps the entries that yield a ``scanno``, and collapses +repeated scan numbers so the subscans of one measurement appear once. That +needs no extra code in the backend and cannot disagree with normal scan +loading, but it has to walk the tree, and it depends on ``parse_h5_node`` +raising for entries that are not scans. + +Nothing here depends on the name of the backend class or on the beamtime id. + Image Objects ------------- diff --git a/doc/source/geometry.rst b/doc/source/geometry.rst index 06483b3..4378995 100644 --- a/doc/source/geometry.rst +++ b/doc/source/geometry.rst @@ -66,6 +66,313 @@ where :math:`\vec{H}_0` and :math:`\vec{H}_1` are vectors in reciprocal lattice units. orGUI calculates the intersections of this line with the Ewald sphere and converts the corresponding detector angles back to detector pixel coordinates. +Momentum Transfer and the Q-Plot +-------------------------------- + +The central plot of orGUI normally shows the detector image in **pixel +coordinates**; the reciprocal-space position of a pixel is reported in +reciprocal lattice units as the cursor is moved. The ``Q-plot`` toolbar action +switches the display to **momentum transfer**, the representation commonly used +for grazing-incidence work. This is a change of the displayed coordinates only, +it does not alter the loaded data. + +Momentum transfer is split into a component parallel to the surface and one +along the surface normal, + +* :math:`q_\parallel` -- in-plane momentum transfer, labelled ``qip`` and + plotted as ``q_par``, +* :math:`q_\perp` -- out-of-plane momentum transfer, labelled ``qoop`` and + plotted as ``q_perp``. + +orGUI can calculate these two quantities in two different ways. Both describe +the same reciprocal space, but they differ in what they return and in when they +are used. + +Per pixel with ``QAlpha`` +~~~~~~~~~~~~~~~~~~~~~~~~~ + +This is the route used everywhere else in orGUI, for example for reflection +positions and for the rod integration. The detector calibration converts pixel +coordinates into the surface angles ``gamma`` and ``delta``, and the Vlieg +diffraction equation converts those angles into a Cartesian momentum transfer +vector in the ``alpha`` frame, + +.. math:: + + \begin{aligned} + Q_x &= K \sin(\delta)\cos(\gamma) \\ + Q_y &= K \left[\cos(\delta)\cos(\gamma) - \cos(\alpha)\right] \\ + Q_z &= K \left[\sin(\gamma) + \sin(\alpha)\right] + \end{aligned} + +with :math:`K = 2\pi/\lambda`. The out-of-plane component is :math:`Q_z`, and +the in-plane component is the radial component in the surface plane: + +.. code-block:: python + + import numpy as np + + detectorCal = ubcalc.detectorCal # DetectorCalibration.Detector2D_SXRD + angles = ubcalc.angles # HKLVlieg.VliegAngles + alpha_i = ubcalc.mu # incidence angle in rad + + gamma, delta = detectorCal.surfaceAngles(alpha_i, shape=image.shape) + Qxyz = angles.QAlpha(alpha_i, delta, gamma) # inverse Angstrom + + q_perp = Qxyz[..., 2] # out-of-plane + q_par = np.hypot(Qxyz[..., 0], Qxyz[..., 1]) # in-plane, unsigned + +The result has the same shape as the detector image: every pixel keeps its own +momentum transfer. The sampling is therefore irregular and curved in +:math:`(q_\parallel, q_\perp)`, which is exactly what is needed to look up +where a reflection falls, but it cannot be handed to an image widget directly. + +On a regular grid with ``FiberIntegrator`` +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +.. warning:: + + The ``Q-plot`` and the :mod:`orgui.app.qconversion` module are + **experimental**. They live in the application layer, not in + ``orgui.datautils.xrayutils``, precisely so that they are not mistaken for + the production reciprocal-space code. The conventions described here may + still change. + +This is the route used by the ``Q-plot`` action. pyFAI's ``FiberIntegrator`` +computes the same per-pixel momentum transfer and additionally **rebins** the +intensities onto a regular grid, which is what makes the result displayable as +an image. orGUI drives it through :mod:`orgui.app.qconversion`, which supplies +pyFAI with orGUI's own angle conventions: + +.. code-block:: python + + from orgui.app import qconversion + + res = qconversion.integrateImage( + detectorCal, + image, + alpha_i, # incidence angle in rad + frame="Q_alpha", # see "Reciprocal-space frames" below + ) + + res.intensity # rebinned image, shape (npt_oop, npt_ip) + res.inplane # q_par axis + res.outofplane # q_perp axis + +``FiberIntegrator`` was introduced in pyFAI 2025.1. With an older pyFAI the +``Q-plot`` action is unavailable and reports an error instead. + +Which one to use +~~~~~~~~~~~~~~~~ + +.. list-table:: + :header-rows: 1 + :widths: 22 39 39 + + * - + - ``QAlpha`` + - ``FiberIntegrator`` + * - Returns + - one :math:`(q_\parallel, q_\perp)` per pixel + - intensities rebinned onto a regular grid + * - Sampling + - irregular, curved + - regular, uniform + * - In-plane sign + - unsigned, :math:`q_\parallel \ge 0` + - signed + * - Requires + - nothing beyond orGUI + - pyFAI >= 2025.1 + * - Use for + - reflection positions, integration, any calculation + - displaying a whole image in reciprocal space + +Use ``QAlpha`` whenever the momentum transfer of a specific pixel or reflection +is needed. Use ``FiberIntegrator`` when a complete image has to be shown or +exported in reciprocal-space coordinates, and accept that the rebinning +quantises positions to the width of one grid cell. + +Reciprocal-space frames +~~~~~~~~~~~~~~~~~~~~~~~ + +``QAlpha`` returns the momentum transfer in the **alpha frame**, the frame of +the sample surface. The drop-down next to the ``Q-plot`` action selects which +frame the image is displayed in; the remaining frames are reached by undoing +the sample rotations, following +:math:`\vec{H} = UB^{-1}\,\Phi^{-1}\,X^{-1}\,\Omega^{-1}\,\vec{Q}_\alpha`. + +.. list-table:: + :header-rows: 1 + :widths: 20 46 34 + + * - Frame + - Reached from the alpha frame by + - Defined for + * - ``Q_alpha`` + - nothing, this is the default + - any image + * - ``Q_lab`` + - the ``alpha`` rotation + - any image + * - ``Q_omega`` + - undoing ``omega`` + - a single image only + * - ``Q_chi`` + - undoing ``omega`` and ``chi`` + - a single image only + * - ``Q_phi`` + - undoing ``omega``, ``chi`` and ``phi`` + - a single image only + * - ``Q_cryst`` + - additionally undoing the orientation matrix ``U`` + - a single image only + +For every frame the out-of-plane component is the ``z`` axis of that frame and +the in-plane component is the radial component in its ``xy`` plane. Maximum and +sum images combine many ``omega`` angles, so the frames that undo ``omega`` are +refused for them. ``Q_cryst`` additionally needs the orientation matrix. + +``Q_cryst`` is the Cartesian reciprocal-space frame of the crystal, so it holds +:math:`B\vec{H}`. Multiplying ``Q_phi`` by :math:`UB^{-1}`, or equivalently +``Q_cryst`` by :math:`B^{-1}`, gives the reciprocal lattice coordinates +:math:`\vec{H} = (h, k, l)^T`, the same result as +:meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.anglesToHkl`. + +Orientation of the Q-plot +~~~~~~~~~~~~~~~~~~~~~~~~~ + +The detector image is displayed with the row index increasing downwards. The +``Q-plot`` keeps that visual orientation, so that a feature seen on the +detector is found in the same place of the reciprocal-space image. + +In the ordinary upward-scattering geometry the out-of-plane coordinate +*de*creases with the row index, and the ordinate points upwards as usual. An +azimuthal reference that points the surface normal the other way — the +inverted, downward-scattering setup used at ID31 — reverses this: ``gamma`` +comes out negative over the whole detector and :math:`q_\perp` grows towards +the bottom of the image. Drawing it upwards anyway would show the +reciprocal-space image mirrored, so orGUI inverts the ordinate in that case. + +The direction is derived from the detector calibration by +:func:`~orgui.app.qconversion.outOfPlaneIncreasesWithRow`, which compares the +out-of-plane coordinate of the top and the bottom of the detector. No beamline +is special-cased. The abscissa always runs from left to right. + +The position readout in the Q-plot +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +While the ``Q-plot`` is on, the abscissa and the ordinate hold momentum +transfer instead of detector pixels. The readout below the plot therefore +labels them ``q_par[]`` and ``q_perp[]``, and obtains the +remaining fields by inverting the conversion, so that a point of the Q-plot +reports the same ``HKL``, ``del``, ``gam`` and ``Q`` as the pixel it was +rebinned from. + +The projection keeps :math:`Q_z` and collapses :math:`Q_x` and :math:`Q_y` +into their radial distance, signed by :math:`Q_x`. The discarded direction is +recovered from the Ewald condition: with :math:`\vec{q} = \vec{Q}/K` and the +incident beam direction :math:`\vec{c}`, the outgoing direction +:math:`\vec{q} + \vec{c}` is a unit vector, so + +.. math:: + + |\vec{q}|^2 + 2\,\vec{q}\cdot\vec{c} = 0, + +which is *linear* in :math:`q_x` and :math:`q_y` once :math:`q_z` and +:math:`|\vec{q}|` are known. Intersecting that line with the circle +:math:`q_x^2 + q_y^2 = (q_\parallel/K)^2` leaves at most two solutions, of +which only those with :math:`\mathrm{sign}(Q_x) = \mathrm{sign}(q_\parallel)` +reproduce the sign convention above. + +In the ``Q_alpha`` frame the beam is :math:`(0, \cos\alpha, -\sin\alpha)` and +has no ``x`` component, so the line is parallel to the ``x`` axis and the sign +condition always leaves a single solution. In the rotated frames the beam +generally does have one and both solutions can carry the same sign; orGUI then +keeps the one that falls onto the detector. That test doubles as the check +whether the position is covered by the image at all, because the rebinned grid +is rectangular and reaches beyond the detector. Positions that no pixel maps +to are shown as ``------``. + +The reported ``Q`` is the reconstructed vector itself, so its in-plane and +out-of-plane components are exactly the coordinates the cursor is at. Unlike +in pixel mode it carries no refraction correction, because the Q-plot axes +carry none either; ``HKL`` is refraction corrected in both modes, as +everywhere else in orGUI. + +Why orGUI supplies its own pyFAI units +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +pyFAI's built-in fiber units cannot express orGUI's geometry, for two +independent reasons. + +**The azimuth is continuous, ``sample_orientation`` is not.** pyFAI's +``sample_orientation`` is the EXIF dihedral group acting on the image array, so +it only spans quarter turns. It describes how the image array is laid out, not +where the sample sits. orGUI's azimuthal reference rotates the ``alpha`` +rotation axis about the incident beam and takes arbitrary values, so it is +passed through ``tilt_angle`` instead, while ``sample_orientation`` stays at the +identity: + +.. math:: + + \mathrm{tilt\_angle} = -\left(\mathrm{azimuth} + \frac{\pi}{2}\right) + +The quarter turn is the offset between orGUI's azimuthal reference and pyFAI's +fiber convention. Encoding the azimuth in ``sample_orientation`` instead would +only be correct for azimuths that are exact multiples of 90 degrees. + +**The rotations are composed in the opposite order.** pyFAI combines the +incidence and tilt rotations extrinsically, about fixed axes, as +:math:`R_x(-\mathrm{tilt})\,R_y(\mathrm{incidence})`. In orGUI the azimuth +rotates the ``alpha`` axis itself, so the incidence rotation has to be applied +about the *already rotated* axis, +:math:`R_y(\mathrm{incidence})\,R_x(-\mathrm{tilt})`. The two agree only when +one of the two angles vanishes. + +:mod:`~orgui.app.qconversion` therefore builds its own +``pyFAI.units.UnitFiber`` objects that implement orGUI's convention, and hands +them to ``FiberIntegrator``. pyFAI still performs the rebinning; only the +coordinate definition is replaced. + +With that in place the two routes agree to numerical precision. This is +verified in ``orgui/app/test/test_q_conversion.py``, which compares them per +pixel for flat and tilted detectors, at several incidence angles and at +azimuths that are deliberately not multiples of 90 degrees. The same file +checks that ``Q_phi`` reproduces ``anglesToHkl``, that ``Q_cryst`` reproduces +:math:`B\vec{H}`, that every frame preserves :math:`|\vec{Q}|`, and that a +single bright pixel is placed by the rebinning within one grid cell of the +position predicted by ``QAlpha``. + +Performance +~~~~~~~~~~~ + +Detector images are large, so the conversion is collapsed into a single affine +relation before any pixel data is touched. With +:math:`n = |(x, y, z)|` every component reduces to + +.. math:: + + q_j = k \left( \frac{G_{j0} x + G_{j1} y + G_{j2} z}{n} - c_j \right) + +where the matrix :math:`G` and the offset :math:`\vec{c}` absorb the sample +orientation map, the beam and incidence rotations, the axis relabelling and the +frame rotation. A compiled kernel, +``orgui/app/cpp/qconversion_cpp.cpp``, evaluates both displayed quantities in a +single pass; a plain numpy implementation is used when the extension has not +been built. The result is cached, because pyFAI evaluates the in-plane and the +out-of-plane unit separately but passes the same pixel positions to both, so an +image is converted once rather than twice. + +On a 6.2 megapixel Pilatus 6M the compiled kernel converts a full image in +roughly 50 ms, and the second unit evaluation is served from the cache. + +.. note:: + + ``numexpr`` is deliberately not used here. Version 2.11.0 returns wrong + results for a small fraction of multi-threaded evaluations of expressions of + this kind, which is not acceptable for a coordinate transform. + Further Reading --------------- diff --git a/doc/source/release_notes.rst b/doc/source/release_notes.rst index f14a668..623797d 100644 --- a/doc/source/release_notes.rst +++ b/doc/source/release_notes.rst @@ -31,6 +31,31 @@ A ***critical bug*** was fixed that affects bulk CTR calculations: - ``UnitCell.F_bulk``'s semi-infinite geometric lattice sum used the raw, untransformed ``l`` index instead of the index converted by ``refHKLTransform`` when computing the out-of-plane attenuation phase. This was only correct for the default case where a component uses its own bulk cell as the reference (no ``reference_uc`` set); any explicit ``reference_uc`` whose out-of-plane reciprocal axis differs from the bulk's — including a plain scale difference between the reference and bulk out-of-plane axis length, not only a rotated or reindexed reference — gave incorrect bulk structure-factor amplitudes. This bug was present in both the accelerated (numba/C++) and plain-Python code paths in all previous released versions, up to and including v1.5.0. See the CTR structure-factor documentation for details. +GUI changes: + +- The position readout of the image plot now shows ``HKL`` as a single bracketed triplet instead of three separate ``H``, ``K`` and ``L`` fields, and reports the momentum transfer of the selected reciprocal-space frame as ``Q[alpha]``, ``Q[lab]``, ``Q[omega]``, ``Q[chi]``, ``Q[phi]`` or ``Q[cryst]``. The field is relabelled when the frame selection changes. Both triplets are shown with five decimals, and the two fields are sized so that the full triplet is visible instead of being elided. Pixel coordinates are shown with two decimals, which is the space this needs. + +GUI fixes: + +- The Q-plot now restores the plot orientation when it is switched off. The x-axis origin was flipped on entering but never recorded or restored, and the y-axis restore was tied to the maximum/sum image, so it was skipped whenever no additional image was displayed. The conversion itself is now guarded as well: a failure reports an error and leaves the displayed image and both axes untouched, instead of leaving the plot half-converted with the scan image hidden. +- The Q-plot ordinate now follows the geometry instead of always pointing upwards. In an inverted, downward-scattering setup such as the one at ID31 — recognised from the azimuthal reference, not hardcoded — the out-of-plane coordinate grows towards the bottom of the detector image, and drawing it upwards showed the reciprocal-space image mirrored with respect to the detector image. +- The position readout is now correct while the Q-plot is switched on. The plot axes hold momentum transfer there and not detector pixels, so the ``X``/``Y`` labels were wrong and ``HKL``, ``del``, ``gam`` and ``Q`` were computed as if the cursor were on a pixel. The two axis fields are now labelled ``q_par[]`` and ``q_perp[]`` and shown with five decimals, and the remaining fields are obtained by inverting the conversion, so hovering a point of the Q-plot reports the same ``HKL`` and the same detector angles as hovering the pixel it came from. The reported ``Q`` is the reconstructed vector, so its components are exactly the coordinates under the cursor; it is unrefracted in the Q-plot because the axes are unrefracted as well, while ``HKL`` stays refraction corrected in both modes. Positions that no pixel maps to — the corners of the rebinned grid, or anything outside the Ewald sphere — are shown as ``------``. +- Cancelling or resetting the machine parameter or the crystal parameter dialog now restores the configuration that was active when the dialog was opened. Both dialogs apply every edit immediately, so restoring the widgets alone left the edited values active, and the discarded configuration stayed in use until it was overwritten or a config file was loaded. + +ESRF ID31 beamline support and reciprocal-space display: + +- Detector images can now be displayed in reciprocal space. The new ``Q-plot`` toolbar action rebins the currently shown image (single image, maximum, or sum) onto a regular grid of in-plane and out-of-plane momentum transfer using pyFAI's ``FiberIntegrator``, which requires pyFAI >= 2025.1. With an older pyFAI the action reports an error and stays disabled. A drop-down next to the action selects the reciprocal-space frame: ``Q_alpha`` (the surface frame returned by ``QAlpha``, the default), ``Q_lab``, ``Q_omega``, ``Q_chi``, ``Q_phi`` and ``Q_cryst``. The frames that undo the ``omega`` rotation are only defined for a single image and are refused for maximum and sum images. **This feature is experimental and its conventions may still change.** +- The reciprocal-space conversion lives in the new module ``orgui.app.qconversion`` and is exact for arbitrary azimuthal references. It is deliberately part of the application layer rather than of ``orgui.datautils.xrayutils``, so that it is not mistaken for production reciprocal-space code. pyFAI's ``sample_orientation`` flag can only express quarter turns and its rotations are composed about fixed axes, neither of which matches orGUI's continuous azimuth convention, so orGUI supplies ``FiberIntegrator`` with its own unit definitions. The result agrees with the per-pixel ``QAlpha`` calculation to numerical precision; see the geometry documentation for when to use which. +- Added the compiled extension ``_qconversion_cpp``, which converts a full detector image in a single pass. A 6.2 megapixel Pilatus 6M image is converted in roughly 50 ms; a numpy implementation is used when the extension has not been built. The pyFAI integrator is reused while the conversion is unchanged, so stepping through images no longer resets it for every image. +- HDF5 files that are still being written can now be refreshed from the GUI. The tree view is rebuilt and the current scan reloaded without restarting the application. +- Added a ``BlissScan_EBS_p4`` backend for the ID31 Pilatus4 detector, an example backend under ``examples/backend/ID31_EBS_p4_backend.py``, and the ``ch8153`` beamtime. +- The default backend is now ``id31_default_p4`` and the fallback geometry describes a Pilatus4 4M CdTe detector. +- Maximum and sum images are now untoggled when the image number changes, and the maximum/sum toolbar icons no longer get out of sync with the displayed image. +- The Q-plot now stays switched on and follows the display instead of being turned off. It is rebuilt when the image number changes, when the maximum or sum image is toggled, when the frame is changed, and when the machine angles, the azimuthal reference, the energy or the orientation matrix are edited. +- The ``Q-plot`` toolbar action now has an icon, a capital Q, instead of being the only text button in that toolbar. The text is kept as the action name and the tooltip is unchanged. +- **BREAKING CHANGE:** HDF5 file locking is now disabled by default (``HDF5_USE_FILE_LOCKING=False``), so that files still open for writing by the acquisition system can be read. A manually set environment variable still wins, and ``--hdflocking`` / ``-l`` restores the previous behavior. +- **BREAKING CHANGE:** the segmented ("interlaced") scan loader no longer guesses a file's layout from a hardcoded list of beamtime ids. Backends can now answer directly with the new optional ``Scan.listScans`` classmethod, which returns the scan identifiers -- numbers such as ``[1, 2, 10]``, names such as ``["ascan_12", "dscan_3"]``, or ``(identifier, label)`` pairs to label the rows of the selection dialog. Backends that do not implement it are handled by applying their own ``parse_h5_node`` to every entry of the file root, which requires no backend change but relies on ``parse_h5_node`` raising for entries that are not scans. Either way the loader can no longer disagree with how the same backend opens a single scan. This fixes the ``id31_default_p4`` backend, which the old list did not cover, and makes custom backends work regardless of what their class is called -- the example backend under ``examples/backend/ID31_EBS_p4_backend.py`` was itself affected. Backends addressed by a name rather than by a number are now supported here as well, which makes segmented scans work for the legacy ``ch5523`` beamtime for the first time: it was listed as ID31-style, so the loader looked for a ``"."`` name its files never had and raised before showing the dialog. Subscans such as ``1.10`` are no longer mistaken for the ``1.1`` fast-counter subscan. + 1.5.0 (2026-06-07) ------------------ diff --git a/examples/backend/CHESS_QM2.py b/examples/backend/CHESS_QM2.py index 21e2a0e..e77d2ab 100644 --- a/examples/backend/CHESS_QM2.py +++ b/examples/backend/CHESS_QM2.py @@ -239,6 +239,43 @@ def parse_h5_node(cls, obj): ddict['name'] = obj.local_name return ddict + @classmethod + def listScans(cls, h5file): + """list every scan in the file, for the segmented scan loader. + + Optional. Without it orGUI calls parse_h5_node on every entry of the + file root instead, which gives the same answer but has to walk the + tree. The file is already open and is handed in, so unlike __init__ + this does not have to deal with a path or a node. + + Every group of the root is a scan, named as in SPEC files, i.e. + 1.1, 1.2, 1.3, ... The identifier is the scan number in front of the + dot, because that is what __init__ takes. It looks up the first group + whose scan number matches, so a scan number is reported once even when + the file holds several subscans for it - offering the others would let + one be selected and a different one be opened. + """ + entries = [] + seen = set() + for name in h5file: + name = str(name) + try: + scanno = int(name.split('.')[0]) + except ValueError: + continue # not a scan, e.g. an instrument or metadata group + if scanno in seen: + continue + seen.add(scanno) + try: + title = h5file[name + '/title'][()] + except Exception: + title = name + if isinstance(title, bytes): + title = title.decode() + entries.append((scanno, str(title))) + entries.sort(key=lambda entry: entry[0]) + return entries + # returns single image, attempts to use the default camera! def get_raw_img(self,img, **kwargs): imagepath = os.path.join(self.imagefolder, self.image_prefix + ("%05i" % img) + self.image_suffix) diff --git a/examples/backend/ID31_EBS_p4_backend.py b/examples/backend/ID31_EBS_p4_backend.py new file mode 100644 index 0000000..46989e6 --- /dev/null +++ b/examples/backend/ID31_EBS_p4_backend.py @@ -0,0 +1,539 @@ +# /*########################################################################## +# +# Copyright (c) 2020-2026 Timo Fuchs, Finn Schröter +# +# Permission is hereby granted, free of charge, to any person obtaining a copy +# of this software and associated documentation files (the "Software"), to deal +# in the Software without restriction, including without limitation the rights +# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +# copies of the Software, and to permit persons to whom the Software is +# furnished to do so, subject to the following conditions: +# +# The above copyright notice and this permission notice shall be included in +# all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +# THE SOFTWARE. +# +# ###########################################################################*/ +"""Example ID31 EBS backend for Pilatus 4 / 4M HDF5 scans. + +This file is meant to be loaded through an orGUI config ``[backend] file = ...`` +entry or through the "Load backend file" GUI action. It intentionally contains +exactly one :class:`orgui.backend.scans.Scan` subclass. +""" + +import copy +import os +import traceback +import warnings + +import numpy as np +import scipy.interpolate +from silx.io import dictdump +import silx.io.h5py_utils + +from orgui.backend.scans import Scan + + +class ID31DiffractLinTilt: + """Convert the ID31 ``linai`` linear translation to incident angle ``mu``. + + :param float muoffset: + Experiment-specific alignment offset in deg. + """ + + def __init__(self, muoffset=-0.06442416994811659): + self.config = { + "a": 938, + "b": 400, + "muoffset": muoffset, + } + + @property + def a(self): + """Return linkage length ``a`` in mm.""" + return self.config.get("a", float) + + @property + def b(self): + """Return linkage length ``b`` in mm.""" + return self.config.get("b", float) + + @property + def muoffset(self): + """Return the experiment-specific ``mu`` offset in deg.""" + return self.config.get("muoffset", float) + + def c(self, a=None, b=None): + """Return the derived linkage length in mm.""" + a = self.a if a is None else a + b = self.b if b is None else b + return np.sqrt(np.square(a) + np.square(b)) + + def d(self, a=None, b=None): + """Return the derived linkage angle in rad.""" + a = self.a if a is None else a + b = self.b if b is None else b + return np.arctan(b / a) + + def calc_from_real(self, positions_dict): + """Convert linear translation in mm to tilt angle in deg. + + :param dict positions_dict: + Must contain ``{"linear": ...}``. + :returns: + ``{"tilt": ...}`` with tilt in deg. + :rtype: dict + """ + a, b = self.a, self.b + c, d = self.c(a, b), self.d(a, b) + a2, c2 = np.square(a), np.square(c) + linear = positions_dict["linear"] + bc2 = np.square(b + linear) + tilt = np.arccos((a2 + c2 - bc2) / (2 * a * c)) - d + return {"tilt": np.rad2deg(tilt)} + + def calc_to_real(self, positions_dict): + """Convert tilt angle in deg to linear translation in mm. + + :param dict positions_dict: + Must contain ``{"tilt": ...}``. + :returns: + ``{"linear": ...}`` with translation in mm. + :rtype: dict + """ + a, b = self.a, self.b + c, d = self.c(a, b), self.d(a, b) + a2, c2 = np.square(a), np.square(c) + tilt = np.deg2rad(positions_dict["tilt"]) + bc = np.sqrt(a2 + c2 - 2 * a * c * np.cos(tilt + d)) + return {"linear": bc - b} + + def linai_to_mu(self, linai): + """Convert ``linai`` in mm to incident angle ``mu`` in deg.""" + tilt = self.calc_from_real({"linear": linai}) + return tilt["tilt"] - self.muoffset + + def mu_to_linai(self, mu): + """Convert incident angle ``mu`` in deg to ``linai`` in mm.""" + tilts = {"tilt": mu + self.muoffset} + return self.calc_to_real(tilts)["linear"] + + +class h5_Image: + """Minimal image container returned by the backend scan API.""" + + def __init__(self, data): + """Store detector image data as a NumPy-like array.""" + self.img = data + self.motors = {} + self.counters = {} + + +class ID31_EBS_p4_2026(Scan): + """Load ID31 EBS Pilatus 4 / 4M scans from BLISS HDF5 files. + + :param hdffilepath_orNode: + Either a path to an HDF5 file or an opened silx HDF5 node. + :param int scanno: + Scan number to load. + :param bool loadimg: + If true, load all detector frames into memory. If false, images are + read lazily from HDF5 in :meth:`get_raw_img`. + :param bool saveh5data: + If true, keep the parsed HDF5 dictionaries on the scan object. + :param float muoffset: + Experiment-specific ``linai`` to ``mu`` offset in deg. The default + value ``-2.31`` matches the CH8153 branch configuration. + """ + + def __init__( + self, + hdffilepath_orNode=None, + scanno=None, + loadimg=True, + saveh5data=False, + muoffset=-2.31, + ): + self.cameras = [] + self.loadimg = loadimg + self.muoffset = muoffset + data_1 = None + data_2 = None + excludenames = None if self.loadimg else ["p4_lima1"] + + if hdffilepath_orNode is None: + return + self.hdffilepath_orNode = hdffilepath_orNode + + if isinstance(hdffilepath_orNode, str): + filepath, _filename = os.path.split(hdffilepath_orNode) + _unused, filename_noext = os.path.split(filepath) + self.filename_base = filename_noext + with silx.io.h5py_utils.File(hdffilepath_orNode) as h5file: + data_1, data_2 = self._read_scan_groups( + h5file, scanno, excludenames + ) + else: + hdffilepath = hdffilepath_orNode.local_filename + filepath, _filename = os.path.split(hdffilepath) + _unused, filename_noext = os.path.split(filepath) + self.filename_base = filename_noext + data_1, data_2 = self._read_scan_groups( + hdffilepath_orNode.file, scanno, excludenames + ) + + if saveh5data: + self.data_1 = data_1 + self.data_2 = data_2 + + if "title" in data_1: + self.title = data_1["title"] + + self.positioners = data_1["instrument"]["positioners"] + self._read_detector_data(data_1) + self._read_scan_axis(data_1) + self._read_metadata_counters(data_1, data_2) + + self.offsetindex = 0 + self.scandatapoints = self.nopoints + self.imageno = np.arange(self.nopoints) + self.omega = -1 * self.th + self.axis = getattr(self, self.axisname) + + def _read_scan_groups(self, h5file, scanno, excludenames): + for group_name in h5file: + scansuffix = group_name.split("_")[-1] + scanname_nosuffix = "_".join(group_name.split("_")[:-1]) + scanno_s, _subscanno = scansuffix.split(".") + if int(scanno_s) == scanno: + break + else: + raise OSError(f"Scan number {scanno} not found in file") + + self.scanno1 = str(scanno) + "." + "1" + self.scanno2 = str(scanno) + "." + "2" + self.scanname_1 = scanname_nosuffix + "_" + self.scanno1 + self.scanname_2 = scanname_nosuffix + "_" + self.scanno2 + self.scanname = self.scanname_1 + self.name = self.scanname + + if self.scanname_1 in h5file: + data_1 = dictdump.h5todict( + h5file, self.scanname_1, exclude_names=excludenames + ) + if self.scanname_2 in h5file: + data_2 = dictdump.h5todict( + h5file, self.scanname_2, exclude_names=excludenames + ) + else: + data_2 = None + else: + data_1 = dictdump.h5todict( + h5file, self.scanno1, exclude_names=excludenames + ) + if self.scanno2 in h5file: + data_2 = dictdump.h5todict( + h5file, self.scanno2, exclude_names=excludenames + ) + else: + data_2 = None + + measurement = data_1["measurement"] + if "p4_lima1" in measurement: + self.cameras.append("p4_lima1") + self.nopoints = measurement["p4_lima1"].shape[0] + if "mpx" in measurement: + self.cameras.append("mpx") + self.nopoints = measurement["mpx"].shape[0] + + return data_1, data_2 + + def _read_detector_data(self, data_1): + measurement = data_1["measurement"] + if "p4_lima1" in measurement: + if "p4_lima1" not in self.cameras: + self.cameras.append("p4_lima1") + self.nopoints = measurement["p4_lima1"].shape[0] + if self.loadimg: + self.p4 = measurement["p4_lima1"][()] + + if "mpx" in measurement: + self.mpx = measurement["mpx"][()] + self.nopoints = measurement["mpx"].shape[0] + if "mpx" not in self.cameras: + self.cameras.append("mpx") + + def _read_scan_axis(self, data_1): + measurement = data_1["measurement"] + if "th" in measurement: + self.th = measurement["th"][: self.nopoints] + if "th_trig" in measurement and "th_delta" in measurement: + self.th = ( + measurement["th_trig"][: self.nopoints] + + measurement["th_delta"][: self.nopoints] / 2 + ) + self.axisname = "th" + self.mu = self.positioners["mu"] + elif "nth" in measurement: + self.th = measurement["nth"][: self.nopoints] + self.axisname = "th" + self.mu = self.positioners["nai"] * -1 + elif "uth" in measurement: + self.th = measurement["uth"][: self.nopoints] + self.axisname = "th" + self.mu = self.positioners["mu"] * -1 + elif "mu" in measurement: + self.mu = measurement["mu"][: self.nopoints] + self.axisname = "mu" + self.th = self.positioners["th"] + elif "nai" in measurement: + self.mu = measurement["nai"][: self.nopoints] + self.axisname = "mu" + self.th = self.positioners["nth"] + elif "linai" in measurement: + lintomu = ID31DiffractLinTilt(muoffset=self.muoffset) + self.linai = measurement["linai"][: self.nopoints] + if "linai_trig" in measurement and "linai_delta" in measurement: + self.linai = ( + measurement["linai_trig"][: self.nopoints] + + measurement["linai_delta"][: self.nopoints] / 2 + ) + self.mu = lintomu.linai_to_mu(self.linai) + self.axisname = "mu" + self.th = self.positioners["th"] + else: + self.axisname = "time" + self.th = self.positioners["th"] + self.mu = self.positioners["mu"] + + def _read_metadata_counters(self, data_1, data_2): + measurement = data_1["measurement"] + if "potv" in measurement: + self.potv = measurement["potv"][: self.nopoints] + if "scaled_potv2f" in measurement: + self.scaled_potv2f = measurement["scaled_potv2f"][: self.nopoints] + + if "srcur" in measurement: + self.srcur = measurement["srcur"][: self.nopoints] + else: + self.srcur = np.ones(self.nopoints) + + if "timer_trig" in measurement: + self.time = measurement["timer_trig"][: self.nopoints] + elif "elapsed_time" in measurement: + self.time = measurement["elapsed_time"][: self.nopoints] + else: + raise OSError("Cannot find time counter in scan") + + if "epoch_trig" in measurement: + self.epoch = measurement["epoch_trig"][: self.nopoints] + elif "epoch" in measurement: + self.epoch = measurement["epoch"][: self.nopoints] + else: + raise OSError("Cannot find epoch counter in scan") + + self._read_slow_counters(data_1, data_2) + + if "mondio" in measurement: + self.mondio = measurement["mondio"][: self.nopoints] + if not hasattr(self, "mondio"): + self.mondio = np.ones(self.nopoints) + + self.mondio = self.mondio / np.mean(self.mondio) + self.srcur = self.srcur / np.mean(self.srcur) + + if "timer_delta" in measurement: + self.exposure_time = measurement["timer_delta"][: self.nopoints] + elif "sec" in measurement: + self.exposure_time = measurement["sec"][: self.nopoints] + else: + raise OSError("Cannot find exposure time in scan") + self.relative_exposure = self.exposure_time / np.mean(self.exposure_time) + + def _read_slow_counters(self, data_1, data_2): + if data_2 is not None: + for counter in ("mondio", "potential", "current"): + if counter in data_2["measurement"]: + try: + interpolator = scipy.interpolate.interp1d( + data_2["measurement"]["epoch"], + data_2["measurement"][counter], + ) + setattr(self, counter, interpolator(self.epoch)) + except ValueError: + warnings.warn( + "Cannot interpolate " + f"{counter} in scan {self.scanno2}\n" + f"{traceback.format_exc()}" + ) + setattr(self, counter, data_2["measurement"][counter]) + else: + measurement = data_1["measurement"] + if "potential" in measurement: + self.potential = measurement["potential"][: self.nopoints] + if "current" in measurement: + self.current = measurement["current"][: self.nopoints] + + @classmethod + def parse_h5_node(cls, obj): + """Parse the scan number from an ID31 BLISS HDF5 node. + + :param obj: + silx HDF5 tree object selected by the GUI. + :returns: + Dictionary with ``scanno`` and ``name`` keys. + :rtype: dict + """ + scanname = obj.local_name + if "_" in scanname: + scansuffix = scanname.split("_")[-1] + elif "/" in scanname: + scansuffix = scanname.split("/")[-1] + else: + scansuffix = scanname + scanno, _subscanno = scansuffix.split(".") + return {"scanno": int(scanno), "name": obj.local_name} + + @classmethod + def listScans(cls, h5file): + """List the scans in the file, for the segmented scan loader. + + Optional. Without it orGUI falls back to calling + :meth:`parse_h5_node` on every entry of the file root, which gives the + same answer but has to walk the tree. Implemented here because the + BLISS naming makes it a few lines: keep the ``.1`` subscan, which is + the one holding the fast counters, and report the scan number in front + of the dot together with the scan title as a label. + + :param h5file: an open h5py-like file object. + :returns: ``(scanno, title)`` pairs, sorted by scan number. + :rtype: list + """ + entries = [] + for name in h5file: + name = str(name) + scan, dot, subscan = name.rpartition(".") + if not dot or subscan != "1": + continue + try: + scanno = int(scan.split("_")[-1]) + except ValueError: + continue + try: + title = h5file[name + "/title"][()] + except Exception: + title = name + if isinstance(title, bytes): + title = title.decode() + entries.append((scanno, str(title))) + entries.sort(key=lambda entry: entry[0]) + return entries + + @property + def auxillary_counters(self): + """Return counters copied to the orGUI integration database.""" + return [ + "current", + "potential", + "exposure_time", + "elapsed_time", + "time", + "srcur", + "mondio", + "epoch", + "scaled_potv2f", + ] + + def get_raw_img(self, img): + """Return detector image ``img`` as an :class:`h5_Image`. + + :param int img: + Zero-based image index. + """ + if not hasattr(self, "p4"): + if isinstance(self.hdffilepath_orNode, str): + with silx.io.h5py_utils.File(self.hdffilepath_orNode, "r") as h5file: + if self.scanname_1 in h5file: + data_1 = h5file[self.scanname_1] + else: + data_1 = h5file[self.scanno1] + image = data_1["measurement"]["p4_lima1"][img][()] + else: + h5file = self.hdffilepath_orNode.file + if self.scanname_1 in h5file: + data_1 = h5file[self.scanname_1] + else: + data_1 = h5file[self.scanno1] + image = data_1["measurement"]["p4_lima1"][img][()] + else: + image = self.p4[img] + return h5_Image(image) + + def slice(self, startno, endno): + """Return a scan object containing images ``startno`` through ``endno``. + + The returned object preserves the original scan metadata and slices + NumPy arrays along their first dimension. + """ + if startno > endno: + startno, endno = endno, startno + + fscan = ID31_EBS_p4_2026() + for key, value in self.__dict__.items(): + if isinstance(value, np.ndarray) and value.ndim > 0: + fscan.__dict__[key] = copy.deepcopy(value[startno:endno]) + else: + fscan.__dict__[key] = copy.deepcopy(value) + + if self.loadimg is False: + with silx.io.h5py_utils.File(self.hdffilepath_orNode, "r") as h5file: + if self.scanname_1 in h5file: + data_1 = h5file[self.scanname_1] + else: + data_1 = h5file[self.scanno1] + fscan.p4 = data_1["measurement"]["p4_lima1"][startno:endno][()] + + fscan.offsetindex = copy.deepcopy(startno) + fscan.nopoints = copy.deepcopy(endno - startno) + return fscan + + def get_p3_img(self, img): + """Return an image with ID31 counters attached for normalization.""" + imgdata = self.get_raw_img(img) + imgdata.counters["Time"] = self.time[img] + imgdata.counters["exposure"] = self.exposure_time[img] + imgdata.counters["TrigTime"] = self.relative_exposure[img] + imgdata.counters["imageno"] = self.imageno[img] + + if self.axisname == "th": + imgdata.counters["th"] = self.th[img] + imgdata.counters["om"] = self.omega[img] + else: + imgdata.counters["th"] = self.th + imgdata.counters["om"] = self.omega + + if self.axisname == "mu": + imgdata.counters["mu"] = self.mu[img] + else: + imgdata.counters["mu"] = self.mu + + if self.srcur is not None: + imgdata.counters["srcur"] = self.srcur[img] + if self.mondio is not None: + imgdata.counters["mondio"] = self.mondio[img] + return imgdata + + def __getitem__(self, key): + """Return image ``key`` with ID31 counters attached.""" + return self.get_p3_img(key) + + def __len__(self): + """Return the number of detector images in the scan.""" + return self.nopoints diff --git a/meson.build b/meson.build index 80e95da..54b800b 100644 --- a/meson.build +++ b/meson.build @@ -225,6 +225,14 @@ py.extension_module( subdir: 'orgui/app' ) +py.extension_module( + '_qconversion_cpp', + files('orgui/app/cpp/qconversion_cpp.cpp'), + dependencies: [pybind11_dep], + install: true, + subdir: 'orgui/app' +) + install_subdir( 'orgui', install_dir: py.get_install_dir(), diff --git a/orgui/app/QScanSelector.py b/orgui/app/QScanSelector.py index 379e322..b0079dd 100644 --- a/orgui/app/QScanSelector.py +++ b/orgui/app/QScanSelector.py @@ -79,6 +79,7 @@ class QScanSelector(qt.QMainWindow): sigROIChanged = qt.pyqtSignal() sigROIintegrate = qt.pyqtSignal() sigSearchHKL = qt.pyqtSignal(list) + sigRefreshH5 = qt.pyqtSignal() def __init__(self, parentmainwindow, parent=None): qt.QMainWindow.__init__(self, parent=None) @@ -107,7 +108,7 @@ def __init__(self, parentmainwindow, parent=None): self.hdfTreeView = silx.gui.hdf5.Hdf5TreeView(self) self.hdfTreeView.setSortingEnabled(True) self.hdfTreeView.addContextMenuCallback(self.nexus_treeview_callback) - self.hdf5model = Hdf5TreeModel(self.hdfTreeView, ownFiles=True) + self.hdf5model = Hdf5TreeModel(self.hdfTreeView, ownFiles=False) self.hdfTreeView.setModel(self.hdf5model) self.hdfTreeView.setExpandsOnDoubleClick(False) self.hdf5model.setFileMoveEnabled(True) @@ -170,7 +171,9 @@ def __init__(self, parentmainwindow, parent=None): qt.QLabel("Backend:", btidsplit) self.btid = qt.QComboBox(btidsplit) [self.btid.addItem(bt) for bt in backends.fscans] - self.btid.setCurrentText("id31_default") + # keep the manual selection in sync with the beamtime autodetect + # fallback, so both paths default to the same backend + self.btid.setCurrentText(backends.default_beamtime) self._selectBackendBtn = qt.QPushButton("...", btidsplit) width = self._selectBackendBtn.fontMetrics().boundingRect(" ... ").width() + 7 @@ -1184,6 +1187,12 @@ def __getRelativePath(self, model, rootIndex, index): raise ValueError("index is not a children of the rootIndex") def _onRefreshFile(self): + # The refresh is driven by orGUI, which has to drop the fscan object + # (to avoid broken links) and possibly scan_image and + # currentAddImageLabel (to avoid a crash) before the tree is rebuilt. + self.sigRefreshH5.emit() + + def _onDoRefresh(self): qt.QApplication.setOverrideCursor(qt.Qt.WaitCursor) selection = self.hdfTreeView.selectionModel() @@ -1224,6 +1233,7 @@ def _onRefreshFile(self): model.clear() """ import gc + gc.collect() # to collect all not anymore used/bound objects for obj in gc.get_objects(): # Browse through ALL objects if isinstance(obj, Hdf5TreeModel): try: @@ -1238,6 +1248,7 @@ def _onRefreshFile(self): pass # Was already closed """ + # why create a new tree view? silx view handles resync differently... self.createTreeView() maintab = self.mainwidget.findChildren(qt.QTabWidget)[0] @@ -1249,6 +1260,7 @@ def _onRefreshFile(self): for h5, filename in h5files: modelnew.insertFile(filename, 0) + # self.hdfTreeView.expandToDepth(0) # this call here can cause the nodes to have broken links!!! # noqa: E501 # self.__expandNodesFromPaths(self.hdfTreeView, index, paths) # self.hdf5model.appendFile(filename) qt.QApplication.restoreOverrideCursor() @@ -1257,12 +1269,12 @@ def createTreeView(self): self.hdfTreeView = silx.gui.hdf5.Hdf5TreeView(self) self.hdfTreeView.setSortingEnabled(True) self.hdfTreeView.addContextMenuCallback(self.nexus_treeview_callback) - self.hdf5model = Hdf5TreeModel(self.hdfTreeView, ownFiles=True) + self.hdf5model = Hdf5TreeModel(self.hdfTreeView, ownFiles=False) self.hdfTreeView.setModel(self.hdf5model) self.hdfTreeView.setExpandsOnDoubleClick(False) self.hdf5model.setFileMoveEnabled(True) - self.hdf5model.sigH5pyObjectLoaded.connect(self.__h5FileLoaded) + # self.hdf5model.sigH5pyObjectLoaded.connect(self.__h5FileLoaded) self.hdf5model.sigH5pyObjectRemoved.connect(self.__h5FileRemoved) self.hdf5model.sigH5pyObjectSynchronized.connect(self.__h5FileSynchonized) diff --git a/orgui/app/QUBCalculator.py b/orgui/app/QUBCalculator.py index 3b746fe..443bde3 100644 --- a/orgui/app/QUBCalculator.py +++ b/orgui/app/QUBCalculator.py @@ -907,11 +907,11 @@ def readConfig(self, configfile): return False def toFallbackConfig(self): - sdd = 0.729 # m - E = 78.0 - pixelsize = 172e-6 - cp = [731.0, 1587.856] - self.mu = np.deg2rad(0.05) + sdd = 1.0 # m + E = 75.0 + pixelsize = 150e-6 + cp = [1100.0, 2020.0] + self.mu = np.deg2rad(0.08) self.chi = 0.0 self.phi = 0.0 self.n = 1 - 1.1415e-06 @@ -926,7 +926,7 @@ def toFallbackConfig(self): self.polfactor = 0 self.azimuth = 0 self.detectorCal = DetectorCalibration.Detector2D_SXRD() - self.detectorCal.detector = pyFAI.detector_factory("Pilatus2m") + self.detectorCal.detector = pyFAI.detector_factory("pilatus44mcdte") self.detectorCal.setFit2D( sdd * 1e3, cp[0], cp[1], pixelX=pixelsize * 1e6, pixelY=pixelsize * 1e6 ) @@ -2485,8 +2485,18 @@ def hideEvent(self, event): super().hideEvent(event) def resetParameters(self): + """Restore the parameters that were active when the dialog was opened. + + The individual editors apply their changes immediately, so restoring + the widgets is not enough: the saved parameters have to be applied + again, otherwise the edited configuration stays active. They are + emitted directly rather than read back from the widgets, so that the + limited precision of the spin boxes cannot alter them. + """ + if self.savedParams is None: + return self.machineparams.setValues(self.savedParams) - # self.machineparams._onAnyValueChanged() + self.machineparams.sigMachineParamsChanged.emit(self.savedParams) def onCancel(self): self.resetParameters() @@ -2540,9 +2550,14 @@ def hideEvent(self, event): super().hideEvent(event) def resetParameters(self): + """Restore the crystal that was active when the dialog was opened. + + As for the machine parameters, the editors apply their changes + immediately, so the saved crystal has to be applied again. + """ if self.savedParams is not None: self.crystalparams.setValues(*self.savedParams) - # self.crystalparams._onAnyValueChanged() + self.crystalparams.sigCrystalParamsChanged.emit(*self.savedParams) def onCancel(self): self.resetParameters() diff --git a/orgui/app/cpp/qconversion_cpp.cpp b/orgui/app/cpp/qconversion_cpp.cpp new file mode 100644 index 0000000..b30538a --- /dev/null +++ b/orgui/app/cpp/qconversion_cpp.cpp @@ -0,0 +1,106 @@ +// Reciprocal-space conversion kernel for the experimental Q-plot. +// +// The whole conversion of a detector image collapses into one affine relation, +// see orgui/app/qconversion.py. With inv = 1 / |(x, y, z)| every Cartesian +// component of the momentum transfer is +// +// q_j = k * ( (G[j] . (x, y, z)) * inv - c[j] ) +// +// so the in-plane and the out-of-plane component can be produced in a single +// pass over the pixel positions. That matters because the arrays are large: +// a Pilatus 6M has more than six million pixels. + +#include +#include +#include +#include + +#include +#include + +namespace py = pybind11; + +using DoubleArray = py::array_t; + +// Returns the signed in-plane and the out-of-plane momentum transfer. +// The sign of the in-plane component follows orGUI's delta direction. +py::tuple q_ip_oop( + const DoubleArray &x, + const DoubleArray &y, + const DoubleArray &z, + const DoubleArray &conversion_matrix, + const DoubleArray &offset, + const double k) +{ + const py::buffer_info x_info = x.request(); + const py::buffer_info y_info = y.request(); + const py::buffer_info z_info = z.request(); + if (y_info.size != x_info.size || z_info.size != x_info.size) { + throw std::invalid_argument( + "x, y and z must have the same number of elements"); + } + + const py::buffer_info matrix_info = conversion_matrix.request(); + if (matrix_info.ndim != 2 || matrix_info.shape[0] != 3 + || matrix_info.shape[1] != 3) { + throw std::invalid_argument("the conversion matrix must be 3x3"); + } + const py::buffer_info offset_info = offset.request(); + if (offset_info.size != 3) { + throw std::invalid_argument("the offset must have three elements"); + } + + const std::vector shape( + x_info.shape.begin(), x_info.shape.end()); + py::array_t q_ip(shape); + py::array_t q_oop(shape); + + const double *const xp = static_cast(x_info.ptr); + const double *const yp = static_cast(y_info.ptr); + const double *const zp = static_cast(z_info.ptr); + const double *const gp = static_cast(matrix_info.ptr); + const double *const cp = static_cast(offset_info.ptr); + double *const ip = static_cast(q_ip.request().ptr); + double *const oop = static_cast(q_oop.request().ptr); + + // Local copies keep the coefficients in registers and let the compiler + // vectorise the loop. + const double g00 = gp[0], g01 = gp[1], g02 = gp[2]; + const double g10 = gp[3], g11 = gp[4], g12 = gp[5]; + const double g20 = gp[6], g21 = gp[7], g22 = gp[8]; + const double c0 = cp[0], c1 = cp[1], c2 = cp[2]; + const py::ssize_t count = x_info.size; + + { + py::gil_scoped_release release; + for (py::ssize_t i = 0; i < count; ++i) { + const double xi = xp[i]; + const double yi = yp[i]; + const double zi = zp[i]; + const double inv = 1.0 / std::sqrt(xi * xi + yi * yi + zi * zi); + const double qx = k * ((g00 * xi + g01 * yi + g02 * zi) * inv - c0); + const double qy = k * ((g10 * xi + g11 * yi + g12 * zi) * inv - c1); + const double qz = k * ((g20 * xi + g21 * yi + g22 * zi) * inv - c2); + const double radial = std::sqrt(qx * qx + qy * qy); + ip[i] = qx >= 0.0 ? radial : -radial; + oop[i] = qz; + } + } + + return py::make_tuple(q_ip, q_oop); +} + +PYBIND11_MODULE(_qconversion_cpp, m) +{ + m.doc() = "Reciprocal-space conversion kernel for the experimental Q-plot."; + m.def( + "q_ip_oop", + &q_ip_oop, + py::arg("x"), + py::arg("y"), + py::arg("z"), + py::arg("conversion_matrix"), + py::arg("offset"), + py::arg("k"), + "Return the signed in-plane and the out-of-plane momentum transfer."); +} diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index 1497146..dd4b87d 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -36,7 +36,7 @@ import copy import sys import os -import re +from contextlib import contextmanager from silx.gui import qt from io import StringIO @@ -56,6 +56,7 @@ from silx.gui.plot.Profile import ProfileToolBar from silx.gui.plot.tools.roi import RegionOfInterestManager from silx.gui.plot.actions import control as control_actions +from silx.gui.widgets.ElidedLabel import ElidedLabel import traceback @@ -78,6 +79,7 @@ import numpy as np from ..datautils.xrayutils import HKLVlieg, CTRcalc +from . import qconversion from ..datautils.xrayutils import ReciprocalNavigation as rn # legacy import: @@ -91,6 +93,10 @@ except Exception: console = False +# Reciprocal-space (Q) conversion of detector images relies on pyFAI's +# FiberIntegrator, which was introduced in pyFAI 2025.1. +HAS_FIBER_INTEGRATOR = qconversion.HAS_FIBER_INTEGRATOR + try: from . import _roi_sum_accel @@ -175,6 +181,30 @@ def __init__(self, configfile, parent=None): toolbar.addAction( control_actions.OpenGLAction(parent=toolbar, plot=self.centralPlot) ) + self.plotAgainstQAct = qt.QAction(resources.getQicon("q-plot"), "Q-plot", self) + self.plotAgainstQAct.setCheckable(True) + self.plotAgainstQAct.setToolTip( + "Experimental: show the image in reciprocal space coordinates. " + "The intensities are rebinned onto a regular grid." + ) + self.plotAgainstQAct.toggled.connect(self._convertImagetoQ) + toolbar.addAction(self.plotAgainstQAct) + + self.qFrameSelector = qt.QComboBox() + for frame in qconversion.FRAMES: + self.qFrameSelector.addItem(qconversion.FRAME_LABELS[frame], frame) + self.qFrameSelector.setCurrentIndex( + qconversion.FRAMES.index(qconversion.DEFAULT_FRAME) + ) + self.qFrameSelector.setToolTip( + "Experimental: reciprocal space frame used by the Q-plot" + ) + self.qFrameSelector.currentIndexChanged.connect(self._onQFrameChanged) + toolbar.addWidget(self.qFrameSelector) + self._qPlotExperimentalWarned = False + self._qPlotPreviousYInverted = None + self._qPlotPreviousXInverted = None + self._qPlotRefreshing = False self.centralPlot.addToolBar(toolbar) self.maskManager = MaskManager() self.maskConfigDialog = MaskConfigDialog(self.maskManager, self) @@ -256,6 +286,7 @@ def __init__(self, configfile, parent=None): self.roPkIntegrTab = RockingPeakIntegrator(self.database) + self.scanSelector.sigRefreshH5.connect(self._onRefreshH5) self.scanSelector.sigImageNoChanged.connect(self._onSliderValueChanged) self.scanSelector.sigImagePathChanged.connect(self._onImagePathChanged) @@ -336,6 +367,9 @@ def __init__(self, configfile, parent=None): self.ubcalc.sigPlottableMachineParamsChanged.connect(self._onPlotMachineParams) self.ubcalc.sigReplotRequest.connect(self.updatePlotItems) + # the Q-plot coordinates depend on the machine angles, the azimuthal + # reference, the energy and, for the crystal frame, on U + self.ubcalc.sigReplotRequest.connect(self._refreshQPlot) self.allimgsum = None self.allimgmax = None self.reflectionSel.setSizePolicy( @@ -595,6 +629,45 @@ def _repair_config_for_image(self, image_shape): col_gaps = np.empty((0, 2), dtype=np.int32) return use_repair, repair, row_gaps, col_gaps + def _onRefreshH5(self): + """GUI/CLI hint: re-read the open HDF5 file and reload the current scan. + + When the file is still being written by the acquisition system, the + tree view has to be rebuilt and the scan reloaded, otherwise the nodes + keep pointing at stale links. For console use the equivalent manual + recipe is:: + + self.fscan = None + for i in self.centralPlot.getAllImages(): + if 'scan' in i.getLegend(): + self.centralPlot.removeImage(i) + self.currentAddImageLabel = None + """ + self.scanSelector._onDoRefresh() # update the treeview + + # Reloading the scan must not regenerate the max/sum images, so the + # auto load and beamtime autodetect options are disabled temporarily. + autoload = self.autoLoadAct.isChecked() + bt_autodetect = self.scanSelector.bt_autodetect_enable.isChecked() + if autoload: + self.autoLoadAct.setChecked(False) + if bt_autodetect: + self.scanSelector.bt_autodetect_enable.setChecked(False) + + # very important! this seems to unlock the file somehow + self.scanSelector._onLoadScan() + + # restore auto load and bt autodetect + if autoload: + self.autoLoadAct.setChecked(True) + if bt_autodetect: + self.scanSelector.bt_autodetect_enable.setChecked(True) + + self.images_loaded = True + # expanding here is safe, in contrast to during the tree rebuild + self.scanSelector.hdfTreeView.expandToDepth(0) + # todo: restore all expanded branches as they were before + def _removeAllIntegrPlotCurves(self): """CLI-capable: clear all integration curves from the plot widget.""" # remove plotted curves @@ -3063,30 +3136,158 @@ def newXyHKLConverter(self): """ def xyToHKL(x, y): - """CLI-safe: convert detector pixels to hkl and detector angles.""" - # print("xytoHKL:") - # print("x,y = %s, %s" % (x,y)) + """CLI-safe: convert plot coordinates to hkl, angles and Q. + + The plot coordinates are detector pixels, or the in-plane and + out-of-plane momentum transfer while the Q-plot is switched on. In + the latter case they are converted back to detector angles first, so + that both modes report the same hkl and the same detector angles for + the same point of the sample. + + The momentum transfer is the one belonging to the displayed + coordinates: it is refraction corrected in pixel mode, like + everywhere else in orGUI, and unrefracted in the Q-plot, where the + axes are unrefracted as well. Reporting the corrected value there + would contradict the position the cursor is at. + + :returns: + ``[h, k, l, delta, gamma, Qx, Qy, Qz]``, with the angles in + degrees and the momentum transfer of the currently selected + reciprocal-space frame in inverse Angstrom. All-NaN when the + position does not correspond to a point on the detector. + """ if self.fscan is None: - return np.array([np.nan, np.nan, np.nan, np.nan, np.nan]) - mu, om = self.getMuOm(self.imageno) - gamma, delta = self.ubcalc.detectorCal.surfaceAnglesPoint( - np.array([y]), np.array([x]), mu - ) - # print(self.ubcalc.detectorCal) - # print(x,y) - # print(self.ubcalc.detectorCal.tth(np.array([y]),np.array([x]))) - pos = [mu, delta[0], gamma[0], om, self.ubcalc.chi, self.ubcalc.phi] + return np.full(8, np.nan) + if self.plotAgainstQAct.isChecked(): + resolved = self._qPlotToPosition(x, y) + if resolved is None: + return np.full(8, np.nan) + mu, om, delta, gamma, Q = resolved + else: + mu, om = self.getMuOm(self.imageno) + gamma_arr, delta_arr = self.ubcalc.detectorCal.surfaceAnglesPoint( + np.array([y]), np.array([x]), mu + ) + delta, gamma = delta_arr[0], gamma_arr[0] + Q = None + pos = [mu, delta, gamma, om, self.ubcalc.chi, self.ubcalc.phi] pos = HKLVlieg.crystalAngles(pos, self.ubcalc.n) - hkl = np.concatenate( + return np.concatenate( ( np.array(self.ubcalc.angles.anglesToHkl(*pos)), - np.rad2deg([delta[0], gamma[0]]), + np.rad2deg([delta, gamma]), + self._qInSelectedFrame(pos) if Q is None else Q, ) ) - return hkl return xyToHKL + def _qPlotOrientationMatrix(self): + """Orientation matrix for the Q-plot, or ``None`` when unavailable.""" + try: + return self.ubcalc.ubCal.getU() + except Exception: + return None + + def _qPlotToPosition(self, q_ip, q_oop): + """Momentum transfer and detector angles of one point of the Q-plot. + + Inverts the projection of :func:`~orgui.app.qconversion.qIpOop` with + :func:`~orgui.app.qconversion.qFromIpOop`. The alpha frame always has a + unique solution; the rotated frames can have two, and the one that falls + onto the detector is kept. That test doubles as the check whether the + position is covered by the image at all: the rebinned grid is + rectangular and reaches beyond the detector. + + The returned momentum transfer is the exact reconstruction, so its + in-plane and out-of-plane components are the requested coordinates. No + refraction correction is applied, because the Q-plot axes carry none + either. + + :returns: + ``(mu, omega, delta, gamma, Q)`` with the angles in rad and ``Q`` + in the selected frame in inverse Angstrom, or ``None`` when the + position is not on the detector. + :rtype: tuple or None + """ + frame = self._selectedQFrame() + # the angles the displayed Q-plot was built with, so that the readout + # inverts the conversion that actually produced the image + alpha, omega = self._qPlotAngles() + try: + solutions = qconversion.qFromIpOop( + frame, + q_ip, + q_oop, + self.ubcalc.ubCal.getK(), + alpha=alpha, + omega=omega, + chi=self.ubcalc.chi, + phi=self.ubcalc.phi, + U=self._qPlotOrientationMatrix(), + ) + except Exception: + return None + + for Q in solutions: + try: + angles = qconversion.detectorAngles( + Q, + self.ubcalc.ubCal.getK(), + frame=frame, + alpha=alpha, + omega=omega, + chi=self.ubcalc.chi, + phi=self.ubcalc.phi, + U=self._qPlotOrientationMatrix(), + ) + except Exception: + continue + if angles is None: + continue + delta, gamma = angles + if self._isOnDetector(delta, gamma, alpha): + return alpha, omega, delta, gamma, Q + return None + + def _isOnDetector(self, delta, gamma, alpha): + """Whether the surface-frame angles hit the active detector.""" + dc = self.ubcalc.detectorCal + try: + yx = dc.pixelsSurfaceAngles( + np.array([gamma]), np.array([delta]), np.array([alpha]) + ) + except Exception: + return False + yx = np.asarray(yx).reshape(-1, 2)[0] + if not np.all(np.isfinite(yx)): + return False + shape = dc.detector.shape + return 0 <= yx[0] < shape[0] and 0 <= yx[1] < shape[1] + + def _qInSelectedFrame(self, pos): + """Momentum transfer of one pixel in the selected reciprocal frame. + + :param pos: + Refracted ``[alpha, delta, gamma, omega, chi, phi]`` angles in rad. + :returns: + ``(Qx, Qy, Qz)`` in inverse Angstrom, or NaN when the frame cannot + be evaluated, for example ``Q_cryst`` without an orientation matrix. + :rtype: numpy.ndarray + """ + try: + rotation = qconversion.frameMatrix( + self._selectedQFrame(), + alpha=pos[0], + omega=pos[3], + chi=pos[4], + phi=pos[5], + U=self.ubcalc.ubCal.getU(), + ) + return rotation @ self.ubcalc.angles.QAlpha(pos[0], pos[1], pos[2]) + except Exception: + return np.full(3, np.nan) + def getMuOm(self, imageno=None): """Return mu and omega for an image or the whole active scan. @@ -3285,6 +3486,84 @@ def _onCreateScan(self): traceback.format_exc(), ) + @staticmethod + def _scanEntriesFromBackend(backendcls, h5file): + """Ask the backend to list the scans of an opened file. + + This is the fast path of the segmented scan loader. The backend may + return bare identifiers or ``(identifier, label)`` pairs, see + :meth:`~orgui.backend.scans.Scan.listScans`. + + The identifier is whatever the backend is opened with. It is reported + as both ``scanno`` and ``name`` because + :func:`~orgui.backend.backends.openScan` reads one or the other + depending on the beamtime, and only the backend knows which of the two + its constructor wants. + + :returns: + List of ``(ddict, label)``, where ``ddict`` holds the keys + identifying the scan and ``label`` may be ``None``. ``None`` when + the backend does not implement ``listScans``. + :rtype: list or None + + .. note:: + CLI-safe. + """ + try: + listed = backendcls.listScans(h5file) + except NotImplementedError: + return None + + entries = [] + for entry in listed: + if isinstance(entry, tuple): + identifier, label = entry + else: + identifier, label = entry, None + entries.append(({"scanno": identifier, "name": str(identifier)}, label)) + return entries + + @staticmethod + def _scanEntriesFromNodes(backendcls, model, rootIndex): + """List the scans of an opened file through ``parse_h5_node``. + + Fallback for backends that do not implement ``listScans``. Every entry + of the file root is wrapped in the same :class:`silx.gui.hdf5.H5Node` + the tree hands to the backend when a single scan is opened, and its own + ``parse_h5_node`` result is kept unchanged. The segmented loader and + the normal loader therefore cannot disagree about what a scan is called + or how it is numbered. + + An entry counts as a scan when ``parse_h5_node`` returns a ``scanno`` + for it. Entries it raises on, or reports without a ``scanno``, are not + scans and are skipped. Scan numbers seen more than once are reported + once, which is what collapses the subscans of one measurement. + + :returns: list of ``(ddict, label)`` with ``label`` always ``None``. + :rtype: list + + .. note:: + GUI-only, it needs the HDF5 tree model. + """ + entries = [] + seen = set() + for row in range(model.rowCount(rootIndex)): + index = model.index(row, 0, rootIndex) + item = model.data(index, role=silx.gui.hdf5.Hdf5TreeModel.H5PY_ITEM_ROLE) + if item is None: + continue + try: + ddict = backendcls.parse_h5_node(silx.gui.hdf5.H5Node(item)) + except Exception: + continue + if not isinstance(ddict, dict) or ddict.get("scanno") is None: + continue + if ddict["scanno"] in seen: + continue + seen.add(ddict["scanno"]) + entries.append((dict(ddict), None)) + return entries + def _onLoadInterlacedScan(self): """GUI-only: build an interlaced scan from selected HDF5 scans.""" # GUI function to concatenate multiple scans into one @@ -3307,61 +3586,37 @@ def _onLoadInterlacedScan(self): # address the root node to correctly get the scan names if rootI.parent().isValid(): - h5file = model.data( - model.parent(rootI), role=silx.gui.hdf5.Hdf5TreeModel.H5PY_OBJECT_ROLE + rootI = model.parent(rootI) + h5file = model.data(rootI, role=silx.gui.hdf5.Hdf5TreeModel.H5PY_OBJECT_ROLE) + + # Which entries of the file are scans, and how they are addressed, is a + # property of the file format, so the backend is asked instead of the + # naming convention being guessed here. Backends that implement + # listScans answer directly; for the others every entry is run through + # their own parse_h5_node. + backendcls = backends.fscans[self.scanSelector.btid.currentText()] + try: + entries = self._scanEntriesFromBackend(backendcls, h5file) + if entries is None: + entries = self._scanEntriesFromNodes(backendcls, model, rootI) + except Exception: + qutils.warning_detailed_message( + self, + "Can not create interlaced scan", + f"Can not list the scans of this file with the " + f"{backendcls.__name__} backend.", + traceback.format_exc(), ) - else: - h5file = model.data( - rootI, role=silx.gui.hdf5.Hdf5TreeModel.H5PY_OBJECT_ROLE + return + if not entries: + qutils.warning_detailed_message( + self, + "Can not create interlaced scan", + f"The {backendcls.__name__} backend did not find any scan in " + "this file.", + "", ) - - isID31 = self.scanSelector.btid.currentText() in [ - "ch5523", - "ch5700", - "ch5918", - "ch6392", - "ch7131", - "ch7149", - "ch7856", - "ch8153", - "id31_default", - ] - kl_full = list(h5file.keys()) - kl = np.empty(0, dtype=int) - for i in kl_full: - if isID31: - pattern = r"\.\d" - result = re.findall(pattern, i)[0][1:] - if result == "1": - # select only scan names which are ending on suffix '.1' (fast counters of id31 hdf5 format) # noqa: E501 - kl = np.append(kl, i) - else: - kl = np.append(kl, i) - - # separate scan nr and delete duplicates suffixes - if isID31: - # try to get the scan nr and '/title' from the hdf5 file - nr = np.empty(0, dtype=int) - name = np.empty(0, dtype=str) - - for i in kl: - pattern = r"\d+\." - result = re.findall(pattern, i) - name = np.append(name, h5file[i + "/title"]) - nr = np.append(nr, int(result[0][:-1])) - - lsort = np.argsort(nr)[::1] - nr = nr[lsort] - name = name[lsort] - else: - nr = np.empty(0, dtype=int) - name = np.empty(0, dtype=str) - for nth, i in enumerate(kl): - name = np.append(name, i) - nr = np.append( - nr, nth + 1 - ) # create scan nr list with ascending integers, starting with 1 - # This will later be used to address the subscans, so check if your scans are handled like this!!! # noqa: E501 + return # open GUI dialog to select which scans to combine interlacedSelectDialog = qt.QDialog() @@ -3373,17 +3628,13 @@ def _onLoadInterlacedScan(self): b = qt.QFormLayout() box = qt.QGroupBox() scanBoxes = [] - # labels = [] - for i, item in enumerate(nr): + for ddict, label in entries: ithScanBox = qt.QCheckBox() scanBoxes.append(ithScanBox) - # labels.append(qt.QLabel('Scan '+str(item)+':'+str(name[i][2:-1]))) - if isID31: - b.addRow( - qt.QLabel("Scan " + str(item) + ": " + name[i].decode()), ithScanBox - ) - else: - b.addRow(qt.QLabel("Scan " + str(item) + ": " + name[i]), ithScanBox) + # the label is the scan title where the backend provides one, and + # the scan name otherwise + description = label if label is not None else ddict.get("name", "") + b.addRow(qt.QLabel(f"Scan {ddict['scanno']}: {description}"), ithScanBox) box.setLayout(b) a.setWidget(box) @@ -3427,21 +3678,19 @@ def _onLoadInterlacedScan(self): return # generate scan objects for selected scans - selectedScans = [] - for i, j in enumerate(scanBoxes): - if j.isChecked(): - selectedScans.append(nr[i]) - # selectedScans.append(name[i]) + selectedScans = [ + ddict for (ddict, _label), box in zip(entries, scanBoxes) if box.isChecked() + ] nodes = list(self.scanSelector.hdfTreeView.selectedH5Nodes()) obj = nodes[0] scansegments = [] - for i in selectedScans: - ddict = dict() - ddict["scanno"] = int(i) + for selected in selectedScans: + # the keys identifying the scan come from the backend, so that a + # backend addressed by name rather than by number also works + ddict = dict(selected) ddict["file"] = obj.local_filename - # ddict['node'] = kl[i] ddict["beamtime"] = IS_btid.currentText() try: scansegments.append(backends.openScan(IS_btid.currentText(), ddict)) @@ -3470,7 +3719,7 @@ def _onLoadInterlacedScan(self): self.scanSelector.setAxis(self.fscan.axis, self.fscan.axisname) self.activescanname = "{}-segmentedScan {} {}-{}".format( self.fscan.axisname, - ",".join(str(itemNr) for itemNr in selectedScans), + ",".join(str(ddict["scanno"]) for ddict in selectedScans), np.amin(self.fscan.axis), np.amax(self.fscan.axis), ) @@ -3823,12 +4072,24 @@ def _onChangeImage(self, imageno): ) return self.scanSelector.slider.setValue(imageno) - self.plotImage(self.scanSelector.slider.value()) + with self._suspendedQPlot(): + if self.scanSelector.showMaxAct.isChecked(): + self.scanSelector.showMaxAct.setChecked(False) + if self.scanSelector.showSumAct.isChecked(): + self.scanSelector.showSumAct.setChecked(False) + self.plotImage(self.scanSelector.slider.value()) + self._refreshQPlot() def _onSliderValueChanged(self, value): """GUI/CLI hint: replot the image selected by the scan slider.""" if self.fscan is not None: - self.plotImage(value) + with self._suspendedQPlot(): + if self.scanSelector.showMaxAct.isChecked(): + self.scanSelector.showMaxAct.setChecked(False) + if self.scanSelector.showSumAct.isChecked(): + self.scanSelector.showSumAct.setChecked(False) + self.plotImage(value) + self._refreshQPlot() # print(self.centralPlot._callback) def _onLoadAll(self): @@ -3958,7 +4219,9 @@ def _onMaxToggled(self, value): z=1, ) self.centralPlot.setActiveImage(self.currentAddImageLabel) - self.scanSelector.alphaslider.setLegend(self.currentAddImageLabel) + self.scanSelector.alphaslider.setLegend("special") + if not self.scanSelector.showMaxAct.isChecked(): + self.scanSelector.showMaxAct.setChecked(True) else: self.scanSelector.showMaxAct.setChecked(False) else: @@ -3966,6 +4229,8 @@ def _onMaxToggled(self, value): self.centralPlot.setActiveImage(self.currentImageLabel) self.centralPlot.removeImage(self.currentAddImageLabel) self.currentAddImageLabel = None + # the displayed data changed, so the Q-plot has to follow + self._refreshQPlot() def _onSumToggled(self, value): """GUI/CLI hint: toggle display of the precomputed summed image.""" @@ -3999,7 +4264,9 @@ def _onSumToggled(self, value): z=1, ) self.centralPlot.setActiveImage(self.currentAddImageLabel) - self.scanSelector.alphaslider.setLegend(self.currentAddImageLabel) + self.scanSelector.alphaslider.setLegend("special") + if not self.scanSelector.showSumAct.isChecked(): + self.scanSelector.showSumAct.setChecked(True) else: self.scanSelector.showSumAct.setChecked(False) else: @@ -4007,6 +4274,8 @@ def _onSumToggled(self, value): self.centralPlot.setActiveImage(self.currentImageLabel) self.centralPlot.removeImage(self.currentAddImageLabel) self.currentAddImageLabel = None + # the displayed data changed, so the Q-plot has to follow + self._refreshQPlot() def _roi_preview_enabled(self): """Return whether fitted-background ROI preview should alter the image.""" @@ -4022,6 +4291,312 @@ def _set_center_roi_preview_color(self, roi): roi.setColor("blue" if self._roi_preview_enabled() else "red") roi.setBgStyle("pink", "-", roi.getLineWidth()) + def _selectedQFrame(self): + """Return the reciprocal space frame selected for the Q-plot.""" + frame = self.qFrameSelector.currentData() + if frame not in qconversion.FRAMES: + return qconversion.DEFAULT_FRAME + return frame + + def _onQFrameChanged(self, index): + """GUI hint: follow the frame selection in the readout and the Q-plot.""" + del index + self.centralPlot.setQFrame(self._selectedQFrame()) + self._refreshQPlot() + + def _refreshQPlot(self, *args): + """GUI hint: rebuild the Q-plot so that it follows the current state. + + Called whenever something the conversion depends on changed: the + displayed image, the machine angles, the azimuthal reference, the + energy or the orientation matrix. The conversion replaces the previous + reciprocal-space image, so the Q-plot stays switched on. Does nothing + while the Q-plot is off, or while a caller batches several changes by + holding ``_qPlotRefreshing``. + """ + del args + if not self.plotAgainstQAct.isChecked() or self._qPlotRefreshing: + return + self._qPlotRefreshing = True + try: + self._convertImagetoQ(True) + finally: + self._qPlotRefreshing = False + + @contextmanager + def _suspendedQPlot(self): + """Batch several changes into a single Q-plot rebuild.""" + previous = self._qPlotRefreshing + self._qPlotRefreshing = True + try: + yield + finally: + self._qPlotRefreshing = previous + + def _qPlotAngles(self): + """Return the alpha and omega angles used for the current Q-plot.""" + try: + return self.getMuOm(self.scanSelector.slider.value()) + except Exception: + return self.ubcalc.mu, 0.0 + + def _convertImagetoQ(self, value): + """GUI/CLI hint: toggle display of the image in reciprocal space. + + **Experimental.** The currently displayed image is rebinned onto a + regular grid of in-plane and out-of-plane momentum transfer by + :func:`~orgui.app.qconversion.integrateImage`, which + drives pyFAI's ``FiberIntegrator`` with orGUI's own angle conventions. + The frame is taken from the frame selector next to the action and + defaults to the alpha (surface) frame, which is what + :meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.QAlpha` returns. + See the geometry section of the documentation. + """ + if value: + # check if conversion possible + if not HAS_FIBER_INTEGRATOR: + logger.error( + "Q conversion requires pyFAI >= 2025.1", + extra={ + "title": "Cannot convert image to Q", + "description": "Conversion of the image to Q requires pyFAI >= 2025.1. Please update pyFAI.", # noqa: E501 + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + if self.fscan is None: + logger.error( + "Q conversion failed: no scan loaded", + extra={ + "title": "Cannot convert image to Q", + "description": "No scan loaded!", + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + if self.fscan.axisname == "mu": + logger.error( + "Q conversion failed: mu scan", + extra={ + "title": "Cannot convert image to Q", + "description": "Conversion not implemented for mu scans!", + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + + # select the data to convert before touching the plot, so that a + # failure leaves the displayed images untouched + showmax = self.scanSelector.showMaxAct.isChecked() + showsum = self.scanSelector.showSumAct.isChecked() + if showmax: + data = self.allimgmax + elif showsum: + data = self.allimgsum + else: + data = None + for i in self.centralPlot.getAllImages(): + if i.getLegend() == "scan_image": + data = i.getData() + break + if data is None: + logger.error( + "Q conversion failed: no image available", + extra={ + "title": "Cannot convert image to Q", + "description": "No image available for the conversion!", + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + + frame = self._selectedQFrame() + if (showmax or showsum) and frame in qconversion.FRAMES_REQUIRING_OMEGA: + logger.error( + "Q conversion failed: frame not defined for max/sum images", + extra={ + "title": "Cannot convert image to Q", + "description": f"The maximum and sum images combine many omega angles, so the {qconversion.FRAME_LABELS[frame]} frame is not defined for them. Use the {qconversion.FRAME_LABELS['Q_alpha']} or {qconversion.FRAME_LABELS['Q_lab']} frame instead.", # noqa: E501 + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + + if not self._qPlotExperimentalWarned: + logger.warning( + "The Q-plot is experimental and its conventions may still change." + ) + self._qPlotExperimentalWarned = True + + # perform conversion into Q coordinates before touching the plot, + # so that a failure leaves the display as it was + alpha, omega = self._qPlotAngles() + try: + res2d = qconversion.integrateImage( + self.ubcalc.detectorCal, + data, + alpha, + frame=frame, + omega=omega, + chi=self.ubcalc.chi, + phi=self.ubcalc.phi, + U=self.ubcalc.ubCal.getU(), + ) + # An inverted geometry, such as the downward scattering setup + # at ID31, has the out-of-plane coordinate growing towards the + # bottom of the detector image. Drawing it upwards anyway would + # show the reciprocal-space image mirrored with respect to the + # detector image, so the ordinate follows the geometry. + flipOrdinate = qconversion.outOfPlaneIncreasesWithRow( + self.ubcalc.detectorCal, + alpha, + frame=frame, + omega=omega, + chi=self.ubcalc.chi, + phi=self.ubcalc.phi, + U=self._qPlotOrientationMatrix(), + ) + except Exception: + logger.error( + "Q conversion failed", + extra={ + "title": "Cannot convert image to Q", + "description": f"The conversion into the {qconversion.FRAME_LABELS[frame]} frame failed:\n{traceback.format_exc()}", # noqa: E501 + "show_dialog": True, + "dialog_level": logging.ERROR, + "parent": self, + }, + ) + self.plotAgainstQAct.setChecked(False) + return + + # hide scan image, remove max/sum image + for i in self.centralPlot.getAllImages(): + if i.getLegend() == "scan_image": + i.setVisible(False) + else: + self.centralPlot.removeImage(i) + + # Reciprocal space is plotted with the ordinate along the + # out-of-plane coordinate, unlike the image convention used for the + # pixel coordinates. The previous state is only recorded when the + # Q-plot is entered, so that rebuilding it in place does not + # overwrite it. + if self._qPlotPreviousYInverted is None: + self._qPlotPreviousYInverted = self.centralPlot.isYAxisInverted() + if self._qPlotPreviousXInverted is None: + self._qPlotPreviousXInverted = self.centralPlot.isXAxisInverted() + self.centralPlot.setYAxisInverted(flipOrdinate) + self.centralPlot.setXAxisInverted(False) + + # plot generated image, q_par on the abscissa and q_perp on the + # ordinate for every frame + oopmin, oopmax = np.min(res2d.outofplane), np.max(res2d.outofplane) + ipmin, ipmax = np.min(res2d.inplane), np.max(res2d.inplane) + scale_ip = (ipmax - ipmin) / res2d.inplane.size + scale_oop = (oopmax - oopmin) / res2d.outofplane.size + self.currentAddImageLabel = self.centralPlot.addImage( + res2d.intensity, + legend="qImage", + replace=False, + resetzoom=False, + copy=True, + z=2, + xlabel=r"q$_\parallel / \, \AA^{-1}$", + ylabel=r"q$_\perp / \, \AA^{-1}$", + scale=(scale_ip, scale_oop), + origin=(res2d.inplane[0], res2d.outofplane[0]), + ) + self.centralPlot.getXAxis().setLimits(ipmin, ipmax) + self.centralPlot.getYAxis().setLimits(oopmin, oopmax) + + # apply active plot settings; the alpha slider expects a legend, + # not the image item that addImage returns + self.centralPlot.setActiveImage(self.currentAddImageLabel) + self.scanSelector.alphaslider.setLegend("qImage") + + # the axes now hold momentum transfer, so the readout has to invert + # the conversion instead of reading detector pixels + self.centralPlot.setQPlotActive(True) + + else: + self.centralPlot.setQPlotActive(False) + + # Restore the axis orientation of the pixel coordinates. This only + # depends on what was recorded when the Q-plot was entered, so it + # must not be tied to the max/sum image below: a conversion that + # failed after the axes were flipped would otherwise leave them + # flipped for good. + restored = False + if self._qPlotPreviousYInverted is not None: + self.centralPlot.setYAxisInverted(self._qPlotPreviousYInverted) + self._qPlotPreviousYInverted = None + restored = True + if self._qPlotPreviousXInverted is not None: + self.centralPlot.setXAxisInverted(self._qPlotPreviousXInverted) + self._qPlotPreviousXInverted = None + restored = True + + # show scan image, delete qImage + hadQImage = False + for i in self.centralPlot.getAllImages(): + if i.getLegend() == "scan_image": + i.setVisible(True) + elif i.getLegend() == "qImage": + self.centralPlot.removeImage(i) + hadQImage = True + + if hadQImage or self.currentAddImageLabel is not None: + # plot max/sum, set correct active image + if ( + self.scanSelector.showMaxAct.isChecked() + and self.allimgmax is not None + ): + self.currentAddImageLabel = self.centralPlot.addImage( + self.allimgmax, + legend="special", + replace=False, + resetzoom=False, + copy=True, + z=1, + ) + self.centralPlot.setActiveImage(self.currentAddImageLabel) + elif ( + self.scanSelector.showSumAct.isChecked() + and self.allimgsum is not None + ): + self.currentAddImageLabel = self.centralPlot.addImage( + self.allimgsum, + legend="special", + replace=False, + resetzoom=False, + copy=True, + z=1, + ) + self.centralPlot.setActiveImage(self.currentAddImageLabel) + else: + self.currentAddImageLabel = None + self.centralPlot.setActiveImage(self.currentImageLabel) + + # apply correct zoom for pixel coordinates + if restored or hadQImage: + self.centralPlot.resetZoom() + def _roi_array_from_key(self, key): """Convert one ROI slice key to the compiled ``(1, 2, 2)`` format.""" roi = np.array( @@ -4880,7 +5455,8 @@ def integrateROI(self): ) else: [ - l.append(np.array([[0, 0], [0, 0]])) for l in roi_lists[1:] # noqa: E741 + l.append(np.array([[0, 0], [0, 0]])) + for l in roi_lists[1:] # noqa: E741 ] # will result in zeros, convert to np.nan later roi_lists[0].append(np.array([[0, 0], [0, 0]])) if hkl_del_gam_2[i, -1]: @@ -4896,11 +5472,13 @@ def integrateROI(self): ) else: [ - l.append(np.array([[0, 0], [0, 0]])) for l in roi_lists[1:] # noqa: E741 + l.append(np.array([[0, 0], [0, 0]])) + for l in roi_lists[1:] # noqa: E741 ] # will result in zeros, convert to np.nan later roi_lists[0].append(np.array([[0, 0], [0, 0]])) roi_lists = [ - np.ascontiguousarray(np.stack(l), dtype=np.int64) for l in roi_lists # noqa: E741 + np.ascontiguousarray(np.stack(l), dtype=np.int64) + for l in roi_lists # noqa: E741 ] roi_lists_accel.append(roi_lists) @@ -5845,19 +6423,108 @@ def closeEvent(self, event): super().closeEvent(event) +def formatPositionTriplet(values): + """Format three numbers as ``[a b c]`` for the plot position readout. + + :param values: three numbers. + :returns: the formatted triplet, or a placeholder if any value is not + finite, which is the case while no scan is loaded. + :rtype: str + """ + values = np.asarray(values, dtype=np.float64) + if not np.all(np.isfinite(values)): + return "------" + return "[{:.5f} {:.5f} {:.5f}]".format(*values) + + +def qFieldName(frame): + """Name of the reciprocal-space field in the plot position readout. + + :param str frame: one of :data:`qconversion.FRAMES`. + :returns: for example ``"Q[alpha]"``. + :rtype: str + """ + return f"Q[{qconversion.frameShortName(frame)}]" + + +def isPositionTripletField(name): + """Whether a position readout field holds a bracketed triplet.""" + return name == "HKL" or name.startswith("Q[") + + +def axisFieldNames(frame, qPlotActive): + """Names of the two plot axis fields in the position readout. + + While the Q-plot is on, the abscissa and the ordinate are the in-plane and + the out-of-plane momentum transfer of the selected frame, not pixels. + + :param str frame: one of :data:`qconversion.FRAMES`. + :param bool qPlotActive: whether the Q-plot is switched on. + :returns: ``(abscissa, ordinate)`` field names. + :rtype: tuple[str, str] + """ + if not qPlotActive: + return ("X", "Y") + short = qconversion.frameShortName(frame) + return (f"q_par[{short}]", f"q_perp[{short}]") + + +# Widest triplet the readout should be able to show without eliding. Five +# decimals are the minimum useful precision for hkl and for Q. +POSITION_TRIPLET_TEMPLATE = "[-99.99999 -99.99999 -99.99999]" + +# The remaining fields hold a pixel coordinate, a detector angle or a data +# value such as "10, Masked". silx reserves about fourteen hash characters for +# each of them, which is wider than any of those values and takes space the +# triplets need. The width of a field is its share of the summed size hints, so +# this is a balance rather than a maximum: made much smaller, the scalar fields +# starve when the window is narrow; made much larger, the triplets are elided +# again. +POSITION_SCALAR_TEMPLATE = "-99999.99999" + + +class _PositionLabel(ElidedLabel): + """Position readout label whose size hint matches its content. + + silx sizes every readout field alike, for about fourteen hash characters. + That elides the last component of a triplet while leaving the scalar fields + wider than they need to be. Only the size *hint* is set here, not the + minimum width, so the window can still be made as narrow as before. The + hint never falls below the current text, and anything elided remains + available as a tooltip, which ``ElidedLabel`` provides by default. + """ + + def __init__(self, template, parent=None): + super().__init__(parent) + self._template = template + + def sizeHint(self): + hint = super().sizeHint() + width = self.fontMetrics().boundingRect(self._template).width() + return qt.QSize(max(hint.width(), width), hint.height()) + + class Plot2DHKL(silx.gui.plot.PlotWindow): sigKeyPressDelete = qt.pyqtSignal() def __init__(self, xyHKLConverter, parent=None, backend=None): """GUI-only: initialize the image plot with hkl position readout.""" self.xyHKLConverter = xyHKLConverter + self._qFrame = qconversion.DEFAULT_FRAME + self._qPlotActive = False + self._qFieldName = qFieldName(self._qFrame) + self._axisFieldNames = axisFieldNames(self._qFrame, self._qPlotActive) posInfo = [ - ("X", lambda x, y: x), - ("Y", lambda x, y: y), - ("H", lambda x, y: self.xyHKLConverter(x, y)[0]), - ("K", lambda x, y: self.xyHKLConverter(x, y)[1]), - ("L", lambda x, y: self.xyHKLConverter(x, y)[2]), + # pixel coordinates need no more than sub-pixel precision, and the + # space is needed for the hkl and Q triplets + (self._axisFieldNames[0], lambda x, y: self._formatAxisValue(x)), + (self._axisFieldNames[1], lambda x, y: self._formatAxisValue(y)), + ("HKL", lambda x, y: formatPositionTriplet(self.xyHKLConverter(x, y)[:3])), + ( + self._qFieldName, + lambda x, y: formatPositionTriplet(self.xyHKLConverter(x, y)[5:8]), + ), ("del", lambda x, y: self.xyHKLConverter(x, y)[3]), ("gam", lambda x, y: self.xyHKLConverter(x, y)[4]), ("Data", WeakMethodProxy(self._getImageValue)), @@ -5883,6 +6550,8 @@ def __init__(self, xyHKLConverter, parent=None, backend=None): mask=True, ) + self._tunePositionFieldWidths() + if parent is None: self.setWindowTitle("Plot2D") self.getXAxis().setLabel("Columns") @@ -5912,6 +6581,109 @@ def keyPressEvent(self, event): self.sigKeyPressDelete.emit() super().keyPressEvent(event) + def _tunePositionFieldWidths(self): + """Size every readout field for the content it actually shows. + + The triplet fields get room for all three components, and the scalar + fields give back the width silx reserves for them but never uses. + + .. note:: + GUI-only. + """ + positionInfo = self.getPositionInfoWidget() + # silx keeps the content widgets in a private list; leave the readout + # untouched if that ever changes + fields = getattr(positionInfo, "_fields", None) + if not fields: + return + layout = positionInfo.layout() + for index, (widget, name, converter) in enumerate(list(fields)): + template = ( + POSITION_TRIPLET_TEMPLATE + if isPositionTripletField(name) + else POSITION_SCALAR_TEMPLATE + ) + replacement = _PositionLabel(template, positionInfo) + replacement.setTextInteractionFlags(qt.Qt.TextSelectableByMouse) + replacement.setText(widget.text()) + layout.replaceWidget(widget, replacement) + widget.setParent(None) + widget.deleteLater() + fields[index] = (replacement, name, converter) + + def _formatAxisValue(self, value): + """Format an abscissa or ordinate value for the position readout. + + Pixel coordinates need no more than sub-pixel precision, while the + momentum transfer shown while the Q-plot is on needs the same five + decimals as the hkl and Q triplets. + + .. note:: + GUI-only. + """ + return f"{value:.5f}" if self._qPlotActive else f"{value:.2f}" + + def _renamePositionField(self, previous, name): + """Rename one position readout field title in place. + + silx builds the field titles once, so the existing title widget is + renamed rather than the readout rebuilt. The displayed values follow + automatically, because the converters read the current state when they + are called. + + .. note:: + GUI-only. + """ + if previous == name: + return + positionInfo = self.getPositionInfoWidget() + if positionInfo is None: + return + old = f"{previous}:" + for widget in positionInfo.findChildren(qt.QLabel): + if widget.text() == old: + widget.setText(f"{name}:") + break + + def _updatePositionFieldNames(self): + """Apply the field titles that match the current frame and plot mode. + + .. note:: + GUI-only. + """ + qFieldNameNew = qFieldName(self._qFrame) + axisNames = axisFieldNames(self._qFrame, self._qPlotActive) + if qFieldNameNew == self._qFieldName and axisNames == self._axisFieldNames: + return + + self._renamePositionField(self._qFieldName, qFieldNameNew) + for previous, name in zip(self._axisFieldNames, axisNames): + self._renamePositionField(previous, name) + self._qFieldName = qFieldNameNew + self._axisFieldNames = axisNames + + positionInfo = self.getPositionInfoWidget() + if positionInfo is not None: + positionInfo.updateInfo() + + def setQFrame(self, frame): + """Select the reciprocal-space frame the readout describes. + + .. note:: + GUI-only. + """ + self._qFrame = frame + self._updatePositionFieldNames() + + def setQPlotActive(self, active): + """Tell the readout whether the axes hold momentum transfer or pixels. + + .. note:: + GUI-only. + """ + self._qPlotActive = bool(active) + self._updatePositionFieldNames() + def setXyHKLconverter(self, xyHKLConverter): """Set the pixel-to-hkl converter used by the plot readout. diff --git a/orgui/app/qconversion.py b/orgui/app/qconversion.py new file mode 100644 index 0000000..b457d63 --- /dev/null +++ b/orgui/app/qconversion.py @@ -0,0 +1,814 @@ +"""Reciprocal-space conversion of detector images for the Q-plot. + +.. warning:: + + **Experimental application code.** This module backs the experimental + ``Q-plot`` action of the GUI. It is deliberately kept out of + :mod:`orgui.datautils.xrayutils` so that it is not mistaken for the + production reciprocal-space code. The authoritative per-pixel calculation is + :meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.QAlpha`; everything + here only exists to display a whole image on a regular grid, and its + conventions may still change. + +The image is rebinned by pyFAI's ``FiberIntegrator``. pyFAI's built-in fiber +units cannot express orGUI's geometry, for two independent reasons: + +* ``sample_orientation`` is the EXIF dihedral group acting on the image array, + so it only spans quarter turns. orGUI's azimuthal reference is continuous and + is therefore passed through ``tilt_angle``, which rotates about the incident + beam. ``sample_orientation`` stays at the identity. +* pyFAI composes ``incident_angle`` and ``tilt_angle`` extrinsically as + ``R_x(-tilt) @ R_y(incidence)``, i.e. about fixed axes. In orGUI the azimuth + rotates the alpha rotation axis itself, so the incidence rotation has to be + applied about the *already rotated* axis, ``R_y(incidence) @ R_x(-tilt)``. + The two agree only when one of the angles vanishes. + +This module therefore builds its own ``pyFAI.units.UnitFiber`` objects; pyFAI +still performs the rebinning. + +Frames +------ + +``QAlpha`` returns the momentum transfer in the alpha frame, the frame of the +sample surface. The other frames are reached with the Vlieg rotation matrices, +following ``hkl = UB^-1 PHI^-1 CHI^-1 OMEGA^-1 Q_alpha``. ``Q_cryst`` +additionally undoes the orientation matrix, so it holds the Cartesian +reciprocal-space vector of the crystal, ``B`` times ``hkl``. For every frame the +out-of-plane component is the ``z`` axis of that frame and the in-plane +component is the radial component in its ``xy`` plane. + +Performance +----------- + +Detector images are large, so the whole conversion is collapsed into a single +affine relation before touching any pixel data. With +``inv = 1 / |(x, y, z)|`` every component reduces to + +.. code-block:: text + + q_j = k * ( (G[j, 0] * x + G[j, 1] * y + G[j, 2] * z) * inv - c[j] ) + +where ``G`` and ``c`` absorb the sample orientation map, the beam and incidence +rotations, the axis relabelling and the frame rotation. The compiled kernel in +``orgui/app/cpp/qconversion_cpp.cpp`` evaluates both displayed quantities in one +pass; a plain numpy fallback is used when the extension is unavailable. + +``numexpr`` is deliberately **not** used here: version 2.11.0 returns wrong +results for a small fraction of multi-threaded evaluations of this kind of +expression, which is not acceptable for a coordinate transform. +""" + +import importlib +import importlib.util +from pathlib import Path + +import numpy as np + +from ..datautils.xrayutils import HKLVlieg + +try: + from pyFAI import units as pyFAI_units + from pyFAI.integrator.fiber import FiberIntegrator + + HAS_FIBER_INTEGRATOR = True +except ImportError: # pyFAI < 2025.1 + HAS_FIBER_INTEGRATOR = False + + +# orGUI's azimuthal reference and pyFAI's fiber convention differ by a quarter +# turn about the incident beam. +AZIMUTH_OFFSET = np.pi / 2 + +# Kept at the identity on purpose, see the module docstring. +SAMPLE_ORIENTATION = 1 + +# Relative slack allowed when inverting the projection, so that a coordinate +# right on the horizon of the accessible region is not rejected by round-off. +_EWALD_TOLERANCE = 1e-9 + +FRAMES = ("Q_alpha", "Q_lab", "Q_omega", "Q_chi", "Q_phi", "Q_cryst") +DEFAULT_FRAME = "Q_alpha" + +# Frames reached through the omega rotation. They are only defined for a single +# image, because a maximum or sum image combines many omega angles. +FRAMES_REQUIRING_OMEGA = ("Q_omega", "Q_chi", "Q_phi", "Q_cryst") + +# Frames that additionally need the orientation matrix U. +FRAMES_REQUIRING_U = ("Q_cryst",) + +FRAME_LABELS = { + "Q_alpha": "Q alpha (surface)", + "Q_lab": "Q lab", + "Q_omega": "Q omega", + "Q_chi": "Q chi", + "Q_phi": "Q phi", + "Q_cryst": "Q crystal", +} + +# (beam, horz, vert) -> orGUI's (Qx, Qy, Qz) +_AXIS_RELABEL = np.array([[0.0, -1.0, 0.0], [1.0, 0.0, 0.0], [0.0, 0.0, 1.0]]) + +_ORIENTATION_BASIS = { + "x": (1.0, 0.0), + "y": (0.0, 1.0), + "(-x)": (-1.0, 0.0), + "(-y)": (0.0, -1.0), +} + +_KERNEL = None +_KERNEL_LOADED = False + + +def _import_kernel(): + """Import the packaged extension or an in-tree Meson build artifact.""" + try: + return importlib.import_module("orgui.app._qconversion_cpp") + except ModuleNotFoundError as package_error: + repo_root = Path(__file__).resolve().parents[2] + candidates = sorted(repo_root.glob("build/**/_qconversion_cpp*.pyd")) + candidates += sorted(repo_root.glob("build/**/_qconversion_cpp*.so")) + if not candidates: + raise package_error + spec = importlib.util.spec_from_file_location( + "_qconversion_cpp", candidates[-1] + ) + if spec is None or spec.loader is None: + raise package_error + module = importlib.util.module_from_spec(spec) + spec.loader.exec_module(module) + return module + + +def kernel(): + """Return the compiled kernel, or ``None`` when it is unavailable.""" + global _KERNEL, _KERNEL_LOADED + if not _KERNEL_LOADED: + try: + _KERNEL = _import_kernel() + except Exception: + _KERNEL = None + _KERNEL_LOADED = True + return _KERNEL + + +def hasKernel(): + """Whether the compiled conversion kernel is available.""" + return kernel() is not None + + +def frameShortName(frame): + """Return the frame name without the leading ``Q_``, for compact labels. + + :param str frame: one of :data:`FRAMES`. + :returns: for example ``"alpha"`` for ``"Q_alpha"``. + """ + if frame not in FRAMES: + raise ValueError(f"Unknown frame {frame!r}, expected one of {FRAMES}") + return frame[len("Q_") :] + + +def tiltAngleFromAzimuth(azimuth): + """Return the pyFAI ``tilt_angle`` reproducing an orGUI azimuthal reference. + + :param float azimuth: azimuthal reference of the detector calibration, rad. + :returns: tilt angle in rad. + """ + return -(azimuth + AZIMUTH_OFFSET) + + +def frameMatrix(frame, alpha=0.0, omega=0.0, chi=0.0, phi=0.0, U=None): + """Rotation taking a vector from the alpha frame into ``frame``. + + :param str frame: one of :data:`FRAMES`. + :param float alpha: incidence angle in rad. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :returns: 3x3 rotation matrix. + """ + if frame not in FRAMES: + raise ValueError(f"Unknown frame {frame!r}, expected one of {FRAMES}") + if frame in FRAMES_REQUIRING_U and U is None: + raise ValueError(f"The {frame} frame requires the orientation matrix U") + + ALPHA, _, _, OMEGA, CHI, PHI = HKLVlieg.createVliegMatrices( + [alpha, None, None, omega, chi, phi] + ) + if frame == "Q_alpha": + return np.identity(3) + if frame == "Q_lab": + return np.asarray(ALPHA) + + matrix = np.asarray(OMEGA).T + if frame == "Q_omega": + return matrix + matrix = np.asarray(CHI).T @ matrix + if frame == "Q_chi": + return matrix + matrix = np.asarray(PHI).T @ matrix + if frame == "Q_phi": + return matrix + # Q_cryst: undo the orientation matrix as well, so that the result is the + # Cartesian reciprocal-space vector of the crystal, i.e. B times hkl. + return np.linalg.inv(np.asarray(U)) @ matrix + + +def _orientationMatrix(sample_orientation): + """Embed the EXIF orientation map as a 3x3 matrix on ``(x, y, z)``.""" + mapping = pyFAI_units.MAPS_SAMPLE_ORIENTATION[sample_orientation] + matrix = np.zeros((3, 3)) + matrix[0, 2] = 1.0 # the beam component picks z + matrix[1, :2] = _ORIENTATION_BASIS[mapping["x"]] + matrix[2, :2] = _ORIENTATION_BASIS[mapping["y"]] + return matrix + + +def conversionCoefficients( + frame, + incident_angle=0.0, + tilt_angle=0.0, + sample_orientation=SAMPLE_ORIENTATION, + omega=0.0, + chi=0.0, + phi=0.0, + U=None, +): + """Collapse the whole conversion into ``(G, c)``. + + The momentum transfer of a pixel at ``(x, y, z)`` is then + ``k * (G @ (x, y, z) / |(x, y, z)| - c)`` in inverse nanometre, with + ``k = 2 pi / wavelength``. + + :returns: ``(G, c)`` with shapes ``(3, 3)`` and ``(3,)``. + """ + incidence = np.array( + [ + [np.cos(incident_angle), 0.0, np.sin(incident_angle)], + [0.0, 1.0, 0.0], + [-np.sin(incident_angle), 0.0, np.cos(incident_angle)], + ] + ) + tilt = np.array( + [ + [1.0, 0.0, 0.0], + [0.0, np.cos(tilt_angle), np.sin(tilt_angle)], + [0.0, -np.sin(tilt_angle), np.cos(tilt_angle)], + ] + ) + combined = ( + frameMatrix(frame, alpha=incident_angle, omega=omega, chi=chi, phi=phi, U=U) + @ _AXIS_RELABEL + @ incidence + @ tilt + ) + return ( + np.ascontiguousarray(combined @ _orientationMatrix(sample_orientation)), + np.ascontiguousarray(combined[:, 0]), + ) + + +def _ipOopNumpy(x, y, z, G, c, k): + """Fallback used when the compiled kernel is unavailable.""" + scratch = np.empty_like(x) + accumulator = np.empty_like(x) + + def component(index): + np.multiply(x, G[index, 0], out=accumulator) + np.multiply(y, G[index, 1], out=scratch) + np.add(accumulator, scratch, out=accumulator) + np.multiply(z, G[index, 2], out=scratch) + np.add(accumulator, scratch, out=accumulator) + np.multiply(accumulator, inverse_norm, out=accumulator) + np.subtract(accumulator, c[index], out=accumulator) + np.multiply(accumulator, k, out=accumulator) + return accumulator.copy() + + inverse_norm = np.empty_like(x) + np.multiply(x, x, out=inverse_norm) + np.multiply(y, y, out=scratch) + np.add(inverse_norm, scratch, out=inverse_norm) + np.multiply(z, z, out=scratch) + np.add(inverse_norm, scratch, out=inverse_norm) + np.sqrt(inverse_norm, out=inverse_norm) + np.reciprocal(inverse_norm, out=inverse_norm) + + q_x = component(0) + q_y = component(1) + q_oop = component(2) + np.multiply(q_x, q_x, out=scratch) + np.multiply(q_y, q_y, out=q_y) + np.add(scratch, q_y, out=scratch) + np.sqrt(scratch, out=scratch) + # signed by the in-plane direction of orGUI's delta angle + return np.copysign(scratch, q_x), q_oop + + +class _ConversionCache: + """Remember the last conversion. + + pyFAI evaluates the in-plane and the out-of-plane unit separately but hands + both the very same position arrays, so the image is converted once instead + of twice. + """ + + def __init__(self): + self._arrays = None + self._signature = None + self._value = None + + def get(self, x, y, z, G, c, k): + signature = (G.tobytes(), c.tobytes(), float(k)) + if self._arrays is not None and self._signature == signature: + cached_x, cached_y, cached_z = self._arrays + if cached_x is x and cached_y is y and cached_z is z: + return self._value + + compiled = kernel() + if compiled is not None: + value = compiled.q_ip_oop(x, y, z, G, c, float(k)) + else: + value = _ipOopNumpy( + np.ascontiguousarray(x, dtype=np.float64), + np.ascontiguousarray(y, dtype=np.float64), + np.ascontiguousarray(z, dtype=np.float64), + G, + c, + k, + ) + self._arrays = (x, y, z) + self._signature = signature + self._value = value + return value + + def clear(self): + self._arrays = None + self._signature = None + self._value = None + + +_CONVERSION = _ConversionCache() + + +class _IntegratorCache: + """Reuse the pyFAI integrator across conversions of the same geometry. + + A fresh ``FiberIntegrator`` starts with empty grazing-incidence parameters, + so passing any non-zero incidence or tilt angle makes pyFAI reset and + recompute all of its cached arrays. Keeping the integrator avoids that on + every image change. + + The cache key contains the collapsed conversion coefficients, so a hit + means the unit equations are byte for byte identical. That matters because + pyFAI caches the per-pixel unit arrays under the unit *name*: reusing an + integrator with a same-named but differently parameterised unit would + silently return the previous image's coordinates. + """ + + def __init__(self): + self._key = None + self._detector = None + self._integrator = None + + def get(self, detectorCal, G, c): + key = ( + float(detectorCal.dist), + float(detectorCal.poni1), + float(detectorCal.poni2), + float(detectorCal.rot1), + float(detectorCal.rot2), + float(detectorCal.rot3), + float(detectorCal.wavelength), + G.tobytes(), + c.tobytes(), + ) + if self._key == key and self._detector is detectorCal.detector: + return self._integrator + + integrator = FiberIntegrator( + dist=detectorCal.dist, + poni1=detectorCal.poni1, + poni2=detectorCal.poni2, + wavelength=detectorCal.wavelength, + rot1=detectorCal.rot1, + rot2=detectorCal.rot2, + rot3=detectorCal.rot3, + detector=detectorCal.detector, + ) + self._key = key + self._detector = detectorCal.detector + self._integrator = integrator + return integrator + + def clear(self): + self._key = None + self._detector = None + self._integrator = None + + +_INTEGRATOR = _IntegratorCache() + + +def qIpOop( + frame, + x, + y, + z, + wavelength, + incident_angle=0.0, + tilt_angle=0.0, + sample_orientation=SAMPLE_ORIENTATION, + omega=0.0, + chi=0.0, + phi=0.0, + U=None, +): + """In-plane and out-of-plane momentum transfer, in inverse nanometre.""" + G, c = conversionCoefficients( + frame, + incident_angle=incident_angle, + tilt_angle=tilt_angle, + sample_orientation=sample_orientation, + omega=omega, + chi=chi, + phi=phi, + U=U, + ) + return _CONVERSION.get(x, y, z, G, c, 2.0e-9 * np.pi / wavelength) + + +def qComponents( + frame, + x, + y, + z, + wavelength, + incident_angle=0.0, + tilt_angle=0.0, + sample_orientation=SAMPLE_ORIENTATION, + omega=0.0, + chi=0.0, + phi=0.0, + U=None, +): + """All three Cartesian components expressed in ``frame``. + + This is the readable reference implementation; the Q-plot itself only needs + the two quantities returned by :func:`qIpOop`. + + :returns: ``(Qx, Qy, Qz)`` of ``frame`` in inverse nanometre. + """ + G, c = conversionCoefficients( + frame, + incident_angle=incident_angle, + tilt_angle=tilt_angle, + sample_orientation=sample_orientation, + omega=omega, + chi=chi, + phi=phi, + U=U, + ) + k = 2.0e-9 * np.pi / wavelength + inverse_norm = 1.0 / np.sqrt(x * x + y * y + z * z) + return tuple( + k * ((G[j, 0] * x + G[j, 1] * y + G[j, 2] * z) * inverse_norm - c[j]) + for j in range(3) + ) + + +def beamDirection(frame, alpha=0.0, omega=0.0, chi=0.0, phi=0.0, U=None): + """Unit vector along the incident beam, expressed in ``frame``. + + In the alpha frame the beam is ``(0, cos alpha, -sin alpha)``, which is the + vector :meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.QAlpha` + subtracts from the outgoing direction. It is the same vector as the ``c`` + returned by :func:`conversionCoefficients`, which is the beam column of the + collapsed conversion. + + :param str frame: one of :data:`FRAMES`. + :param float alpha: incidence angle in rad. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :returns: unit vector of shape ``(3,)``. + """ + rotation = frameMatrix(frame, alpha=alpha, omega=omega, chi=chi, phi=phi, U=U) + return rotation @ np.array([0.0, np.cos(alpha), -np.sin(alpha)]) + + +def qFromIpOop( + frame, + q_ip, + q_oop, + k, + alpha=0.0, + omega=0.0, + chi=0.0, + phi=0.0, + U=None, +): + """Momentum transfer vectors that project onto one Q-plot coordinate. + + This inverts the projection performed by :func:`qIpOop`, which keeps the + out-of-plane component and collapses the two in-plane components into their + radial distance, signed by ``Qx``. The missing information is recovered from + the Ewald condition: with ``q = Q / k`` and the beam direction ``c`` of + :func:`beamDirection`, the outgoing direction ``q + c`` is a unit vector, so + + .. code-block:: text + + |q|^2 + 2 q.c = 0 + + which is *linear* in ``qx`` and ``qy`` once ``qz`` and ``|q|`` are known. + Intersecting that line with the circle ``qx^2 + qy^2 = (q_ip / k)^2`` leaves + at most two solutions, of which only those with ``sign(Qx) == sign(q_ip)`` + reproduce the sign convention of :func:`qIpOop`. + + In the alpha frame the beam has no ``x`` component, the line is parallel to + the ``x`` axis and the sign condition always leaves a single solution. In + the rotated frames the beam generally does have one, and both solutions can + carry the same sign; the caller has to decide between them, for example by + keeping the one that falls onto the detector. + + ``k``, ``q_ip`` and ``q_oop`` must share the same unit, and the result is + returned in that unit. Note that :func:`qIpOop` works in inverse nanometre + while orGUI's own reciprocal space is in inverse Angstrom. + + :param str frame: one of :data:`FRAMES`. + :param float q_ip: in-plane coordinate, signed as in :func:`qIpOop`. + :param float q_oop: out-of-plane coordinate. + :param float k: magnitude of the wave vector. + :param float alpha: incidence angle in rad. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :returns: + Tuple of zero, one or two momentum transfer vectors of shape ``(3,)``, + expressed in ``frame``. Empty when the coordinate is outside the part + of reciprocal space the Ewald sphere can reach. + :rtype: tuple[numpy.ndarray, ...] + """ + q_ip = float(q_ip) + q_oop = float(q_oop) + k = float(k) + if not np.isfinite([q_ip, q_oop, k]).all() or k == 0.0: + return () + + c = beamDirection(frame, alpha=alpha, omega=omega, chi=chi, phi=phi, U=U) + + # everything below is in units of k, so that the Ewald sphere is the unit + # sphere around the tip of the incident wave vector + radius = abs(q_ip) / k + q_z = q_oop / k + squared_length = radius * radius + q_z * q_z + + # c[0] * q_x + c[1] * q_y = offset, from the Ewald condition + offset = -0.5 * squared_length - c[2] * q_z + normal_length = c[0] * c[0] + c[1] * c[1] + if normal_length <= 0.0: + # the beam is along z, so the condition no longer constrains the plane + return () + + # distance from the origin to the line, compared with the circle radius + half_chord_squared = radius * radius - offset * offset / normal_length + if half_chord_squared < 0.0: + # tolerate round-off on the horizon of the accessible region + if half_chord_squared < -_EWALD_TOLERANCE * max(radius * radius, 1.0): + return () + half_chord_squared = 0.0 + half_chord = np.sqrt(half_chord_squared / normal_length) + + # foot of the perpendicular from the origin onto the line + foot_x = c[0] * offset / normal_length + foot_y = c[1] * offset / normal_length + + solutions = [] + for sign in (1.0, -1.0): + q_x = foot_x - sign * c[1] * half_chord + q_y = foot_y + sign * c[0] * half_chord + # np.copysign in qIpOop signs the radius by Qx, so a solution is only + # valid when its Qx carries the sign of the requested in-plane value + if q_ip != 0.0 and np.sign(q_x) != np.sign(q_ip): + continue + solutions.append(k * np.array([q_x, q_y, q_z])) + if half_chord == 0.0: + # the line touches the circle, both signs give the same point + break + return tuple(solutions) + + +def detectorAngles( + Q, k, frame=DEFAULT_FRAME, alpha=0.0, omega=0.0, chi=0.0, phi=0.0, U=None +): + """Detector angles producing a momentum transfer expressed in ``frame``. + + The momentum transfer is rotated back into the alpha frame, where + ``Q / k + c`` is the outgoing unit direction + ``(sin delta cos gamma, cos delta cos gamma, sin gamma)`` of + :meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.QAlpha`. + + :param Q: momentum transfer of ``frame``, shape ``(3,)``. + :param float k: magnitude of the wave vector, in the unit of ``Q``. + :param str frame: one of :data:`FRAMES`. + :param float alpha: incidence angle in rad. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :returns: ``(delta, gamma)`` in rad, or ``None`` when ``Q`` is not on the + Ewald sphere and the direction cannot be normalised. + :rtype: tuple[float, float] or None + """ + rotation = frameMatrix(frame, alpha=alpha, omega=omega, chi=chi, phi=phi, U=U) + q_alpha = np.linalg.inv(rotation) @ np.asarray(Q, dtype=np.float64) + direction = q_alpha / float(k) + np.array([0.0, np.cos(alpha), -np.sin(alpha)]) + + norm = np.linalg.norm(direction) + if not np.isfinite(norm) or norm == 0.0: + return None + # the caller may have built Q from a rounded plot coordinate, so normalise + # instead of insisting on an exactly unit direction + direction = direction / norm + + gamma = np.arcsin(np.clip(direction[2], -1.0, 1.0)) + delta = np.arctan2(direction[0], direction[1]) + return float(delta), float(gamma) + + +def fiberUnits(frame=DEFAULT_FRAME, omega=0.0, chi=0.0, phi=0.0, U=None): + """Build the in-plane and out-of-plane pyFAI units for ``frame``. + + :returns: ``(unit_ip, unit_oop)`` as ``pyFAI.units.UnitFiber``. + """ + if not HAS_FIBER_INTEGRATOR: + raise RuntimeError("pyFAI >= 2025.1 with FiberIntegrator is required") + if frame not in FRAMES: + raise ValueError(f"Unknown frame {frame!r}, expected one of {FRAMES}") + if frame in FRAMES_REQUIRING_U and U is None: + raise ValueError(f"The {frame} frame requires the orientation matrix U") + + def evaluate(index, x, y, z, wavelength, incident_angle, tilt, orientation): + return qIpOop( + frame, + x, + y, + z, + wavelength, + incident_angle=incident_angle, + tilt_angle=tilt, + sample_orientation=orientation, + omega=omega, + chi=chi, + phi=phi, + U=U, + )[index] + + def equation_ip( + x, + y, + z, + wavelength, + incident_angle=0.0, + tilt_angle=0.0, + sample_orientation=SAMPLE_ORIENTATION, + ): + return evaluate( + 0, x, y, z, wavelength, incident_angle, tilt_angle, sample_orientation + ) + + def equation_oop( + x, + y, + z, + wavelength, + incident_angle=0.0, + tilt_angle=0.0, + sample_orientation=SAMPLE_ORIENTATION, + ): + return evaluate( + 1, x, y, z, wavelength, incident_angle, tilt_angle, sample_orientation + ) + + label = FRAME_LABELS[frame] + unit_ip = pyFAI_units.UnitFiber( + f"qip_{frame}_A^-1", + scale=0.1, + label=rf"{label} $q_\parallel$ ($\AA^{{-1}}$)", + equation=equation_ip, + short_name=f"qip_{frame}", + unit_symbol=r"\AA^{-1}", + positive=False, + ) + unit_oop = pyFAI_units.UnitFiber( + f"qoop_{frame}_A^-1", + scale=0.1, + label=rf"{label} $q_\perp$ ($\AA^{{-1}}$)", + equation=equation_oop, + short_name=f"qoop_{frame}", + unit_symbol=r"\AA^{-1}", + positive=False, + ) + return unit_ip, unit_oop + + +def outOfPlaneIncreasesWithRow( + detectorCal, alpha, frame=DEFAULT_FRAME, omega=0.0, chi=0.0, phi=0.0, U=None +): + """Whether the out-of-plane coordinate grows with the detector row index. + + The detector image is displayed with the row index increasing downwards, + while the Q-plot draws the out-of-plane coordinate upwards. The two only + agree visually when the out-of-plane coordinate *de*creases with the row + index, which is the ordinary upward-scattering case. + + An azimuthal reference that points the surface normal the other way -- the + inverted, downward-scattering geometry used at ID31 -- reverses this, and + the reciprocal-space image then appears vertically mirrored with respect to + the detector image. The caller can invert the Q-plot's ordinate to + compensate. + + Only the sign matters here, so the wave vector and the beam offset drop + out of the comparison and no wavelength is needed. + + :param detectorCal: :class:`DetectorCalibration.Detector2D_SXRD`. + :param float alpha: incidence angle in rad. + :param str frame: one of :data:`FRAMES`. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :rtype: bool + """ + shape = detectorCal.detector.shape + rows = np.array([0.0, float(shape[0] - 1)]) + cols = np.full(2, 0.5 * (shape[1] - 1)) + param = np.array( + [ + detectorCal.dist, + detectorCal.poni1, + detectorCal.poni2, + detectorCal.rot1, + detectorCal.rot2, + detectorCal.rot3, + ], + dtype=np.float64, + ) + z, y, x = detectorCal.calc_pos_zyx(d1=rows, d2=cols, param=param, do_parallax=True) + + G, _ = conversionCoefficients( + frame, + incident_angle=alpha, + tilt_angle=tiltAngleFromAzimuth(detectorCal.getAzimuthalReference()), + omega=omega, + chi=chi, + phi=phi, + U=U, + ) + # the out-of-plane component up to the positive factor k and a constant + # offset, both of which cancel in the comparison + direction = np.stack([x, y, z], axis=-1) + direction = direction / np.linalg.norm(direction, axis=-1, keepdims=True) + out_of_plane = direction @ G[2] + return bool(out_of_plane[1] > out_of_plane[0]) + + +def integrateImage( + detectorCal, + image, + alpha, + frame=DEFAULT_FRAME, + omega=0.0, + chi=0.0, + phi=0.0, + U=None, + npt=1000, +): + """Rebin ``image`` onto a regular grid of momentum transfer. Experimental. + + :param detectorCal: :class:`DetectorCalibration.Detector2D_SXRD`. + :param numpy.ndarray image: detector image. + :param float alpha: incidence angle in rad. + :param str frame: one of :data:`FRAMES`. + :param float omega, chi, phi: sample circle angles in rad. + :param U: orientation matrix, required for :data:`FRAMES_REQUIRING_U`. + :param int npt: number of bins along each axis. + :returns: the pyFAI 2D integration result. + """ + if not HAS_FIBER_INTEGRATOR: + raise RuntimeError("pyFAI >= 2025.1 with FiberIntegrator is required") + + unit_ip, unit_oop = fiberUnits(frame, omega=omega, chi=chi, phi=phi, U=U) + tilt = tiltAngleFromAzimuth(detectorCal.getAzimuthalReference()) + # the coefficients identify the conversion completely, see _IntegratorCache + G, c = conversionCoefficients( + frame, + incident_angle=alpha, + tilt_angle=tilt, + omega=omega, + chi=chi, + phi=phi, + U=U, + ) + integrator = _INTEGRATOR.get(detectorCal, G, c) + return integrator.integrate2d_grazing_incidence( + image, + npt_ip=npt, + npt_oop=npt, + sample_orientation=SAMPLE_ORIENTATION, + incident_angle=alpha, + tilt_angle=tilt, + unit_ip=unit_ip, + unit_oop=unit_oop, + ) diff --git a/orgui/app/test/test_parameter_dialog_cancel.py b/orgui/app/test/test_parameter_dialog_cancel.py new file mode 100644 index 0000000..2c9642c --- /dev/null +++ b/orgui/app/test/test_parameter_dialog_cancel.py @@ -0,0 +1,110 @@ +"""Cancelling a parameter dialog must not leave the edited configuration active. + +The machine and crystal parameter editors apply every change immediately, so +restoring the widget values is not sufficient: the saved parameters have to be +applied again when the dialog is cancelled or reset. +""" + +import pytest + +QUBCalculator = pytest.importorskip( + "orgui.app.QUBCalculator", reason="Qt bindings are required" +) + + +class _StubSignal: + def __init__(self): + self.emitted = [] + + def emit(self, *args): + self.emitted.append(args) + + +class _StubMachineParams: + def __init__(self): + self.restored = [] + self.sigMachineParamsChanged = _StubSignal() + + def setValues(self, params): + self.restored.append(params) + + +class _StubMachineDialog: + def __init__(self, saved): + self.savedParams = saved + self.machineparams = _StubMachineParams() + + +class _StubCrystalParams: + def __init__(self): + self.restored = [] + self.sigCrystalParamsChanged = _StubSignal() + + def setValues(self, crystal, n): + self.restored.append((crystal, n)) + + +class _StubCrystalDialog: + def __init__(self, saved): + self.savedParams = saved + self.crystalparams = _StubCrystalParams() + + +def _machine_params(): + return { + "diffractometer": {"mu": 0.1, "chi": 0.2, "phi": 0.3}, + "source": {"E": 78.0}, + "SXRD_geometry": object(), + } + + +def test_reset_reapplies_the_saved_machine_parameters(): + saved = _machine_params() + dialog = _StubMachineDialog(saved) + + QUBCalculator.QMachineParametersDialog.resetParameters(dialog) + + assert dialog.machineparams.restored == [saved] + # applied again, otherwise the edited configuration would stay active + assert dialog.machineparams.sigMachineParamsChanged.emitted == [(saved,)] + + +def test_machine_parameters_are_reapplied_unrounded(): + """The saved snapshot is emitted, not the value read back from a spin box.""" + saved = _machine_params() + dialog = _StubMachineDialog(saved) + + QUBCalculator.QMachineParametersDialog.resetParameters(dialog) + + (emitted,) = dialog.machineparams.sigMachineParamsChanged.emitted[0] + assert emitted is saved + + +def test_reset_without_a_snapshot_does_nothing(): + dialog = _StubMachineDialog(None) + + QUBCalculator.QMachineParametersDialog.resetParameters(dialog) + + assert dialog.machineparams.restored == [] + assert dialog.machineparams.sigMachineParamsChanged.emitted == [] + + +def test_reset_reapplies_the_saved_crystal(): + crystal, refraction_index = object(), 0.9999 + dialog = _StubCrystalDialog((crystal, refraction_index)) + + QUBCalculator.QCrystalParameterDialog.resetParameters(dialog) + + assert dialog.crystalparams.restored == [(crystal, refraction_index)] + assert dialog.crystalparams.sigCrystalParamsChanged.emitted == [ + (crystal, refraction_index) + ] + + +def test_crystal_reset_without_a_snapshot_does_nothing(): + dialog = _StubCrystalDialog(None) + + QUBCalculator.QCrystalParameterDialog.resetParameters(dialog) + + assert dialog.crystalparams.restored == [] + assert dialog.crystalparams.sigCrystalParamsChanged.emitted == [] diff --git a/orgui/app/test/test_position_readout.py b/orgui/app/test/test_position_readout.py new file mode 100644 index 0000000..7640452 --- /dev/null +++ b/orgui/app/test/test_position_readout.py @@ -0,0 +1,255 @@ +"""The plot status line reports HKL and the selected reciprocal-space frame. + +The readout shows ``HKL`` and ``Q[]`` as bracketed triplets, and the +``Q`` field is relabelled when a different frame is selected. +""" + +import numpy as np +import pytest + +from orgui.app import qconversion +from orgui.datautils.xrayutils import HKLVlieg + +orGUI = pytest.importorskip("orgui.app.orGUI", reason="Qt bindings are required") + + +ENERGY = 78.0 # keV + +# refracted [alpha, delta, gamma, omega, chi, phi] in rad +POSITION = [ + np.deg2rad(0.6), + np.deg2rad(12.0), + np.deg2rad(4.0), + np.deg2rad(23.0), + np.deg2rad(3.0), + np.deg2rad(-5.0), +] + +EXPECTED_FIELD_NAMES = { + "Q_alpha": "Q[alpha]", + "Q_lab": "Q[lab]", + "Q_omega": "Q[omega]", + "Q_chi": "Q[chi]", + "Q_phi": "Q[phi]", + "Q_cryst": "Q[cryst]", +} + + +def _ub_calculator(with_orientation=True): + lattice = HKLVlieg.Lattice(np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0])) + ub = HKLVlieg.UBCalculator(lattice, ENERGY) + if with_orientation: + ub.defaultU() + return ub + + +class _StubUBCalculatorWidget: + def __init__(self, ub): + self.ubCal = ub + self.angles = HKLVlieg.VliegAngles(ub) + + +class _StubOrGUI: + """The readout helper only needs the frame and the UB calculator.""" + + def __init__(self, frame, ub): + self._frame = frame + self.ubcalc = _StubUBCalculatorWidget(ub) + + def _selectedQFrame(self): + return self._frame + + +def test_triplet_is_formatted_as_a_bracketed_vector(): + """Five decimals are the minimum useful precision for hkl and for Q.""" + formatted = orGUI.formatPositionTriplet([1.23456, -0.5, 2.0]) + + assert formatted == "[1.23456 -0.50000 2.00000]" + + +def test_triplet_template_is_wide_enough_for_the_format(): + """The size hint must cover the widest value the format can produce.""" + widest = orGUI.formatPositionTriplet([-99.999994, -99.999994, -99.999994]) + + assert len(widest) <= len(orGUI.POSITION_TRIPLET_TEMPLATE) + + +@pytest.mark.parametrize("values", [[np.nan, 1.0, 2.0], [np.inf, 1.0, 2.0]]) +def test_triplet_without_a_scan_shows_a_placeholder(values): + assert orGUI.formatPositionTriplet(values) == "------" + + +def test_every_frame_has_an_expected_field_name(): + assert set(EXPECTED_FIELD_NAMES) == set(qconversion.FRAMES) + + +@pytest.mark.parametrize("frame, expected", sorted(EXPECTED_FIELD_NAMES.items())) +def test_field_name_follows_the_frame(frame, expected): + assert orGUI.qFieldName(frame) == expected + + +def test_alpha_frame_reports_qalpha(): + stub = _StubOrGUI("Q_alpha", _ub_calculator()) + + q = orGUI.orGUI._qInSelectedFrame(stub, POSITION) + + expected = stub.ubcalc.angles.QAlpha(*POSITION[:3]) + np.testing.assert_allclose(q, expected, atol=1e-12) + + +def test_crystal_frame_reports_b_times_hkl(): + """Q_cryst is the crystal Cartesian frame, so it must equal B times hkl.""" + ub = _ub_calculator() + stub = _StubOrGUI("Q_cryst", ub) + + q = orGUI.orGUI._qInSelectedFrame(stub, POSITION) + + hkl = np.array(stub.ubcalc.angles.anglesToHkl(*POSITION)) + np.testing.assert_allclose(q, ub.lattice.B_mat @ hkl, atol=1e-9) + + +def test_crystal_frame_without_an_orientation_matrix_is_unavailable(): + stub = _StubOrGUI("Q_cryst", _ub_calculator(with_orientation=False)) + + q = orGUI.orGUI._qInSelectedFrame(stub, POSITION) + + assert np.all(np.isnan(q)) + assert orGUI.formatPositionTriplet(q) == "------" + + +class _StubTitle: + """Stand-in for the bold field title that silx creates once.""" + + def __init__(self, text): + self._text = text + + def text(self): + return self._text + + def setText(self, text): + self._text = text + + +class _StubPositionInfo: + def __init__(self, titles): + self._titles = titles + self.updates = 0 + + def findChildren(self, cls): + del cls + return self._titles + + def updateInfo(self): + self.updates += 1 + + +class _StubPlot: + def __init__(self): + self._qFrame = "Q_alpha" + self._qPlotActive = False + self._qFieldName = "Q[alpha]" + self._axisFieldNames = ("X", "Y") + self.titles = [ + _StubTitle("X:"), + _StubTitle("Y:"), + _StubTitle("HKL:"), + _StubTitle("Q[alpha]:"), + ] + self.positionInfo = _StubPositionInfo(self.titles) + + def getPositionInfoWidget(self): + return self.positionInfo + + # exercised through the real implementations + _renamePositionField = orGUI.Plot2DHKL._renamePositionField + _updatePositionFieldNames = orGUI.Plot2DHKL._updatePositionFieldNames + _formatAxisValue = orGUI.Plot2DHKL._formatAxisValue + + +@pytest.mark.parametrize("frame, name", sorted(EXPECTED_FIELD_NAMES.items())) +def test_selecting_a_frame_relabels_only_the_q_field(frame, name): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQFrame(plot, frame) + + assert plot.titles[3].text() == f"{name}:" + assert plot.titles[2].text() == "HKL:" # untouched + assert [t.text() for t in plot.titles[:2]] == ["X:", "Y:"] + assert plot._qFieldName == name + + +def test_reselecting_the_same_frame_is_a_noop(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQFrame(plot, "Q_alpha") + + assert plot.positionInfo.updates == 0 + assert plot.titles[3].text() == "Q[alpha]:" + + +def test_relabelling_refreshes_the_displayed_values(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQFrame(plot, "Q_phi") + + assert plot.positionInfo.updates == 1 + + +def test_axis_fields_are_pixels_while_the_q_plot_is_off(): + assert orGUI.axisFieldNames("Q_alpha", False) == ("X", "Y") + + +@pytest.mark.parametrize("frame", qconversion.FRAMES) +def test_axis_fields_name_the_momentum_transfer_in_the_q_plot(frame): + """The Q-plot axes are momentum transfer, so ``X``/``Y`` would be wrong.""" + short = qconversion.frameShortName(frame) + + assert orGUI.axisFieldNames(frame, True) == ( + f"q_par[{short}]", + f"q_perp[{short}]", + ) + + +def test_switching_on_the_q_plot_relabels_the_axis_fields(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQPlotActive(plot, True) + + assert [t.text() for t in plot.titles[:2]] == [ + "q_par[alpha]:", + "q_perp[alpha]:", + ] + assert plot.titles[3].text() == "Q[alpha]:" # untouched + + +def test_switching_off_the_q_plot_restores_the_pixel_fields(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQPlotActive(plot, True) + orGUI.Plot2DHKL.setQPlotActive(plot, False) + + assert [t.text() for t in plot.titles[:2]] == ["X:", "Y:"] + + +def test_changing_the_frame_follows_through_to_the_axis_fields(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQPlotActive(plot, True) + orGUI.Plot2DHKL.setQFrame(plot, "Q_phi") + + assert [t.text() for t in plot.titles[:2]] == [ + "q_par[phi]:", + "q_perp[phi]:", + ] + assert plot.titles[3].text() == "Q[phi]:" + + +def test_axis_values_gain_precision_in_the_q_plot(): + """Sub-pixel precision is enough for pixels, but not for inverse Angstrom.""" + plot = _StubPlot() + + assert orGUI.Plot2DHKL._formatAxisValue(plot, 1234.5678) == "1234.57" + + orGUI.Plot2DHKL.setQPlotActive(plot, True) + + assert orGUI.Plot2DHKL._formatAxisValue(plot, 1.2345678) == "1.23457" diff --git a/orgui/app/test/test_q_conversion.py b/orgui/app/test/test_q_conversion.py new file mode 100644 index 0000000..8848b99 --- /dev/null +++ b/orgui/app/test/test_q_conversion.py @@ -0,0 +1,568 @@ +"""Equivalence of the Q-plot conversion with orGUI's own QAlpha calculation. + +orGUI describes reciprocal space through the surface angles ``gamma`` and +``delta`` (:meth:`DetectorCalibration.Detector2D_SXRD.surfaceAngles`) combined +with the Vlieg diffraction equation (:meth:`HKLVlieg.VliegAngles.QAlpha`). The +experimental Q-plot instead rebins a whole image onto a regular grid using +pyFAI's ``FiberIntegrator``, driven by :mod:`orgui.app.qconversion`. + +The azimuthal reference is deliberately tested at values that are *not* +multiples of 90 degrees. pyFAI's ``sample_orientation`` flag can only express +quarter turns, so those cases are exactly the ones an implementation based on +that flag alone gets wrong. +""" + +import unittest + +import numpy as np +import pyFAI + +from orgui.app import qconversion +from orgui.datautils.xrayutils import DetectorCalibration, HKLVlieg + +requires_fiber = unittest.skipUnless( + qconversion.HAS_FIBER_INTEGRATOR, + "pyFAI >= 2025.1 with FiberIntegrator is required", +) + +ENERGY = 78.0 # keV + +FLAT_DETECTOR = (0.0, 0.0, 0.0) +TILTED_DETECTOR = (np.deg2rad(1.7), np.deg2rad(-0.9), np.deg2rad(3.4)) + +# not restricted to multiples of 90 degrees on purpose +AZIMUTHS = (0.0, 17.3, 45.0, 90.0, 135.0, 233.7, 341.9) +INCIDENCE_ANGLES = (0.0, 0.08, 0.6, 2.5) + + +def _geometry(azimuth_deg, rotations=FLAT_DETECTOR): + """Return a Detector2D_SXRD with a fully specified geometry.""" + det = DetectorCalibration.Detector2D_SXRD() + det.detector = pyFAI.detector_factory("Pilatus2m") + det.wavelength = (12.39842 / ENERGY) * 1e-10 + det.dist = 0.729 + det.poni1 = 0.21 + det.poni2 = 0.14 + det.rot1, det.rot2, det.rot3 = rotations + det.setAzimuthalReference(np.deg2rad(azimuth_deg)) + return det + + +def _ub_calculator(): + lattice = HKLVlieg.Lattice(np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0])) + ub = HKLVlieg.UBCalculator(lattice, ENERGY) + ub.defaultU() + return ub + + +def _pixel_grid(det, stride=251): + """Flattened pixel coordinates covering the whole detector.""" + shape = det.detector.shape + rows, cols = np.meshgrid( + np.arange(0, shape[0], stride, dtype=np.float64), + np.arange(0, shape[1], stride, dtype=np.float64), + indexing="ij", + ) + return rows.ravel(), cols.ravel() + + +def _pixel_positions(det, rows, cols): + """pyFAI pixel positions as ``(x, y, z)`` = (horizontal, vertical, beam).""" + param = np.array( + [det.dist, det.poni1, det.poni2, det.rot1, det.rot2, det.rot3], + dtype=np.float64, + ) + t3, t1, t2 = det.calc_pos_zyx(d1=rows, d2=cols, param=param, do_parallax=True) + return t2, t1, t3 + + +def _q_from_conversion(det, rows, cols, alpha_i, frame="Q_alpha", **frame_angles): + """Momentum transfer in inverse Angstrom from :mod:`qconversion`.""" + x, y, z = _pixel_positions(det, rows, cols) + components = qconversion.qComponents( + frame, + x, + y, + z, + det.wavelength, + incident_angle=alpha_i, + tilt_angle=qconversion.tiltAngleFromAzimuth(det.getAzimuthalReference()), + sample_orientation=qconversion.SAMPLE_ORIENTATION, + **frame_angles, + ) + return np.stack([component / 10.0 for component in components], axis=-1) + + +@requires_fiber +class TestQAlphaEquivalence(unittest.TestCase): + """The alpha frame must reproduce QAlpha component by component.""" + + def test_matches_qalpha_for_arbitrary_azimuth(self): + angles = HKLVlieg.VliegAngles(_ub_calculator()) + for azimuth_deg in AZIMUTHS: + for rotations in (FLAT_DETECTOR, TILTED_DETECTOR): + for alpha_deg in INCIDENCE_ANGLES: + with self.subTest( + azimuth=azimuth_deg, + tilted=rotations is TILTED_DETECTOR, + alpha_i=alpha_deg, + ): + det = _geometry(azimuth_deg, rotations) + rows, cols = _pixel_grid(det) + alpha_i = np.deg2rad(alpha_deg) + + actual = _q_from_conversion(det, rows, cols, alpha_i) + gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha_i) + expected = angles.QAlpha(alpha_i, delta, gamma) + + np.testing.assert_allclose(actual, expected, atol=1e-9) + + def test_builtin_pyfai_units_cannot_express_a_general_azimuth(self): + """Motivation for the custom units, see :mod:`qconversion`. + + ``sample_orientation`` only spans quarter turns, so no value of it can + reproduce an azimuthal reference of 45 degrees. + """ + from pyFAI import units as pyFAI_units + + angles = HKLVlieg.VliegAngles(_ub_calculator()) + det = _geometry(45.0) + rows, cols = _pixel_grid(det) + alpha_i = np.deg2rad(0.6) + + gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha_i) + expected = angles.QAlpha(alpha_i, delta, gamma) + expected_oop = expected[..., 2] + x, y, z = _pixel_positions(det, rows, cols) + + best = np.inf + for orientation in range(1, 9): + actual_oop = ( + pyFAI_units.eq_qoop( + x=x, + y=y, + z=z, + wavelength=det.wavelength, + incident_angle=alpha_i, + tilt_angle=0.0, + sample_orientation=orientation, + ) + / 10.0 + ) + best = min(best, np.max(np.abs(expected_oop - actual_oop))) + + # the closest flag is still off by inverse Angstrom, not by rounding + self.assertGreater(best, 1.0) + + +@requires_fiber +class TestFrames(unittest.TestCase): + """The frame rotations must follow the Vlieg convention.""" + + def test_q_phi_reproduces_hkl(self): + ub = _ub_calculator() + angles = HKLVlieg.VliegAngles(ub) + det = _geometry(45.0) + rows, cols = _pixel_grid(det, stride=401) + alpha_i = np.deg2rad(0.6) + omega, chi, phi = np.deg2rad(23.0), np.deg2rad(4.0), np.deg2rad(-7.0) + + q_phi = _q_from_conversion( + det, rows, cols, alpha_i, frame="Q_phi", omega=omega, chi=chi, phi=phi + ) + actual = np.einsum("ij,...j->...i", np.linalg.inv(ub.getUB()), q_phi) + + gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha_i) + h, k, ll = angles.anglesToHkl(alpha_i, delta, gamma, omega, chi, phi) + + np.testing.assert_allclose(actual[..., 0], h, atol=1e-9) + np.testing.assert_allclose(actual[..., 1], k, atol=1e-9) + np.testing.assert_allclose(actual[..., 2], ll, atol=1e-9) + + def test_q_cryst_is_b_times_hkl(self): + ub = _ub_calculator() + angles = HKLVlieg.VliegAngles(ub) + det = _geometry(45.0) + rows, cols = _pixel_grid(det, stride=401) + alpha_i = np.deg2rad(0.6) + omega, chi, phi = np.deg2rad(23.0), np.deg2rad(4.0), np.deg2rad(-7.0) + + actual = _q_from_conversion( + det, + rows, + cols, + alpha_i, + frame="Q_cryst", + omega=omega, + chi=chi, + phi=phi, + U=ub.getU(), + ) + + gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha_i) + hkl = np.stack( + angles.anglesToHkl(alpha_i, delta, gamma, omega, chi, phi), axis=-1 + ) + expected = np.einsum("ij,...j->...i", ub.lattice.B_mat, hkl) + + np.testing.assert_allclose(actual, expected, atol=1e-9) + + def test_frames_requiring_u_are_rejected_without_it(self): + for frame in qconversion.FRAMES_REQUIRING_U: + with self.subTest(frame=frame): + with self.assertRaises(ValueError): + qconversion.frameMatrix(frame, alpha=0.1, omega=0.2) + with self.assertRaises(ValueError): + qconversion.fiberUnits(frame) + + def test_every_frame_preserves_the_momentum_transfer(self): + det = _geometry(45.0) + rows, cols = _pixel_grid(det, stride=401) + alpha_i = np.deg2rad(0.6) + frame_angles = { + "omega": np.deg2rad(23.0), + "chi": np.deg2rad(4.0), + "phi": np.deg2rad(-7.0), + "U": _ub_calculator().getU(), + } + + reference = np.linalg.norm( + _q_from_conversion(det, rows, cols, alpha_i), axis=-1 + ) + for frame in qconversion.FRAMES: + with self.subTest(frame=frame): + q = _q_from_conversion( + det, rows, cols, alpha_i, frame=frame, **frame_angles + ) + np.testing.assert_allclose( + np.linalg.norm(q, axis=-1), reference, atol=1e-9 + ) + + def test_alpha_frame_matrix_is_the_identity(self): + np.testing.assert_allclose( + qconversion.frameMatrix("Q_alpha", alpha=0.3), np.identity(3), atol=1e-15 + ) + + def test_unknown_frame_is_rejected(self): + with self.assertRaises(ValueError): + qconversion.frameMatrix("Q_nonsense") + with self.assertRaises(ValueError): + qconversion.fiberUnits("Q_nonsense") + + def test_frame_metadata_is_consistent(self): + self.assertTrue( + set(qconversion.FRAMES_REQUIRING_OMEGA) <= set(qconversion.FRAMES) + ) + self.assertTrue(set(qconversion.FRAMES_REQUIRING_U) <= set(qconversion.FRAMES)) + # undoing U also requires the sample rotations to be undone first + self.assertTrue( + set(qconversion.FRAMES_REQUIRING_U) + <= set(qconversion.FRAMES_REQUIRING_OMEGA) + ) + self.assertNotIn(qconversion.DEFAULT_FRAME, qconversion.FRAMES_REQUIRING_OMEGA) + self.assertEqual(set(qconversion.FRAME_LABELS), set(qconversion.FRAMES)) + + +class TestInverseConversion(unittest.TestCase): + """The Q-plot readout inverts the projection, so it must round-trip. + + The forward conversion keeps the out-of-plane component and collapses the + two in-plane components into their radial distance. The inverse recovers + the discarded direction from the Ewald condition, so a coordinate taken + from the forward conversion has to lead back to the same detector angles. + + These tests do not need pyFAI's rebinning, only the geometry. + """ + + SAMPLE_ANGLES = dict( + omega=np.deg2rad(23.0), chi=np.deg2rad(4.0), phi=np.deg2rad(-7.0) + ) + + def _forward(self, frame, alpha_i, delta, gamma, ub): + """``(Q, q_ip, q_oop)`` of one detector angle pair, in the frame.""" + angles = HKLVlieg.VliegAngles(ub) + q_alpha = angles.QAlpha( + np.asarray(alpha_i), np.asarray(delta), np.asarray(gamma) + ) + rotation = qconversion.frameMatrix( + frame, alpha=alpha_i, U=ub.getU(), **self.SAMPLE_ANGLES + ) + Q = rotation @ q_alpha + # the signing convention of qIpOop + return Q, np.copysign(np.hypot(Q[0], Q[1]), Q[0]), Q[2] + + def test_round_trip_recovers_the_detector_angles(self): + ub = _ub_calculator() + k = ub.getK() + for frame in qconversion.FRAMES: + for alpha_deg in INCIDENCE_ANGLES: + for delta_deg, gamma_deg in ((12.0, 4.0), (-31.5, 22.0), (48.0, -3.0)): + with self.subTest( + frame=frame, alpha_i=alpha_deg, delta=delta_deg, gamma=gamma_deg + ): + alpha_i = np.deg2rad(alpha_deg) + delta, gamma = np.deg2rad([delta_deg, gamma_deg]) + Q, q_ip, q_oop = self._forward(frame, alpha_i, delta, gamma, ub) + + solutions = qconversion.qFromIpOop( + frame, + q_ip, + q_oop, + k, + alpha=alpha_i, + U=ub.getU(), + **self.SAMPLE_ANGLES, + ) + + self.assertTrue(solutions) + # the true momentum transfer is among the solutions + best = min(solutions, key=lambda s: np.linalg.norm(s - Q)) + np.testing.assert_allclose(best, Q, atol=1e-9) + + recovered = qconversion.detectorAngles( + best, + k, + frame=frame, + alpha=alpha_i, + U=ub.getU(), + **self.SAMPLE_ANGLES, + ) + np.testing.assert_allclose(recovered, (delta, gamma), atol=1e-9) + + def test_the_alpha_frame_is_never_ambiguous(self): + """The beam has no x component there, so the sign fixes the solution. + + This is why the default frame needs no disambiguation at all. + """ + ub = _ub_calculator() + rng = np.random.default_rng(20260803) + for _ in range(200): + alpha_i = np.deg2rad(rng.uniform(0.0, 3.0)) + delta = np.deg2rad(rng.uniform(-70.0, 70.0)) + gamma = np.deg2rad(rng.uniform(-20.0, 60.0)) + _, q_ip, q_oop = self._forward("Q_alpha", alpha_i, delta, gamma, ub) + + solutions = qconversion.qFromIpOop( + "Q_alpha", q_ip, q_oop, ub.getK(), alpha=alpha_i + ) + + self.assertEqual(len(solutions), 1) + + def test_solutions_keep_the_sign_of_the_in_plane_coordinate(self): + ub = _ub_calculator() + for frame in qconversion.FRAMES: + for delta_deg in (-40.0, -5.0, 5.0, 40.0): + with self.subTest(frame=frame, delta=delta_deg): + alpha_i = np.deg2rad(0.6) + delta = np.deg2rad(delta_deg) + gamma = np.deg2rad(14.0) + _, q_ip, q_oop = self._forward(frame, alpha_i, delta, gamma, ub) + + solutions = qconversion.qFromIpOop( + frame, + q_ip, + q_oop, + ub.getK(), + alpha=alpha_i, + U=ub.getU(), + **self.SAMPLE_ANGLES, + ) + + self.assertTrue(solutions) + for Q in solutions: + self.assertEqual(np.sign(Q[0]), np.sign(q_ip)) + + def test_unreachable_coordinates_have_no_solution(self): + """Beyond the Ewald sphere nothing can be reported.""" + ub = _ub_calculator() + + solutions = qconversion.qFromIpOop( + "Q_alpha", 10.0 * ub.getK(), 0.0, ub.getK(), alpha=np.deg2rad(0.6) + ) + + self.assertEqual(solutions, ()) + + def test_non_finite_coordinates_have_no_solution(self): + ub = _ub_calculator() + for q_ip, q_oop in ((np.nan, 0.0), (0.0, np.inf)): + with self.subTest(q_ip=q_ip, q_oop=q_oop): + self.assertEqual( + qconversion.qFromIpOop("Q_alpha", q_ip, q_oop, ub.getK()), () + ) + + def test_beam_direction_matches_the_collapsed_conversion(self): + """``beamDirection`` must be the ``c`` of the forward conversion.""" + for frame in qconversion.FRAMES: + for alpha_deg in INCIDENCE_ANGLES: + with self.subTest(frame=frame, alpha_i=alpha_deg): + ub = _ub_calculator() + alpha_i = np.deg2rad(alpha_deg) + _, c = qconversion.conversionCoefficients( + frame, + incident_angle=alpha_i, + tilt_angle=0.83, + U=ub.getU(), + **self.SAMPLE_ANGLES, + ) + + beam = qconversion.beamDirection( + frame, alpha=alpha_i, U=ub.getU(), **self.SAMPLE_ANGLES + ) + + np.testing.assert_allclose(beam, c, atol=1e-12) + + +@requires_fiber +class TestIntegrateImage(unittest.TestCase): + """The rebinned map places intensity where QAlpha predicts.""" + + def test_single_pixel_peak_lands_at_the_qalpha_position(self): + angles = HKLVlieg.VliegAngles(_ub_calculator()) + alpha_i = np.deg2rad(0.6) + + for azimuth_deg in (45.0, 233.7): + det = _geometry(azimuth_deg) + for row, col in ((900, 700), (600, 1000)): + with self.subTest(azimuth=azimuth_deg, pixel=(row, col)): + image = np.zeros(det.detector.shape, dtype=np.float64) + image[row, col] = 1000.0 + + result = qconversion.integrateImage(det, image, alpha_i) + + peak = np.unravel_index( + np.argmax(result.intensity), result.intensity.shape + ) + q_ip_binned = result.inplane[peak[1]] + q_oop_binned = result.outofplane[peak[0]] + + gamma, delta = det.surfaceAnglesPoint( + np.array([float(row)]), np.array([float(col)]), alpha_i + ) + expected = angles.QAlpha(alpha_i, delta, gamma) + expected_ip = np.hypot(expected[..., 0], expected[..., 1])[0] + expected_oop = expected[..., 2][0] + + # the rebinning quantises the position to one grid cell + bin_ip = ( + result.inplane.max() - result.inplane.min() + ) / result.inplane.size + bin_oop = ( + result.outofplane.max() - result.outofplane.min() + ) / result.outofplane.size + + self.assertLess(abs(abs(q_ip_binned) - expected_ip), bin_ip) + self.assertLess(abs(q_oop_binned - expected_oop), bin_oop) + + +@requires_fiber +class TestIntegratorCache(unittest.TestCase): + """The pyFAI integrator is reused, but never at the cost of stale results. + + A fresh integrator resets on every call, which is slow and floods the log. + Reusing one is only safe while the conversion is unchanged, because pyFAI + caches the per-pixel unit arrays under the unit name. + """ + + def setUp(self): + self.detector = _geometry(45.0) + self.image = np.zeros(self.detector.detector.shape, dtype=np.float64) + self.image[900, 700] = 1000.0 + self.alpha = np.deg2rad(0.6) + qconversion._INTEGRATOR.clear() + self.addCleanup(qconversion._INTEGRATOR.clear) + self.addCleanup(qconversion._CONVERSION.clear) + + def _integrate(self, **keywords): + qconversion._CONVERSION.clear() + return qconversion.integrateImage( + self.detector, self.image, self.alpha, npt=200, **keywords + ) + + def test_integrator_is_reused_when_nothing_changed(self): + self._integrate() + first = qconversion._INTEGRATOR._integrator + self._integrate() + + self.assertIs(qconversion._INTEGRATOR._integrator, first) + + def test_integrator_is_rebuilt_when_the_conversion_changes(self): + self._integrate(frame="Q_alpha") + first = qconversion._INTEGRATOR._integrator + self._integrate(frame="Q_chi", omega=np.deg2rad(10.0), chi=np.deg2rad(5.0)) + + self.assertIsNot(qconversion._INTEGRATOR._integrator, first) + + def test_reuse_does_not_serve_stale_coordinates(self): + """chi tilts the surface normal, so the axes really have to change.""" + + def peak(chi_deg): + result = self._integrate( + frame="Q_chi", omega=np.deg2rad(10.0), chi=np.deg2rad(chi_deg) + ) + index = np.unravel_index( + np.argmax(result.intensity), result.intensity.shape + ) + return result.inplane[index[1]], result.outofplane[index[0]] + + straight = peak(0.0) + tilted = peak(9.0) + again = peak(0.0) + + self.assertNotAlmostEqual(straight[1], tilted[1], places=3) + np.testing.assert_allclose(straight, again, atol=1e-9) + + +class TestKernelBackends(unittest.TestCase): + """The compiled kernel and the numpy fallback must agree.""" + + def setUp(self): + det = _geometry(45.0) + rows, cols = _pixel_grid(det, stride=53) + self.arguments = ( + "Q_alpha", + *_pixel_positions(det, rows, cols), + det.wavelength, + ) + self.keywords = { + "incident_angle": np.deg2rad(0.6), + "tilt_angle": qconversion.tiltAngleFromAzimuth(np.deg2rad(45.0)), + } + self.addCleanup(qconversion._CONVERSION.clear) + + def _q_ip_oop(self): + qconversion._CONVERSION.clear() + return qconversion.qIpOop(*self.arguments, **self.keywords) + + def test_q_ip_oop_matches_the_reference_components(self): + q_ip, q_oop = self._q_ip_oop() + qx, qy, qz = qconversion.qComponents(*self.arguments, **self.keywords) + + np.testing.assert_allclose(q_ip, np.copysign(np.hypot(qx, qy), qx), atol=1e-9) + np.testing.assert_allclose(q_oop, qz, atol=1e-9) + + def test_compiled_kernel_matches_the_numpy_fallback(self): + if not qconversion.hasKernel(): + self.skipTest("the compiled conversion kernel is not available") + + compiled_ip, compiled_oop = self._q_ip_oop() + + kernel = qconversion._KERNEL + qconversion._KERNEL = None # force the numpy fallback + try: + fallback_ip, fallback_oop = self._q_ip_oop() + finally: + qconversion._KERNEL = kernel + + np.testing.assert_allclose(compiled_ip, fallback_ip, atol=1e-9) + np.testing.assert_allclose(compiled_oop, fallback_oop, atol=1e-9) + + def test_repeated_call_is_served_from_the_cache(self): + first = self._q_ip_oop() + second = qconversion.qIpOop(*self.arguments, **self.keywords) + # pyFAI evaluates the two units separately; the image is converted once + self.assertIs(first[0], second[0]) + self.assertIs(first[1], second[1]) + + +if __name__ == "__main__": + unittest.main() diff --git a/orgui/app/test/test_q_plot_orientation.py b/orgui/app/test/test_q_plot_orientation.py new file mode 100644 index 0000000..a74218d --- /dev/null +++ b/orgui/app/test/test_q_plot_orientation.py @@ -0,0 +1,448 @@ +"""The Q-plot must be wired to the validated conversion in :mod:`qconversion`. + +The numerical equivalence of the Q-plot with orGUI's own QAlpha calculation is +established in ``orgui/app/test/test_q_conversion.py``. These tests only cover +the GUI wiring: which frames are offered, which one is preselected, and how the +selection is read back. +""" + +import os + +import numpy as np +import pyFAI +import pytest + +from orgui.app import qconversion +from orgui.datautils.xrayutils import DetectorCalibration, HKLVlieg + +orGUI = pytest.importorskip("orgui.app.orGUI", reason="Qt bindings are required") + + +class _StubComboBox: + """Minimal stand-in for the frame selector combo box.""" + + def __init__(self, data): + self._data = data + + def currentData(self): + return self._data + + +class _StubOrGUI: + """The frame accessor only reads the combo box.""" + + def __init__(self, data): + self.qFrameSelector = _StubComboBox(data) + + +@pytest.mark.parametrize("frame", qconversion.FRAMES) +def test_selected_frame_is_read_from_the_combo_box(frame): + assert orGUI.orGUI._selectedQFrame(_StubOrGUI(frame)) == frame + + +@pytest.mark.parametrize("data", [None, "Q_nonsense", ""]) +def test_unknown_selection_falls_back_to_the_default_frame(data): + selected = orGUI.orGUI._selectedQFrame(_StubOrGUI(data)) + assert selected == qconversion.DEFAULT_FRAME + + +def test_default_frame_is_the_alpha_frame(): + """QAlpha is what the rest of orGUI uses, so it is the natural default.""" + assert qconversion.DEFAULT_FRAME == "Q_alpha" + assert qconversion.DEFAULT_FRAME in qconversion.FRAMES + + +def test_every_frame_has_a_label_for_the_selector(): + for frame in qconversion.FRAMES: + assert qconversion.FRAME_LABELS[frame] + + +def test_fiber_integrator_flag_is_taken_from_qconversion(): + """The action is only usable when FiberIntegrator can be imported.""" + assert orGUI.HAS_FIBER_INTEGRATOR == qconversion.HAS_FIBER_INTEGRATOR + + +@pytest.mark.parametrize("suffix", [".svg", ".png"]) +def test_the_toolbar_icon_ships_with_the_package(suffix): + """The action loads its icon by name, which silently yields a blank icon. + + ``QIcon`` never reports a missing file, so a packaging mistake would only + show up as an empty toolbar button. + """ + from orgui import resources + + assert os.path.isfile(resources.getPath("q-plot" + suffix)) + + +class _StubAction: + def __init__(self, checked): + self._checked = checked + + def isChecked(self): + return self._checked + + +class _RefreshStub: + """Records the toggles the refresh performs.""" + + def __init__(self, checked, reentrant=False): + self.plotAgainstQAct = _StubAction(checked) + self._qPlotRefreshing = False + self._reentrant = reentrant + self.calls = [] + + def _convertImagetoQ(self, value): + self.calls.append(value) + if self._reentrant: + # emulate a signal arriving while the plot is being rebuilt + orGUI.orGUI._refreshQPlot(self) + + +def test_refresh_does_nothing_while_the_q_plot_is_off(): + stub = _RefreshStub(checked=False) + + orGUI.orGUI._refreshQPlot(stub) + + assert stub.calls == [] + + +def test_refresh_rebuilds_in_place_and_stays_in_q_mode(): + """Changing the image must not drop out of the reciprocal-space view.""" + stub = _RefreshStub(checked=True) + + orGUI.orGUI._refreshQPlot(stub) + + # a single conversion replaces the previous one; the action stays checked + assert stub.calls == [True] + assert stub.plotAgainstQAct.isChecked() + + +def test_refresh_accepts_the_replot_signal_argument(): + """sigReplotRequest carries a bool, which the slot must tolerate.""" + stub = _RefreshStub(checked=True) + + orGUI.orGUI._refreshQPlot(stub, True) + + assert stub.calls == [True] + + +def test_refresh_is_not_reentrant(): + stub = _RefreshStub(checked=True, reentrant=True) + + orGUI.orGUI._refreshQPlot(stub) + + assert stub.calls == [True] + assert stub._qPlotRefreshing is False + + +def test_suspended_q_plot_batches_into_one_rebuild(): + """Untoggling max/sum before replotting must not rebuild several times.""" + stub = _RefreshStub(checked=True) + + with orGUI.orGUI._suspendedQPlot(stub): + orGUI.orGUI._refreshQPlot(stub) + orGUI.orGUI._refreshQPlot(stub) + assert stub.calls == [] + orGUI.orGUI._refreshQPlot(stub) + + assert stub.calls == [True] + assert stub._qPlotRefreshing is False + + +# --- the position readout while the Q-plot is on ------------------------- +# +# The plot axes are momentum transfer there, not pixels, so the readout has to +# invert the conversion. These tests drive the real methods against a fully +# specified geometry. + +ENERGY = 78.0 # keV +ALPHA = np.deg2rad(0.6) +OMEGA = np.deg2rad(23.0) + + +def _geometry(azimuth_deg=45.0): + det = DetectorCalibration.Detector2D_SXRD() + det.detector = pyFAI.detector_factory("Pilatus2m") + det.wavelength = (12.39842 / ENERGY) * 1e-10 + det.dist = 0.729 + det.poni1 = 0.21 + det.poni2 = 0.14 + det.setAzimuthalReference(np.deg2rad(azimuth_deg)) + return det + + +class _StubUBCalculatorWidget: + def __init__(self, det): + lattice = HKLVlieg.Lattice( + np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0]) + ) + self.ubCal = HKLVlieg.UBCalculator(lattice, ENERGY) + self.ubCal.defaultU() + self.angles = HKLVlieg.VliegAngles(self.ubCal) + self.detectorCal = det + self.chi = 0.0 + self.phi = 0.0 + self.n = 1.0 + + +class _ReadoutStub: + """Everything the Q-plot readout inversion reads.""" + + def __init__(self, frame="Q_alpha"): + self._frame = frame + self.ubcalc = _StubUBCalculatorWidget(_geometry()) + + def _selectedQFrame(self): + return self._frame + + def _qPlotAngles(self): + return ALPHA, OMEGA + + _qPlotOrientationMatrix = orGUI.orGUI._qPlotOrientationMatrix + _qPlotToPosition = orGUI.orGUI._qPlotToPosition + _isOnDetector = orGUI.orGUI._isOnDetector + + +def _plot_coordinates(stub, row, col): + """Q-plot coordinates of one detector pixel, in inverse Angstrom.""" + det = stub.ubcalc.detectorCal + gamma, delta = det.surfaceAnglesPoint( + np.array([float(row)]), np.array([float(col)]), ALPHA + ) + q_alpha = stub.ubcalc.angles.QAlpha(ALPHA, delta, gamma)[0] + rotation = qconversion.frameMatrix( + stub._selectedQFrame(), + alpha=ALPHA, + omega=OMEGA, + chi=stub.ubcalc.chi, + phi=stub.ubcalc.phi, + U=stub.ubcalc.ubCal.getU(), + ) + Q = rotation @ q_alpha + return np.copysign(np.hypot(Q[0], Q[1]), Q[0]), Q[2], delta[0], gamma[0] + + +@pytest.mark.parametrize("frame", qconversion.FRAMES) +@pytest.mark.parametrize("pixel", [(900, 700), (600, 1000), (150, 1400)]) +def test_readout_inverts_the_conversion_back_to_the_same_pixel(frame, pixel): + """Hovering a Q-plot point must report the angles of the pixel it came from.""" + stub = _ReadoutStub(frame) + q_ip, q_oop, delta, gamma = _plot_coordinates(stub, *pixel) + + resolved = orGUI.orGUI._qPlotToPosition(stub, q_ip, q_oop) + + assert resolved is not None + alpha, omega, recovered_delta, recovered_gamma, _ = resolved + assert (alpha, omega) == (ALPHA, OMEGA) + np.testing.assert_allclose( + [recovered_delta, recovered_gamma], [delta, gamma], atol=1e-9 + ) + + +@pytest.mark.parametrize("frame", qconversion.FRAMES) +@pytest.mark.parametrize("pixel", [(900, 700), (600, 1000), (150, 1400)]) +def test_reported_q_reproduces_the_cursor_coordinates(frame, pixel): + """The reported Q must be the point the cursor is on, not a corrected one. + + The Q-plot axes carry no refraction correction, so applying one to the + reported momentum transfer would contradict the displayed position. + """ + stub = _ReadoutStub(frame) + q_ip, q_oop, _, _ = _plot_coordinates(stub, *pixel) + + _, _, _, _, Q = orGUI.orGUI._qPlotToPosition(stub, q_ip, q_oop) + + np.testing.assert_allclose(np.copysign(np.hypot(Q[0], Q[1]), Q[0]), q_ip, atol=1e-9) + np.testing.assert_allclose(Q[2], q_oop, atol=1e-9) + + +def test_readout_reports_nothing_outside_the_detector(): + """The rebinned grid is rectangular and reaches past the detector.""" + stub = _ReadoutStub() + + assert orGUI.orGUI._qPlotToPosition(stub, 0.0, -5.0) is None + + +def test_readout_reports_nothing_beyond_the_ewald_sphere(): + stub = _ReadoutStub() + k = stub.ubcalc.ubCal.getK() + + assert orGUI.orGUI._qPlotToPosition(stub, 10.0 * k, 0.0) is None + + +# --- which way the Q-plot ordinate runs ---------------------------------- +# +# The detector image is drawn with the row index increasing downwards. The +# Q-plot should keep that visual orientation, which means its ordinate has to +# follow the geometry rather than always pointing upwards. + + +def _out_of_plane_at_top_and_bottom(det, alpha): + """``q_perp`` at the top-centre and bottom-centre pixel, in A^-1.""" + lattice = HKLVlieg.Lattice(np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0])) + ub = HKLVlieg.UBCalculator(lattice, ENERGY) + ub.defaultU() + shape = det.detector.shape + rows = np.array([0.0, float(shape[0] - 1)]) + cols = np.full(2, 0.5 * (shape[1] - 1)) + gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha) + return HKLVlieg.VliegAngles(ub).QAlpha(alpha, delta, gamma)[..., 2] + + +@pytest.mark.parametrize("azimuth_deg", [0.0, 45.0, 90.0, 135.0, 180.0, 270.0]) +@pytest.mark.parametrize("alpha_deg", [0.6, -0.6]) +def test_ordinate_direction_matches_the_detector_image(azimuth_deg, alpha_deg): + """The flag must say where q_perp really grows on the displayed image.""" + det = _geometry(azimuth_deg) + alpha = np.deg2rad(alpha_deg) + top, bottom = _out_of_plane_at_top_and_bottom(det, alpha) + + flipped = qconversion.outOfPlaneIncreasesWithRow(det, alpha) + + assert flipped == bool(bottom > top) + + +def test_upward_scattering_keeps_the_ordinate_pointing_up(): + """The ordinary geometry must not be inverted.""" + det = _geometry(90.0) + + assert not qconversion.outOfPlaneIncreasesWithRow(det, np.deg2rad(0.6)) + + +def test_inverted_geometry_inverts_the_ordinate(): + """ID31's downward scattering setup, reached by rotating the azimuth.""" + det = _geometry(270.0) + + assert qconversion.outOfPlaneIncreasesWithRow(det, np.deg2rad(0.6)) + + +@pytest.mark.parametrize("frame", qconversion.FRAMES) +def test_ordinate_direction_is_defined_for_every_frame(frame): + lattice = HKLVlieg.Lattice(np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0])) + ub = HKLVlieg.UBCalculator(lattice, ENERGY) + ub.defaultU() + det = _geometry(270.0) + + flipped = qconversion.outOfPlaneIncreasesWithRow( + det, + np.deg2rad(0.6), + frame=frame, + omega=OMEGA, + chi=0.05, + phi=-0.1, + U=ub.getU(), + ) + + assert isinstance(flipped, bool) + + +# --- restoring the axis orientation when the Q-plot is switched off ------- + + +class _AxisPlotStub: + """Records the axis orientation the Q-plot toggle applies.""" + + def __init__(self, xInverted, yInverted): + self.xInverted = xInverted + self.yInverted = yInverted + self.zoomResets = 0 + + def isXAxisInverted(self): + return self.xInverted + + def isYAxisInverted(self): + return self.yInverted + + def setXAxisInverted(self, flag): + self.xInverted = flag + + def setYAxisInverted(self, flag): + self.yInverted = flag + + def resetZoom(self): + self.zoomResets += 1 + + # not exercised by the restore path + def setQPlotActive(self, active): + pass + + def getAllImages(self): + return [] + + def setActiveImage(self, legend): + pass + + +class _RestoreStub: + """The part of orGUI the axis restore touches.""" + + def __init__(self, xInverted, yInverted): + self.centralPlot = _AxisPlotStub(xInverted, yInverted) + self._qPlotPreviousXInverted = None + self._qPlotPreviousYInverted = None + self.currentAddImageLabel = None + self.currentImageLabel = "scan_image" + self.allimgmax = None + self.allimgsum = None + + def _enter(self, flipOrdinate=False): + """The axis part of switching the Q-plot on.""" + if self._qPlotPreviousYInverted is None: + self._qPlotPreviousYInverted = self.centralPlot.isYAxisInverted() + if self._qPlotPreviousXInverted is None: + self._qPlotPreviousXInverted = self.centralPlot.isXAxisInverted() + self.centralPlot.setYAxisInverted(flipOrdinate) + self.centralPlot.setXAxisInverted(False) + + +@pytest.mark.parametrize("xInverted", [False, True]) +@pytest.mark.parametrize("yInverted", [False, True]) +def test_leaving_the_q_plot_restores_both_axes(xInverted, yInverted): + """Both origins must come back, whatever they were set to.""" + stub = _RestoreStub(xInverted, yInverted) + stub._enter() + + orGUI.orGUI._convertImagetoQ(stub, False) + + assert stub.centralPlot.xInverted == xInverted + assert stub.centralPlot.yInverted == yInverted + assert stub._qPlotPreviousXInverted is None + assert stub._qPlotPreviousYInverted is None + + +def test_restore_does_not_depend_on_the_max_sum_image(): + """A conversion that failed after flipping the axes must still roll back. + + The restore used to sit inside the block that rebuilds the maximum or sum + image, so it was skipped whenever no additional image was displayed. + """ + stub = _RestoreStub(xInverted=True, yInverted=True) + stub._enter(flipOrdinate=True) + assert stub.currentAddImageLabel is None # nothing was ever added + + orGUI.orGUI._convertImagetoQ(stub, False) + + assert stub.centralPlot.xInverted is True + assert stub.centralPlot.yInverted is True + assert stub.centralPlot.zoomResets == 1 + + +def test_leaving_without_ever_entering_touches_nothing(): + stub = _RestoreStub(xInverted=True, yInverted=False) + + orGUI.orGUI._convertImagetoQ(stub, False) + + assert stub.centralPlot.xInverted is True + assert stub.centralPlot.yInverted is False + assert stub.centralPlot.zoomResets == 0 + + +def test_rebuilding_in_place_keeps_the_recorded_orientation(): + """Stepping through images re-enters the Q-plot; the record must survive.""" + stub = _RestoreStub(xInverted=False, yInverted=True) + stub._enter() + stub._enter() # a rebuild, not a fresh entry + + orGUI.orGUI._convertImagetoQ(stub, False) + + assert stub.centralPlot.yInverted is True diff --git a/orgui/app/test/test_segmented_scan_entries.py b/orgui/app/test/test_segmented_scan_entries.py new file mode 100644 index 0000000..f6b9f72 --- /dev/null +++ b/orgui/app/test/test_segmented_scan_entries.py @@ -0,0 +1,242 @@ +"""The segmented scan loader asks the backend, it does not guess. + +Two routes lead to the same answer: a backend implementing ``listScans`` is +asked directly, and any other backend has its own ``parse_h5_node`` applied to +every entry of the file root. Neither route may disagree with how a single scan +is opened, because the identifiers produced here are handed to ``openScan``. +""" + +import pytest + +from orgui.backend import scans + +orGUI = pytest.importorskip("orgui.app.orGUI", reason="Qt bindings are required") + + +# --- fast path ----------------------------------------------------------- + + +class _ListingBackend: + """A backend that implements the fast path.""" + + returns = [] + + @classmethod + def listScans(cls, h5file): + del h5file + return cls.returns + + +def _fastPath(returns): + backend = type("_Backend", (_ListingBackend,), {"returns": returns}) + return orGUI.orGUI._scanEntriesFromBackend(backend, object()) + + +def test_bare_numbers_are_accepted(): + """The simplest thing a backend can return.""" + entries = _fastPath([1, 2, 10]) + + assert [ddict["scanno"] for ddict, _ in entries] == [1, 2, 10] + assert [label for _, label in entries] == [None, None, None] + + +def test_bare_names_are_accepted(): + """A backend may be addressed by a string rather than by a number.""" + entries = _fastPath(["ascan_12", "dscan_3"]) + + assert [ddict["scanno"] for ddict, _ in entries] == ["ascan_12", "dscan_3"] + + +def test_an_identifier_is_reported_under_both_keys(): + """openScan reads 'name' for some beamtimes and 'scanno' for others.""" + ((ddict, _label),) = _fastPath(["ascan_12"]) + + assert ddict["scanno"] == "ascan_12" + assert ddict["name"] == "ascan_12" + + +def test_pairs_carry_a_label_for_the_dialog(): + entries = _fastPath([(1, "ascan th 0 90 90 1"), (2, "dscan mu 0 1 20 1")]) + + assert [ddict["scanno"] for ddict, _ in entries] == [1, 2] + assert [label for _, label in entries] == [ + "ascan th 0 90 90 1", + "dscan mu 0 1 20 1", + ] + + +def test_the_legacy_id31_listing_fills_the_key_openscan_reads(): + """openScan opens ch5523 with ddict['name'], not with a scan number. + + The segmented loader never set that key, so this beamtime raised KeyError + before it could open anything. Its listing reports group names, which the + normalisation has to put where openScan looks. + """ + from orgui.backend.beamline import id31_tools + + class _LegacyFile: + names = {"ascan_12": b"ascan th 0 90 90 1"} + + def __iter__(self): + return iter(self.names) + + def __getitem__(self, path): + name, _, field = path.partition("/") + if name not in self.names or field != "title": + raise KeyError(path) + return type("_D", (), {"__getitem__": lambda s, k: self.names[name]})() + + listed = id31_tools.BlissScan.listScans(_LegacyFile()) + backend = type("_Backend", (_ListingBackend,), {"returns": listed}) + + ((ddict, label),) = orGUI.orGUI._scanEntriesFromBackend(backend, object()) + + assert ddict["name"] == "ascan_12" + assert label == "ascan th 0 90 90 1" + + +def test_a_backend_without_the_fast_path_asks_for_the_fallback(): + """None is the signal to fall back, and must not be confused with [].""" + + class _NoFastPath(scans.Scan): + def __init__(self): + pass + + @classmethod + def parse_h5_node(cls, node): + return {} + + def __len__(self): + return 0 + + def get_raw_img(self, i): + return None + + assert orGUI.orGUI._scanEntriesFromBackend(_NoFastPath, object()) is None + assert _fastPath([]) == [] + + +# --- fallback ------------------------------------------------------------ + + +class _StubItem: + def __init__(self, name): + self.name = name + + +class _StubModel: + """Enough of the silx HDF5 tree model to walk the root's children.""" + + def __init__(self, names): + self._names = list(names) + + def rowCount(self, parent): + del parent + return len(self._names) + + def index(self, row, column, parent): + del column, parent + return row + + def data(self, index, role): + del role + return _StubItem(self._names[index]) + + +class _StubNode: + """What silx H5Node reports for a direct child of the file root.""" + + def __init__(self, item): + self.basename = item.name + self.local_name = "/" + item.name + + +@pytest.fixture +def h5node(monkeypatch): + """Wrap tree items in the stub node instead of the real silx one.""" + monkeypatch.setattr(orGUI.silx.gui.hdf5, "H5Node", _StubNode) + + +def _fallback(backend, names): + return orGUI.orGUI._scanEntriesFromNodes(backend, _StubModel(names), None) + + +class _BlissLikeBackend: + """Numbered by the integer in front of the dot, as the ID31 backends are.""" + + @classmethod + def parse_h5_node(cls, node): + scanno, _subscan = node.basename.split(".") + return {"scanno": int(scanno), "name": node.local_name} + + +class _NamedBackend: + """Addressed by name, as the legacy ch5523 backend is.""" + + @classmethod + def parse_h5_node(cls, node): + return { + "scanno": int(node.local_name.split("_")[-1]), + "name": node.local_name.strip("/"), + } + + +def test_fallback_numbers_scans_the_way_the_backend_does(h5node): + entries = _fallback(_BlissLikeBackend, ["1.1", "2.1", "10.1"]) + + assert [ddict["scanno"] for ddict, _ in entries] == [1, 2, 10] + + +def test_fallback_collapses_the_subscans_of_one_measurement(h5node): + """1.1 and 1.2 are the same scan, and must appear once.""" + entries = _fallback(_BlissLikeBackend, ["1.1", "1.2", "2.1", "2.2"]) + + assert [ddict["scanno"] for ddict, _ in entries] == [1, 2] + # the first occurrence wins, which is the fast counter subscan + assert entries[0][0]["name"] == "/1.1" + + +@pytest.mark.parametrize("name", ["instrument", "title", "start_time", "entry"]) +def test_fallback_skips_entries_the_backend_cannot_parse(h5node, name): + """Raising from parse_h5_node is how a backend says "not a scan".""" + entries = _fallback(_BlissLikeBackend, [name, "3.1"]) + + assert [ddict["scanno"] for ddict, _ in entries] == [3] + + +def test_fallback_skips_a_backend_that_reports_no_scan_number(h5node): + """InterlacedScan and ImportImagesScan return an empty dict.""" + + class _NoScanNo: + @classmethod + def parse_h5_node(cls, node): + return {} + + assert _fallback(_NoScanNo, ["1.1", "2.1"]) == [] + + +def test_fallback_keeps_the_name_a_backend_is_opened_with(h5node): + """ch5523 is opened by name; that name has to survive to openScan.""" + entries = _fallback(_NamedBackend, ["ascan_12", "dscan_3"]) + + assert [ddict["name"] for ddict, _ in entries] == ["ascan_12", "dscan_3"] + assert [ddict["scanno"] for ddict, _ in entries] == [12, 3] + + +def test_fallback_reports_no_label(h5node): + """parse_h5_node has no title to offer, so the dialog shows the name.""" + entries = _fallback(_BlissLikeBackend, ["1.1"]) + + assert entries[0][1] is None + + +def test_both_routes_agree_on_the_same_file(h5node): + """The whole point: the fast path may not disagree with parse_h5_node.""" + names = ["1.1", "1.2", "2.1", "2.2", "10.1"] + + fallback = _fallback(_BlissLikeBackend, names) + fast = _fastPath([1, 2, 10]) + + assert [ddict["scanno"] for ddict, _ in fallback] == [ + ddict["scanno"] for ddict, _ in fast + ] diff --git a/orgui/backend/backends.py b/orgui/backend/backends.py index 99ce782..8517d12 100644 --- a/orgui/backend/backends.py +++ b/orgui/backend/backends.py @@ -47,7 +47,7 @@ # the box in the lower left corner in the gui and manually selecting # a backend -default_beamtime = "id31_default" +default_beamtime = "id31_default_p4" beamtimes = { "ch5523": (datetime(2018, 9, 22), datetime(2018, 10, 5)), @@ -93,6 +93,7 @@ def getBeamtimeId(dt): "P212_default": H5Fastsweep, "ch7149": BlissScan_EBS, "ch7856": BlissScan_EBS, + "ch8153": BlissScan_EBS_p4, "id31_default": BlissScan_EBS, "id31_default_p4": BlissScan_EBS_p4, } @@ -197,6 +198,16 @@ def openScan(btid, ddict): else: fscan = fscancls(ddict["file"], ddict["scanno"], loadimg=False, muoffset=0) + elif btid == "ch8153": + if "node" in ddict: + fscan = fscancls( + ddict["node"], ddict["scanno"], loadimg=False, muoffset=-2.31 + ) + else: + fscan = fscancls( + ddict["file"], ddict["scanno"], loadimg=False, muoffset=-2.31 + ) + else: try: fscan = fscancls(ddict["file"], ddict["scanno"]) diff --git a/orgui/backend/beamline/id31_tools.py b/orgui/backend/beamline/id31_tools.py index 4a50cb5..4cb4d62 100644 --- a/orgui/backend/beamline/id31_tools.py +++ b/orgui/backend/beamline/id31_tools.py @@ -40,6 +40,7 @@ from silx.io import dictdump from silx.io import specfile import silx.io +import silx.io.h5py_utils # noqa: F401 (registers silx.io.h5py_utils for attribute access) from .ID31DiffractLinTilt import ID31DiffractLinTilt from ..scans import Scan @@ -120,6 +121,17 @@ def __init__(self, filename): self.img = np.rot90(self.img, 3) +def _scanTitle(h5file, name): + """Title of one scan, falling back to its name when there is none.""" + try: + title = h5file[name + "/title"][()] + except Exception: + return name + if isinstance(title, bytes): + return title.decode() + return str(title) + + # currently only set up for th scans in TOMO session, cbf fileformat class Fastscan(Scan): def __init__(self, fastscan_specfile, scanno): @@ -219,6 +231,29 @@ def __getitem__(self, key): def __len__(self): return self.nopoints + @classmethod + def listScans(cls, h5file): + """List the scans of a BLISS ``"."`` file. + + Only the ``.1`` subscan is a selectable scan, the others belong to the + same measurement. The scan number is the integer in front of the dot, + the same one :meth:`parse_h5_node` reports, and the label is the scan + title stored in the file. + """ + entries = [] + for name in h5file: + name = str(name) + scan, dot, subscan = name.rpartition(".") + if not dot or subscan != "1": + continue + try: + scanno = int(scan.split("_")[-1]) + except ValueError: + continue + entries.append((scanno, _scanTitle(h5file, name))) + entries.sort(key=lambda entry: entry[0]) + return entries + class BlissScan_EBS(Fastscan): def __init__( @@ -247,7 +282,7 @@ def __init__( self.filename_base = ( filename_noext # filename_noext[:filename_noext.rfind('_')] ) - with silx.io.open(hdffilepath_orNode) as f: + with silx.io.h5py_utils.File(hdffilepath_orNode) as f: # print([d for d in f]) for d in f: scansuffix = d.split("_")[-1] @@ -635,7 +670,7 @@ def __init__( self.filename_base = ( filename_noext # filename_noext[:filename_noext.rfind('_')] ) - with silx.io.open(hdffilepath_orNode) as f: + with silx.io.h5py_utils.File(hdffilepath_orNode) as f: # print([d for d in f]) for d in f: scansuffix = d.split("_")[-1] @@ -1137,6 +1172,32 @@ def parse_h5_node(cls, obj): ddict["name"] = obj.local_name.strip("/") return ddict + @classmethod + def listScans(cls, h5file): + """List the scans of a pre-EBS BLISS file, named ``"_"``. + + This overrides the ``"."`` implementation of + :class:`Fastscan`, which does not apply here: the scans of this legacy + format are named after the command that recorded them, ``ascan_12``, + ``dscan_3`` and so on, and there are no subscans. + + The identifier is the group name rather than the scan number, because + that is what ``__init__`` takes and what + :func:`~orgui.backend.backends.openScan` passes it for this beamtime. + The entries are ordered by the trailing scan number, the way + :meth:`parse_h5_node` reads it. + """ + entries = [] + for name in h5file: + name = str(name) + try: + scanno = int(name.split("_")[-1]) + except ValueError: + continue + entries.append((scanno, name)) + entries.sort(key=lambda entry: entry[0]) + return [(name, _scanTitle(h5file, name)) for _scanno, name in entries] + @property def auxillary_counters(self): """Optional: provide a list of counters or motor names, that should diff --git a/orgui/backend/scans.py b/orgui/backend/scans.py index 1c0d715..c746890 100644 --- a/orgui/backend/scans.py +++ b/orgui/backend/scans.py @@ -95,6 +95,44 @@ def auxillary_counters(self): def parse_h5_node(cls, node): pass + @classmethod + def listScans(cls, h5file): + """Optional: the scans an opened HDF5 file contains. + + Used by the segmented ("interlaced") scan loader to fill its selection + dialog, which has to know about every scan in the file before any scan + object exists. + + Return whatever your ``__init__`` accepts as its scan identifier, which + is the same thing :meth:`parse_h5_node` reports for a single node: + + .. code-block:: python + + return [1, 2, 10] # scan numbers + return ["ascan_12", "dscan_3"] # scan names + + Optionally pair each identifier with a label shown in the dialog: + + .. code-block:: python + + return [(1, "ascan th 0 90 90 1"), (2, "dscan mu 0 1 20 1")] + + Do not use ``float`` for BLISS style ``"."`` names: + ``1.1`` and ``1.10`` are the same float. Integers, strings and numpy + integers are all fine, the identifier is passed on untouched and is + only converted to ``str`` for display. + + Implementing this is optional. Without it orGUI calls + :meth:`parse_h5_node` on every entry of the file root instead, which + needs no extra code in the backend but is slower, and which relies on + ``parse_h5_node`` raising for entries that are not scans. + + :param h5file: an open h5py-like file object. + :returns: identifiers, or ``(identifier, label)`` pairs. + :rtype: list + """ + raise NotImplementedError(f"{cls.__name__} does not implement listScans") + @abstractmethod def __len__(self): """returns the number of entries in the scan.""" diff --git a/orgui/backend/test/test_scan_listing.py b/orgui/backend/test/test_scan_listing.py new file mode 100644 index 0000000..5c7137f --- /dev/null +++ b/orgui/backend/test/test_scan_listing.py @@ -0,0 +1,241 @@ +"""Backends describe their own file layout to the segmented scan loader. + +The loader used to guess the layout from the backend's class name, which +silently mishandled every backend not called ``BlissScan*`` -- including the +example backend shipped with orGUI. A backend now either implements +``Scan.listScans``, or lets the loader fall back to its own ``parse_h5_node``. +""" + +import importlib.util +import pathlib + +import pytest + +from orgui.backend import scans +from orgui.backend.beamline import id31_tools + + +class _StubDataset: + def __init__(self, value): + self._value = value + + def __getitem__(self, key): + assert key == () + return self._value + + +class _StubH5File: + """Just enough of an h5py file: iteration over names and item access.""" + + def __init__(self, entries): + self._entries = dict(entries) + + def __iter__(self): + return iter(self._entries) + + def __getitem__(self, path): + name, _, field = path.partition("/") + if name not in self._entries: + raise KeyError(path) + title = self._entries[name] + if field != "title" or title is None: + raise KeyError(path) + return _StubDataset(title) + + +BLISS_FILE = _StubH5File( + [ + ("2.1", b"ascan th 0 90 90 1"), + ("2.2", b"subscan of 2"), + ("1.1", b"ascan th 0 45 45 1"), + ("1.2", b"subscan of 1"), + ("10.1", b"ascan th 0 10 10 1"), + ] +) + +BLISS_BACKENDS = (id31_tools.BlissScan_EBS, id31_tools.BlissScan_EBS_p4) + + +def test_the_base_class_declares_the_fast_path_as_unavailable(): + """The template declares listScans; not implementing it is legal.""" + with pytest.raises(NotImplementedError): + scans.Scan.listScans(BLISS_FILE) + + +def test_the_declaration_is_not_abstract(): + """An @abstractmethod would break every existing backend on construction.""" + assert "listScans" not in scans.Scan.__abstractmethods__ + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +def test_only_the_fast_counter_subscan_is_listed(backend): + entries = backend.listScans(BLISS_FILE) + + assert [scanno for scanno, _ in entries] == [1, 2, 10] + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +def test_scans_are_listed_in_scan_order(backend): + """The dialog lists them by scan number, not in file order.""" + entries = backend.listScans(BLISS_FILE) + + assert entries == sorted(entries) + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +def test_the_label_is_the_scan_title(backend): + entries = dict(backend.listScans(BLISS_FILE)) + + assert entries[1] == "ascan th 0 45 45 1" + assert entries[2] == "ascan th 0 90 90 1" + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +def test_a_double_digit_subscan_is_not_the_fast_counter_subscan(backend): + """``1.10`` must not be mistaken for ``1.1``.""" + entries = backend.listScans( + _StubH5File([("1.1", b"real"), ("1.10", b"other subscan")]) + ) + + assert [scanno for scanno, _ in entries] == [1] + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +@pytest.mark.parametrize("name", ["notascan", "instrument", "1.x", "x.1"]) +def test_entries_that_are_not_scans_are_skipped(backend, name): + entries = backend.listScans(_StubH5File([(name, b"title"), ("3.1", b"ok")])) + + assert [scanno for scanno, _ in entries] == [3] + + +@pytest.mark.parametrize("backend", BLISS_BACKENDS) +def test_a_scan_without_a_title_falls_back_to_its_name(backend): + entries = backend.listScans(_StubH5File([("7.1", None)])) + + assert entries == [(7, "7.1")] + + +# --- the legacy ch5523 format, named "_" ---------------- + +LEGACY_FILE = _StubH5File( + [ + ("dscan_3", b"dscan mu 0 1 20 1"), + ("ascan_12", b"ascan th 0 90 90 1"), + ("ascan_7", b"ascan th 0 45 45 1"), + ("instrument", b"not a scan"), + ] +) + + +def test_the_legacy_backend_lists_its_own_naming_convention(): + """It has no subscans, so Fastscan's implementation does not apply.""" + entries = id31_tools.BlissScan.listScans(LEGACY_FILE) + + assert [name for name, _ in entries] == ["dscan_3", "ascan_7", "ascan_12"] + + +def test_the_legacy_backend_reports_the_name_it_is_opened_with(): + """openScan passes ddict['name'] to BlissScan, not the scan number. + + Reporting the number would make orGUI look up a group that does not exist. + """ + entries = dict(id31_tools.BlissScan.listScans(LEGACY_FILE)) + + assert "ascan_12" in entries + assert entries["ascan_12"] == "ascan th 0 90 90 1" + + +def test_the_legacy_backend_ignores_entries_that_are_not_scans(): + entries = id31_tools.BlissScan.listScans(LEGACY_FILE) + + assert all(name != "instrument" for name, _ in entries) + + +def test_the_legacy_backend_does_not_list_subscan_style_names(): + """A "1.1" name is not a scan of this format and must not be offered.""" + assert id31_tools.BlissScan.listScans(BLISS_FILE) == [] + + +def _example_backend(filename, classname): + """Load a shipped example backend, or skip when it is not in the tree.""" + path = ( + pathlib.Path(__file__).resolve().parents[3] / "examples" / "backend" / filename + ) + if not path.is_file(): + pytest.skip("examples are not part of an installed orGUI") + spec = importlib.util.spec_from_file_location(path.stem + "_example", path) + module = importlib.util.module_from_spec(spec) + spec.loader.exec_module(module) + return getattr(module, classname) + + +def test_the_example_backend_implements_the_fast_path(): + """The reviewed case: a backend with no "BlissScan" in its name.""" + backend = _example_backend("ID31_EBS_p4_backend.py", "ID31_EBS_p4_2026") + + assert "BlissScan" not in backend.__name__ + assert backend.listScans(BLISS_FILE) == id31_tools.BlissScan_EBS_p4.listScans( + BLISS_FILE + ) + + +# --- the CHESS QM2 example, SPEC style names with one scan per group ------ + + +def _qm2(): + return _example_backend("CHESS_QM2.py", "QM2_backend_2026") + + +def _qm2_resolves(h5file, scanno): + """The group QM2_backend_2026.__init__ would open for this identifier.""" + for name in h5file: + head, _, tail = str(name).partition(".") + if int(head) == scanno: + return ".".join((head, tail)) + return None + + +QM2_FILE = _StubH5File( + [ + ("2.1", b"ascan th 0 90 90 1"), + ("1.1", b"ascan th 0 45 45 1"), + ("10.1", b"dscan mu 0 1 20 1"), + ] +) + + +def test_qm2_lists_the_groups_of_the_file_root(): + entries = _qm2().listScans(QM2_FILE) + + assert [scanno for scanno, _ in entries] == [1, 2, 10] + assert dict(entries)[2] == "ascan th 0 90 90 1" + + +def test_qm2_reports_the_scan_number_its_constructor_takes(): + """__init__ matches int(name before the dot), so the number is the id.""" + backend = _qm2() + + for scanno, _title in backend.listScans(QM2_FILE): + assert _qm2_resolves(QM2_FILE, scanno) is not None + + +def test_qm2_offers_a_scan_number_once_per_file(): + """__init__ opens the first group with a matching number. + + Offering the later subscans would let one be selected and another opened. + """ + file = _StubH5File([("1.1", b"first"), ("1.2", b"second"), ("2.1", b"other")]) + + entries = _qm2().listScans(file) + + assert [scanno for scanno, _ in entries] == [1, 2] + assert dict(entries)[1] == "first" + + +def test_qm2_skips_groups_that_are_not_scans(): + """__init__ would raise on these, so they must not be offered.""" + file = _StubH5File( + [("1.1", b"real scan"), ("instrument", b"meta"), ("title", None)] + ) + + assert [scanno for scanno, _ in _qm2().listScans(file)] == [1] diff --git a/orgui/main.py b/orgui/main.py index 1384527..674ffc3 100644 --- a/orgui/main.py +++ b/orgui/main.py @@ -150,17 +150,17 @@ def main(): "-l", dest="locking", action="store_true", - help="HDF5_USE_FILE_LOCKING=True (default)", + help="HDF5_USE_FILE_LOCKING=True", ) locking_parser.add_argument( "--no-hdflocking", "-nl", dest="locking", action="store_false", - help="HDF5_USE_FILE_LOCKING=False", + help="HDF5_USE_FILE_LOCKING=False (default)", ) - parser.set_defaults(locking=True) + parser.set_defaults(locking=False) options = parser.parse_args() diff --git a/orgui/resources/icons/q-plot.png b/orgui/resources/icons/q-plot.png new file mode 100644 index 0000000..aef7abe Binary files /dev/null and b/orgui/resources/icons/q-plot.png differ diff --git a/orgui/resources/icons/q-plot.svg b/orgui/resources/icons/q-plot.svg new file mode 100644 index 0000000..fd79a32 --- /dev/null +++ b/orgui/resources/icons/q-plot.svg @@ -0,0 +1,13 @@ + + + Q-plot + Momentum transfer: a capital Q. + + + + + + +