From 47664234562c5723024fa052110c0ce9e8423846 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 16:31:13 +0200 Subject: [PATCH 01/10] feat!: change the default hdf5 locking setting for live data analysis and working with hdf5 files where new data is still added, the hdf5 locking must be turned off --- orgui/main.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/orgui/main.py b/orgui/main.py index c080113..ea31b35 100644 --- a/orgui/main.py +++ b/orgui/main.py @@ -117,7 +117,7 @@ def main(): locking_parser.add_argument('--hdflocking', '-l', dest='locking', action='store_true', help="HDF5_USE_FILE_LOCKING=True (default)") locking_parser.add_argument('--no-hdflocking', '-nl', dest='locking', action='store_false', help="HDF5_USE_FILE_LOCKING=False") - parser.set_defaults(locking=True) + parser.set_defaults(locking=False) options = parser.parse_args() From 409c4d7f4cd0db213c9ddcc80ab4fd4cf605bd80 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 16:39:50 +0200 Subject: [PATCH 02/10] feat!: change how id31 scan objects are recognized As ESRF fastscans contain two separate hdf5 nodes, the data must be treated separately when loading a segmented scan. Requirement: class name includes 'BlissScan' --- orgui/app/orGUI.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index 7503ce5..ac6cacf 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -2570,11 +2570,11 @@ def _onLoadInterlacedScan(self): else: h5file = model.data(rootI, role=silx.gui.hdf5.Hdf5TreeModel.H5PY_OBJECT_ROLE) - isID31 = self.scanSelector.btid.currentText() in ['ch5523','ch5700','ch5918','ch6392','ch7131','ch7149','ch7856','ch8153','id31_default'] + backendcls = backends.fscans[self.scanSelector.btid.currentText()] kl_full = list(h5file.keys()) kl = np.empty(0,dtype=int) for i in kl_full: - if isID31: + if 'BlissScan' in backendcls.__name__: pattern = r'\.\d' result = re.findall(pattern, i)[0][1:] if result == '1': @@ -2584,7 +2584,7 @@ def _onLoadInterlacedScan(self): kl = np.append(kl,i) # separate scan nr and delete duplicates suffixes - if isID31: + if 'BlissScan' in backendcls.__name__: # try to get the scan nr and '/title' from the hdf5 file nr = np.empty(0,dtype=int) name = np.empty(0,dtype=str) @@ -2621,7 +2621,7 @@ def _onLoadInterlacedScan(self): ithScanBox = qt.QCheckBox() scanBoxes.append(ithScanBox) #labels.append(qt.QLabel('Scan '+str(item)+':'+str(name[i][2:-1]))) - if isID31: + if 'BlissScan' in backendcls.__name__: b.addRow(qt.QLabel('Scan '+str(item)+': '+name[i].decode()),ithScanBox) else: b.addRow(qt.QLabel('Scan '+str(item)+': '+name[i]),ithScanBox) From 4afbeb323bb0ffe8e9502f080fcce123e307ab4e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 16:43:56 +0200 Subject: [PATCH 03/10] feat: untoggle max/sum images when changing image nr the max image or sum image are automatically deactivated when the user changes the image using the slider below the central plot --- orgui/app/orGUI.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index ac6cacf..befc783 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -2986,11 +2986,19 @@ def _onChangeImage(self,imageno): ) return self.scanSelector.slider.setValue(imageno) + 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()) def _onSliderValueChanged(self,value): """GUI/CLI hint: replot the image selected by the scan slider.""" if self.fscan is not None: + if self.scanSelector.showMaxAct.isChecked(): + self.scanSelector.showMaxAct.setChecked(False) + if self.scanSelector.showSumAct.isChecked(): + self.scanSelector.showSumAct.setChecked(False) self.plotImage(value) #print(self.centralPlot._callback) From c40b4439176b215c541fc8ed02e7cedc1cba01b3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 16:48:24 +0200 Subject: [PATCH 04/10] fix: max/sum icon sometimes not correct fix a bug that causes the Qt icon to be deactivated even though the max/sum image is shown in the central plot --- orgui/app/orGUI.py | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index befc783..bb35296 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -3104,6 +3104,8 @@ def _onMaxToggled(self,value): replace=False,resetzoom=False,copy=True,z=1) self.centralPlot.setActiveImage(self.currentAddImageLabel) self.scanSelector.alphaslider.setLegend(self.currentAddImageLabel) + if not self.scanSelector.showMaxAct.isChecked(): + self.scanSelector.showMaxAct.setChecked(True) else: self.scanSelector.showMaxAct.setChecked(False) else: @@ -3133,6 +3135,8 @@ def _onSumToggled(self,value): replace=False,resetzoom=False,copy=True,z=1) self.centralPlot.setActiveImage(self.currentAddImageLabel) self.scanSelector.alphaslider.setLegend(self.currentAddImageLabel) + if not self.scanSelector.showSumAct.isChecked(): + self.scanSelector.showSumAct.setChecked(True) else: self.scanSelector.showSumAct.setChecked(False) else: From cd2722c02c86aa2e22cada6171af4269a9c3abd5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 17:08:05 +0200 Subject: [PATCH 05/10] feat: refreshing of id31 hdf5 files now possible this commit allows the refreshing of hdf5 files while new scans in the same dataset are still performed - without causing broken link errors at ESRF beamline id31 --- orgui/app/QScanSelector.py | 15 ++++++++++---- orgui/app/orGUI.py | 31 ++++++++++++++++++++++++++++ orgui/backend/beamline/id31_tools.py | 4 ++-- 3 files changed, 44 insertions(+), 6 deletions(-) diff --git a/orgui/app/QScanSelector.py b/orgui/app/QScanSelector.py index 1968603..10dbd8a 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) self.parentmainwindow = parentmainwindow @@ -103,7 +104,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) @@ -1013,6 +1014,10 @@ def __getRelativePath(self, model, rootIndex, index): def _onRefreshFile(self): + self.sigRefreshH5.emit() + # delete fscan object (to avoid broken links), and potentially scan_image and currentAddImageLabel (to avoid crash) + + def _onDoRefresh(self): qt.QApplication.setOverrideCursor(qt.Qt.WaitCursor) selection = self.hdfTreeView.selectionModel() @@ -1052,6 +1057,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: @@ -1066,7 +1072,7 @@ def _onRefreshFile(self): pass # Was already closed ''' - self.createTreeView() + self.createTreeView() # why create new tree view? silx view handles resync different... maintab = self.mainwidget.findChildren(qt.QTabWidget)[0] maintab.removeTab(0) @@ -1077,6 +1083,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!!! #self.__expandNodesFromPaths(self.hdfTreeView, index, paths) #self.hdf5model.appendFile(filename) qt.QApplication.restoreOverrideCursor() @@ -1085,12 +1092,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/orGUI.py b/orgui/app/orGUI.py index bb35296..ea03e73 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -217,6 +217,7 @@ def __init__(self,configfile,parent=None): + self.scanSelector.sigRefreshH5.connect(self._onRefreshH5) self.scanSelector.sigImageNoChanged.connect(self._onSliderValueChanged) self.scanSelector.sigImagePathChanged.connect(self._onImagePathChanged) @@ -463,6 +464,36 @@ def __init__(self,configfile,parent=None): if configfile is not None: self.ubcalc.readConfig(configfile) + def _onRefreshH5(self): + ''' + # manual way for console use to avoid broken nodes after update: + 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 + + autoload = self.autoLoadAct.isChecked() + bt_autodetect = self.scanSelector.bt_autodetect_enable.isChecked() + if autoload: + self.autoLoadAct.setChecked(False) # no other way to avoid re-generation of max/sum when loading scan + if bt_autodetect: + self.scanSelector.bt_autodetect_enable.setChecked(False) + + self.scanSelector._onLoadScan() # very important! this seems to unlock the file somehow + + # 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 + self.scanSelector.hdfTreeView.expandToDepth(0) # here the call is possible + # 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 diff --git a/orgui/backend/beamline/id31_tools.py b/orgui/backend/beamline/id31_tools.py index 3e0bd1e..cfc493d 100644 --- a/orgui/backend/beamline/id31_tools.py +++ b/orgui/backend/beamline/id31_tools.py @@ -229,7 +229,7 @@ def __init__(self, hdffilepath_orNode=None, scanno=None, loadimg=True, saveh5dat _,filename_noext = os.path.split(filepath) #filename_noext = filename.split('.')[0] 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] @@ -558,7 +558,7 @@ def __init__(self, hdffilepath_orNode=None, scanno=None, loadimg=True, saveh5dat _,filename_noext = os.path.split(filepath) #filename_noext = filename.split('.')[0] 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] From b658f14f915e7c3a82900612801c17c6bc3e5491 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 17:13:35 +0200 Subject: [PATCH 06/10] feat: change default detector, backend and geometry settings The Pilatus 4 4M detector is now set as default, the fallback geometry is more universal and the ch8153 beamtime has been added to the backend selection --- orgui/app/QUBCalculator.py | 12 ++++++------ orgui/backend/backends.py | 9 ++++++++- 2 files changed, 14 insertions(+), 7 deletions(-) diff --git a/orgui/app/QUBCalculator.py b/orgui/app/QUBCalculator.py index 6f212b8..57673cc 100644 --- a/orgui/app/QUBCalculator.py +++ b/orgui/app/QUBCalculator.py @@ -647,11 +647,11 @@ def readConfig(self,configfile): return False def toFallbackConfig(self): - sdd = 0.729 #m - E = 78. - pixelsize = 172e-6 - cp = [731.0,1587.856] - self.mu = np.deg2rad(0.05) + sdd = 1.0 #m + E = 75. + pixelsize = 150e-6 + cp = [1100.0,2020.0] + self.mu = np.deg2rad(0.08) self.chi = 0. self.phi = 0. self.n = 1 - 1.1415e-06 @@ -666,7 +666,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) self.detectorCal.wavelength = self.ubCal.getLambda()*1e-10 self.detectorCal.setAzimuthalReference(np.deg2rad(90.)) diff --git a/orgui/backend/backends.py b/orgui/backend/backends.py index 82584ef..4d9f700 100644 --- a/orgui/backend/backends.py +++ b/orgui/backend/backends.py @@ -44,7 +44,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)), '20190017': (datetime(2019, 12, 8), datetime(2019, 12, 24)), @@ -89,6 +89,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 } @@ -144,6 +145,12 @@ 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']) From d5a9c87996f9590b8872ac865120c5765855d9b3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Finn=20Schr=C3=B6ter?= Date: Mon, 29 Jun 2026 17:19:03 +0200 Subject: [PATCH 07/10] feat: add functionality to convert detector images into Q coordinates using the pyFAI FiberIntegrator package, grazing incidence detector images can now easily be converted from pixel coordinates to their Q values in units of inverse Angstrom --- orgui/app/orGUI.py | 169 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 169 insertions(+) diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index ea03e73..c61ba06 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -64,6 +64,11 @@ except: console = False +from packaging.version import Version +import pyFAI.version +if Version(pyFAI.version) >= Version('2025.01'): + from pyFAI.integrator.fiber import FiberIntegrator + import traceback from . import qutils, ROIutils @@ -162,6 +167,10 @@ def __init__(self,configfile,parent=None): self.centralPlot.setCallback(self._graphCallback) toolbar = qt.QToolBar() toolbar.addAction(control_actions.OpenGLAction(parent=toolbar, plot=self.centralPlot)) + self.plotAgainstQAct = qt.QAction("Q-plot",self) + self.plotAgainstQAct.setCheckable(True) + self.plotAgainstQAct.toggled.connect(self._convertImagetoQ) + toolbar.addAction(self.plotAgainstQAct) self.centralPlot.addToolBar(toolbar) self.currentImageLabel = None @@ -3017,6 +3026,8 @@ def _onChangeImage(self,imageno): ) return self.scanSelector.slider.setValue(imageno) + if self.plotAgainstQAct.isChecked(): + self.plotAgainstQAct.setChecked(False) if self.scanSelector.showMaxAct.isChecked(): self.scanSelector.showMaxAct.setChecked(False) if self.scanSelector.showSumAct.isChecked(): @@ -3026,6 +3037,8 @@ def _onChangeImage(self,imageno): def _onSliderValueChanged(self,value): """GUI/CLI hint: replot the image selected by the scan slider.""" if self.fscan is not None: + if self.plotAgainstQAct.isChecked(): + self.plotAgainstQAct.setChecked(False) if self.scanSelector.showMaxAct.isChecked(): self.scanSelector.showMaxAct.setChecked(False) if self.scanSelector.showSumAct.isChecked(): @@ -3117,6 +3130,8 @@ def readfile_max(imgno): def _onMaxToggled(self,value): """GUI/CLI hint: toggle display of the precomputed maximum image.""" + if self.plotAgainstQAct.isChecked(): + self.plotAgainstQAct.setChecked(False) if self.scanSelector.showSumAct.isChecked(): self.scanSelector.showSumAct.setChecked(False) if value: @@ -3148,6 +3163,8 @@ def _onMaxToggled(self,value): def _onSumToggled(self,value): """GUI/CLI hint: toggle display of the precomputed summed image.""" + if self.plotAgainstQAct.isChecked(): + self.plotAgainstQAct.setChecked(False) if self.scanSelector.showMaxAct.isChecked(): self.scanSelector.showMaxAct.setChecked(False) if value: @@ -3177,6 +3194,158 @@ def _onSumToggled(self,value): self.currentAddImageLabel = None + def _convertImagetoQ(self,value): + if value: + # check if conversion possible + if Version(pyFAI.version) < Version('2025.1'): + print('conversion of image to Q not possible! Install pyFAI > 2025.1 !') + self.plotAgainstQAct.setChecked(False) + return + if self.fscan is None: + if logger_utils.get_logging_context() == 'gui': + qt.QMessageBox.critical(self, + "Q conversion failed", + "Cannot convert image to Q\n"\ + "No scan loaded!") + else: + logger.exception('Q conversion failed', + extra={'title' : 'Cannot convert image to Q', + 'description' : 'No scan loaded!', + 'show_dialog' : False, + "dialog_level" : logging.ERROR, + 'parent' : self}) + self.plotAgainstQAct.setChecked(False) + return + if self.fscan.axisname == 'mu': + if logger_utils.get_logging_context() == 'gui': + qt.QMessageBox.critical(self, + "Q conversion failed", + "Cannot convert image to Q\n"\ + "Conversion not implemented for mu scans!") + else: + logger.exception('Q conversion failed', + extra={'title' : 'Cannot convert image to Q', + 'description' : 'Conversion not implemented for mu scans!', + 'show_dialog' : False, + "dialog_level" : logging.ERROR, + 'parent' : self}) + self.plotAgainstQAct.setChecked(False) + return + + # hide scan image, remove max/sum image, get data for conversion + for i in self.centralPlot.getAllImages(): + if i.getLegend() == "scan_image": + self.centralPlot.hideCurve(i) + if not self.scanSelector.showMaxAct.isChecked() and not self.scanSelector.showSumAct.isChecked(): + data = i.getData() + elif i.getLegend() != "scan_image": + self.centralPlot.removeImage(i) + if self.scanSelector.showMaxAct.isChecked(): + data = self.allimgmax + elif self.scanSelector.showSumAct.isChecked(): + data = self.allimgsum + + # determine sample orientation + azim = self.ubcalc.detectorCal.getAzimuthalReference() + if -np.pi/4 < azim < np.pi/4: + orientation = 7 + elif np.pi/4 < azim < np.pi*3/4: + orientation = 4 + if self.centralPlot.isYAxisInverted(): + self.centralPlot.setYAxisInverted(False) + else: + self.centralPlot.setYAxisInverted(True) + elif np.pi*3/4 < azim < np.pi*5/4: + orientation = 8 + if self.centralPlot.isXAxisInverted(): + self.centralPlot.setXAxisInverted(False) + else: + self.centralPlot.setXAxisInverted(True) + elif np.pi*5/4 < azim < np.pi*7/4: + orientation = 1 + else: + orientation = 1 + + # perform conversion into Q coordinates + fi = FiberIntegrator(dist=self.ubcalc.detectorCal.dist, poni1=self.ubcalc.detectorCal.poni1, poni2=self.ubcalc.detectorCal.poni2, + wavelength=self.ubcalc.detectorCal.wavelength, + rot1=self.ubcalc.detectorCal.rot1, rot2=self.ubcalc.detectorCal.rot2, rot3=self.ubcalc.detectorCal.rot3, + detector=self.ubcalc.detectorCal.detector, + ) + res2d = fi.integrate2d_grazing_incidence(data, sample_orientation=orientation, + incident_angle=self.ubcalc.mu, + tilt_angle=0, + unit_oop="qoop_A^-1",unit_ip="qip_A^-1") + + # plot generated image + oopmin, oopmax = np.min(res2d.outofplane), np.max(res2d.outofplane) + ipmin, ipmax = np.min(res2d.inplane), np.max(res2d.inplane) + orig_ip, orig_oop = res2d.inplane[0], res2d.outofplane[0] + if orientation in [1,4]: # specular axis on vertical detector axis + 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=((ipmax-ipmin)/1000,(oopmax-oopmin)/1000), + origin=(orig_ip,orig_oop), + ) + # apply correct zoom + self.centralPlot.getXAxis().setLimits(ipmin,ipmax) + self.centralPlot.getYAxis().setLimits(oopmin,oopmax) + else: # specular axis on horizontal detector axis + self.currentAddImageLabel = self.centralPlot.addImage(res2d.intensity.T,legend='qImage',replace=False, + resetzoom=False,copy=True,z=2, + xlabel=r"q$_\perp / \, \AA^{-1}$", + ylabel=r"q$_\parallel / \, \AA^{-1}$", + scale=((oopmax-oopmin)/1000,(ipmax-ipmin)/1000), + origin=(orig_oop,orig_ip), + ) + # apply correct zoom + self.centralPlot.getYAxis().setLimits(ipmin,ipmax) + self.centralPlot.getXAxis().setLimits(oopmin,oopmax) + + # apply active plot settings + self.centralPlot.setActiveImage(self.currentAddImageLabel) + self.scanSelector.alphaslider.setLegend(self.currentAddImageLabel) + + else: + # restore status from before + if self.currentAddImageLabel is not None: + # show scan image, delete qImage + for i in self.centralPlot.getAllImages(): + if i.getLegend() == "scan_image": + self.centralPlot.hideCurve(i,False) + elif i.getLegend() == 'qImage': + self.centralPlot.removeImage(i) + + # 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) + + # restore status of xaxis and yaxis + if self.fscan is not None: + if self.fscan.axisname != 'mu' and Version(pyFAI.version) > Version('2025.1'): + azim = self.ubcalc.detectorCal.getAzimuthalReference() + if np.pi/4 < azim < np.pi*3/4: + if not self.centralPlot.isYAxisInverted(): + self.centralPlot.setYAxisInverted(True) + else: + self.centralPlot.setYAxisInverted(False) + elif np.pi*3/4 < azim < np.pi*5/4: + if not self.centralPlot.isXAxisInverted(): + self.centralPlot.setXAxisInverted(True) + else: + self.centralPlot.setXAxisInverted(False) + + # apply correct zoom for pixel coordinates + self.centralPlot.resetZoom() def plotImage(self,key=0): From 9c99984a580b1d689381d315b8e6041a95091b6f Mon Sep 17 00:00:00 2001 From: "T. Fuchs" Date: Tue, 30 Jun 2026 21:41:44 -0400 Subject: [PATCH 08/10] feat(backend): Update ID31 backend with p4 to new API --- examples/backend/ID31_EBS_p4_backend.py | 504 ++++++++++++++++++++++++ 1 file changed, 504 insertions(+) create mode 100644 examples/backend/ID31_EBS_p4_backend.py diff --git a/examples/backend/ID31_EBS_p4_backend.py b/examples/backend/ID31_EBS_p4_backend.py new file mode 100644 index 0000000..a9cfde2 --- /dev/null +++ b/examples/backend/ID31_EBS_p4_backend.py @@ -0,0 +1,504 @@ +# /*########################################################################## +# +# Copyright (c) 2020-2026 Timo Fuchs +# +# 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} + + @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 From 57a3971a9e817219651739375161080d90f020ac Mon Sep 17 00:00:00 2001 From: "T. Fuchs" Date: Wed, 29 Jul 2026 18:17:14 -0400 Subject: [PATCH 09/10] feat(app): rework the experimental Q-plot with selectable frames Make the reciprocal-space conversion exact for arbitrary azimuthal references, use orGUI diffractometer conventions for single crystals instead of the 2d FiberIntegrator. - add the Q_alpha, Q_lab, Q_omega, Q_chi, Q_phi and Q_cryst frames with a selector next to the Q-plot action, defaulting to Q_alpha - add the _qconversion_cpp extension, which converts a 6.2 megapixel Pilatus 6M image in ~50 ms instead of ~740 ms, with a numpy fallback - reuse the pyFAI integrator while the conversion is unchanged, keyed on the collapsed conversion coefficients, so stepping through images no longer resets the integrator for every image - keep the Q-plot switched on and rebuild it in place when the image, the max/sum selection, the frame, the machine angles, the energy or the orientation matrix change - pass a legend rather than an image item to the alpha slider - take the manual backend selection from backends.default_beamtime --- CHANGELOG.md | 37 +- doc/source/geometry.rst | 186 ++++-- doc/source/release_notes.rst | 9 +- meson.build | 8 + orgui/app/QScanSelector.py | 4 +- orgui/app/QUBCalculator.py | 19 +- orgui/app/cpp/qconversion_cpp.cpp | 106 ++++ orgui/app/orGUI.py | 297 +++++---- orgui/app/qconversion.py | 585 ++++++++++++++++++ .../app/test/test_parameter_dialog_cancel.py | 110 ++++ orgui/app/test/test_q_conversion.py | 419 +++++++++++++ orgui/app/test/test_q_plot_orientation.py | 145 +++-- .../xrayutils/test/test_q_conversion.py | 238 ------- 13 files changed, 1696 insertions(+), 467 deletions(-) create mode 100644 orgui/app/cpp/qconversion_cpp.cpp create mode 100644 orgui/app/qconversion.py create mode 100644 orgui/app/test/test_parameter_dialog_cancel.py create mode 100644 orgui/app/test/test_q_conversion.py delete mode 100644 orgui/datautils/xrayutils/test/test_q_conversion.py diff --git a/CHANGELOG.md b/CHANGELOG.md index 74c3dfc..f16b011 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -66,15 +66,42 @@ 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 fixes: + +- 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. The conversion - is equivalent to orGUI's own per-pixel ``QAlpha`` calculation; see the - geometry documentation for when to use which. + 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. @@ -86,6 +113,10 @@ ESRF ID31 beamline support and reciprocal-space display: - 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. - **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 diff --git a/doc/source/geometry.rst b/doc/source/geometry.rst index bd5c59a..7fbf572 100644 --- a/doc/source/geometry.rst +++ b/doc/source/geometry.rst @@ -130,32 +130,29 @@ 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: +an image. orGUI drives it through :mod:`orgui.app.qconversion`, which supplies +pyFAI with orGUI's own angle conventions: .. code-block:: python - from pyFAI.integrator.fiber import FiberIntegrator - - fi = 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, - ) - res = fi.integrate2d_grazing_incidence( + from orgui.app import qconversion + + res = qconversion.integrateImage( + detectorCal, image, - sample_orientation=4, # depends on the azimuthal reference, see below - incident_angle=alpha_i, - tilt_angle=0, - unit_ip="qip_A^-1", - unit_oop="qoop_A^-1", + alpha_i, # incidence angle in rad + frame="Q_alpha", # see "Reciprocal-space frames" below ) res.intensity # rebinned image, shape (npt_oop, npt_ip) @@ -196,47 +193,124 @@ 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. -Sample orientation and sign conventions -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +Reciprocal-space frames +~~~~~~~~~~~~~~~~~~~~~~~ -The two routes place the surface normal on the detector differently. orGUI -derives it from the azimuthal reference of the detector calibration, while -pyFAI expects the EXIF-style ``sample_orientation`` argument. The ``Q-plot`` -action maps one onto the other: +``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: 34 33 33 - - * - Azimuthal reference - - ``sample_orientation`` - - Flipped plot axis - * - 0 degrees - - 7 - - none - * - 90 degrees (default) - - 4 - - vertical - * - 180 degrees - - 8 - - horizontal - * - 270 degrees - - 1 - - none - -pyFAI pairs each orientation with an in-plane mirrored partner -- (6, 7), -(3, 4), (5, 8) and (1, 2). Partners return an identical :math:`q_\perp` and -differ only in the **sign** of :math:`q_\parallel`. Because ``QAlpha`` returns -the unsigned radial in-plane component, both members of a pair are equally -correct for it; the sign is a display convention, which is why the ``Q-plot`` -action flips a plot axis for some orientations. - -With that convention fixed, the two routes agree to numerical precision. This -is verified in ``orgui/datautils/xrayutils/test/test_q_conversion.py``, which -compares them per pixel for flat and tilted detectors at several incidence -angles, and additionally checks that a single bright pixel is placed by -``integrate2d_grazing_incidence`` within one grid cell of the position -predicted by ``QAlpha``. + :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`. + +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 38ed200..0b3266c 100644 --- a/doc/source/release_notes.rst +++ b/doc/source/release_notes.rst @@ -28,13 +28,20 @@ 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 fixes: + +- 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. The conversion is equivalent to orGUI's own per-pixel ``QAlpha`` calculation; see the geometry documentation for when to use which. +- 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. - **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:** ID31-style scan objects are now recognized from the configured backend class instead of a hardcoded list of beamtime ids. This fixes the ``id31_default_p4`` backend, which the old list did not cover, and makes new ID31 beamtimes work without code changes. Custom backends whose class name does not contain ``BlissScan`` are no longer treated as ID31-style, even if their beamtime id was previously listed. 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 ca945ca..b0079dd 100644 --- a/orgui/app/QScanSelector.py +++ b/orgui/app/QScanSelector.py @@ -171,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 diff --git a/orgui/app/QUBCalculator.py b/orgui/app/QUBCalculator.py index 6bd1aca..443bde3 100644 --- a/orgui/app/QUBCalculator.py +++ b/orgui/app/QUBCalculator.py @@ -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 1cbb510..9c26a16 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -37,6 +37,7 @@ import sys import os import re +from contextlib import contextmanager from silx.gui import qt from io import StringIO @@ -57,10 +58,6 @@ from silx.gui.plot.tools.roi import RegionOfInterestManager from silx.gui.plot.actions import control as control_actions -import pyFAI -import pyFAI.version -from packaging.version import Version - import traceback from . import qutils, ROIutils, autoBraggWorkflow @@ -82,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: @@ -97,9 +95,7 @@ # Reciprocal-space (Q) conversion of detector images relies on pyFAI's # FiberIntegrator, which was introduced in pyFAI 2025.1. -HAS_FIBER_INTEGRATOR = Version(pyFAI.version) >= Version("2025.1") -if HAS_FIBER_INTEGRATOR: - from pyFAI.integrator.fiber import FiberIntegrator +HAS_FIBER_INTEGRATOR = qconversion.HAS_FIBER_INTEGRATOR try: from . import _roi_sum_accel @@ -187,8 +183,27 @@ def __init__(self, configfile, parent=None): ) self.plotAgainstQAct = qt.QAction("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._qPlotRefreshing = False self.centralPlot.addToolBar(toolbar) self.maskManager = MaskManager() self.maskConfigDialog = MaskConfigDialog(self.maskManager, self) @@ -351,6 +366,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( @@ -3872,24 +3890,24 @@ def _onChangeImage(self, imageno): ) return self.scanSelector.slider.setValue(imageno) - if self.plotAgainstQAct.isChecked(): - self.plotAgainstQAct.setChecked(False) - 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()) + 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: - if self.plotAgainstQAct.isChecked(): - self.plotAgainstQAct.setChecked(False) - if self.scanSelector.showMaxAct.isChecked(): - self.scanSelector.showMaxAct.setChecked(False) - if self.scanSelector.showSumAct.isChecked(): - self.scanSelector.showSumAct.setChecked(False) - 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): @@ -3989,8 +4007,6 @@ def readfile_max(imgno): def _onMaxToggled(self, value): """GUI/CLI hint: toggle display of the precomputed maximum image.""" - if self.plotAgainstQAct.isChecked(): - self.plotAgainstQAct.setChecked(False) if self.scanSelector.showSumAct.isChecked(): self.scanSelector.showSumAct.setChecked(False) if value: @@ -4021,7 +4037,7 @@ 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: @@ -4031,11 +4047,11 @@ 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.""" - if self.plotAgainstQAct.isChecked(): - self.plotAgainstQAct.setChecked(False) if self.scanSelector.showMaxAct.isChecked(): self.scanSelector.showMaxAct.setChecked(False) if value: @@ -4066,7 +4082,7 @@ 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: @@ -4076,6 +4092,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.""" @@ -4091,34 +4109,65 @@ def _set_center_roi_preview_color(self, roi): roi.setColor("blue" if self._roi_preview_enabled() else "red") roi.setBgStyle("pink", "-", roi.getLineWidth()) - def _qConversionSampleOrientation(self): - """Return the pyFAI sample orientation for the current azimuthal reference. - - The azimuthal reference fixes where the surface normal points on the - detector, which in pyFAI terms is the EXIF-style ``sample_orientation``. - The second return value tells whether the vertical (``"y"``), the - horizontal (``"x"``) or no (``None``) plot axis has to be flipped so - that the converted image keeps the same handedness as the pixel image. + 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: re-render the Q-plot after the frame selection changed.""" + del index + 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``. """ - azim = self.ubcalc.detectorCal.getAzimuthalReference() - if -np.pi / 4 < azim < np.pi / 4: - return 7, None - elif np.pi / 4 < azim < np.pi * 3 / 4: - return 4, "y" - elif np.pi * 3 / 4 < azim < np.pi * 5 / 4: - return 8, "x" - else: - return 1, None + 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. - The currently displayed image is rebinned onto a regular grid of - in-plane and out-of-plane momentum transfer using pyFAI's - ``FiberIntegrator``. The per-pixel momentum transfer used everywhere - else in orGUI is calculated by - :meth:`~orgui.datautils.xrayutils.HKLVlieg.VliegAngles.QAlpha`; both - routes are equivalent, see the geometry section of the documentation. + **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 @@ -4190,6 +4239,27 @@ def _convertImagetoQ(self, value): 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 + # hide scan image, remove max/sum image for i in self.centralPlot.getAllImages(): if i.getLegend() == "scan_image": @@ -4197,80 +4267,53 @@ def _convertImagetoQ(self, value): else: self.centralPlot.removeImage(i) - # determine sample orientation - orientation, flipaxis = self._qConversionSampleOrientation() - if flipaxis == "y": - self.centralPlot.setYAxisInverted( - not self.centralPlot.isYAxisInverted() - ) - elif flipaxis == "x": - self.centralPlot.setXAxisInverted( - not self.centralPlot.isXAxisInverted() - ) + # Reciprocal space is plotted with the ordinate pointing upwards, + # 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() + self.centralPlot.setYAxisInverted(False) + self.centralPlot.setXAxisInverted(False) # perform conversion into Q coordinates - detectorCal = self.ubcalc.detectorCal - fi = 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, - ) - res2d = fi.integrate2d_grazing_incidence( + alpha, omega = self._qPlotAngles() + res2d = qconversion.integrateImage( + self.ubcalc.detectorCal, data, - sample_orientation=orientation, - incident_angle=self.ubcalc.mu, - tilt_angle=0, - unit_oop="qoop_A^-1", - unit_ip="qip_A^-1", + alpha, + frame=frame, + omega=omega, + chi=self.ubcalc.chi, + phi=self.ubcalc.phi, + U=self.ubcalc.ubCal.getU(), ) - # plot generated image + # 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) - orig_ip, orig_oop = res2d.inplane[0], res2d.outofplane[0] scale_ip = (ipmax - ipmin) / res2d.inplane.size scale_oop = (oopmax - oopmin) / res2d.outofplane.size - if orientation in [1, 4]: # specular axis on vertical detector axis - 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=(orig_ip, orig_oop), - ) - # apply correct zoom - self.centralPlot.getXAxis().setLimits(ipmin, ipmax) - self.centralPlot.getYAxis().setLimits(oopmin, oopmax) - else: # specular axis on horizontal detector axis - self.currentAddImageLabel = self.centralPlot.addImage( - res2d.intensity.T, - legend="qImage", - replace=False, - resetzoom=False, - copy=True, - z=2, - xlabel=r"q$_\perp / \, \AA^{-1}$", - ylabel=r"q$_\parallel / \, \AA^{-1}$", - scale=(scale_oop, scale_ip), - origin=(orig_oop, orig_ip), - ) - # apply correct zoom - self.centralPlot.getYAxis().setLimits(ipmin, ipmax) - self.centralPlot.getXAxis().setLimits(oopmin, oopmax) + 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 + # 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(self.currentAddImageLabel) + self.scanSelector.alphaslider.setLegend("qImage") else: # restore status from before @@ -4313,21 +4356,10 @@ def _convertImagetoQ(self, value): self.currentAddImageLabel = None self.centralPlot.setActiveImage(self.currentImageLabel) - # restore status of xaxis and yaxis - if ( - self.fscan is not None - and self.fscan.axisname != "mu" - and HAS_FIBER_INTEGRATOR - ): - _, flipaxis = self._qConversionSampleOrientation() - if flipaxis == "y": - self.centralPlot.setYAxisInverted( - not self.centralPlot.isYAxisInverted() - ) - elif flipaxis == "x": - self.centralPlot.setXAxisInverted( - not self.centralPlot.isXAxisInverted() - ) + # restore the image convention of the pixel coordinates + if self._qPlotPreviousYInverted is not None: + self.centralPlot.setYAxisInverted(self._qPlotPreviousYInverted) + self._qPlotPreviousYInverted = None # apply correct zoom for pixel coordinates self.centralPlot.resetZoom() @@ -5190,7 +5222,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]: @@ -5206,11 +5239,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) diff --git a/orgui/app/qconversion.py b/orgui/app/qconversion.py new file mode 100644 index 0000000..15a08b3 --- /dev/null +++ b/orgui/app/qconversion.py @@ -0,0 +1,585 @@ +"""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 + +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 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 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 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_q_conversion.py b/orgui/app/test/test_q_conversion.py new file mode 100644 index 0000000..12f959f --- /dev/null +++ b/orgui/app/test/test_q_conversion.py @@ -0,0 +1,419 @@ +"""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)) + + +@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 index c2913b2..593ba51 100644 --- a/orgui/app/test/test_q_plot_orientation.py +++ b/orgui/app/test/test_q_plot_orientation.py @@ -1,57 +1,132 @@ -"""The Q-plot action must use the sample orientations validated against QAlpha. +"""The Q-plot must be wired to the validated conversion in :mod:`qconversion`. -The numerical equivalence of the Q-plot conversion with orGUI's own QAlpha -calculation is established in -``orgui/datautils/xrayutils/test/test_q_conversion.py`` for a fixed mapping of -the azimuthal reference onto pyFAI's EXIF sample orientation. This test keeps -the GUI implementation tied to that validated mapping. +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 numpy as np import pytest -from orgui.datautils.xrayutils.test.test_q_conversion import AZIMUTH_ORIENTATION +from orgui.app import qconversion orGUI = pytest.importorskip("orgui.app.orGUI", reason="Qt bindings are required") -# Which plot axis has to be flipped so that the converted image keeps the same -# handedness as the pixel image, see orGUI._qConversionSampleOrientation. -EXPECTED_FLIP = {0.0: None, 90.0: "y", 180.0: "x", 270.0: None} +class _StubComboBox: + """Minimal stand-in for the frame selector combo box.""" -class _StubDetectorCal: - def __init__(self, azimuth): - self._azimuth = azimuth + def __init__(self, data): + self._data = data - def getAzimuthalReference(self): - return self._azimuth + def currentData(self): + return self._data -class _StubUBCalculator: - def __init__(self, azimuth): - self.detectorCal = _StubDetectorCal(azimuth) +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 + + +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] -class _StubOrGUI: - """Minimal stand-in: the method only reads the azimuthal reference.""" - def __init__(self, azimuth): - self.ubcalc = _StubUBCalculator(azimuth) +def test_refresh_is_not_reentrant(): + stub = _RefreshStub(checked=True, reentrant=True) + orGUI.orGUI._refreshQPlot(stub) -@pytest.mark.parametrize("azimuth_deg, orientation", AZIMUTH_ORIENTATION) -def test_orientation_matches_validated_mapping(azimuth_deg, orientation): - stub = _StubOrGUI(np.deg2rad(azimuth_deg)) + assert stub.calls == [True] + assert stub._qPlotRefreshing is False - result, flipaxis = orGUI.orGUI._qConversionSampleOrientation(stub) - assert result == orientation - assert flipaxis == EXPECTED_FLIP[azimuth_deg] +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) -def test_fiber_integrator_flag_matches_pyfai_version(): - """The Q-plot action is only offered when FiberIntegrator can be imported.""" - if orGUI.HAS_FIBER_INTEGRATOR: - assert hasattr(orGUI, "FiberIntegrator") - else: - assert not hasattr(orGUI, "FiberIntegrator") + assert stub.calls == [True] + assert stub._qPlotRefreshing is False diff --git a/orgui/datautils/xrayutils/test/test_q_conversion.py b/orgui/datautils/xrayutils/test/test_q_conversion.py deleted file mode 100644 index 2f1533c..0000000 --- a/orgui/datautils/xrayutils/test/test_q_conversion.py +++ /dev/null @@ -1,238 +0,0 @@ -"""Equivalence of orGUI's QAlpha conversion with pyFAI's FiberIntegrator. - -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 "Q-plot" action of the GUI instead -rebins the displayed image onto a regular grid of in-plane and out-of-plane -momentum transfer using pyFAI's ``FiberIntegrator``. - -Both routes have to describe the same reciprocal space. These tests pin that -down numerically: first per pixel, which is an exact identity, and then through -the full ``integrate2d_grazing_incidence`` call actually used by the GUI, where -the agreement is limited by the width of the rebinning grid. -""" - -import unittest - -import numpy as np -import pyFAI - -from .. import DetectorCalibration, HKLVlieg - -try: - from pyFAI import units as pyFAI_units - from pyFAI.integrator.fiber import FiberIntegrator - - fiber_available = True -except ImportError: # pyFAI < 2025.1 - fiber_available = False - -requires_fiber = unittest.skipUnless( - fiber_available, "pyFAI >= 2025.1 with FiberIntegrator is required" -) - -ENERGY = 78.0 # keV - -# Azimuthal reference (in degrees) -> pyFAI EXIF sample orientation, as -# implemented by orgui.app.orGUI.orGUI._qConversionSampleOrientation. -# -# pyFAI pairs every orientation with an in-plane mirrored partner -- -# (6, 7), (3, 4), (5, 8) and (1, 2). Partners share q_oop and differ only in -# the sign of q_ip, see TestSampleOrientationMirrorPairs. QAlpha returns the -# in-plane momentum transfer as the unsigned radial component -# sqrt(Qx**2 + Qy**2), so it cannot distinguish the partners and the tests -# below compare |q_ip|. -AZIMUTH_ORIENTATION = ((0.0, 7), (90.0, 4), (180.0, 8), (270.0, 1)) - -FLAT_DETECTOR = (0.0, 0.0, 0.0) -TILTED_DETECTOR = (np.deg2rad(1.7), np.deg2rad(-0.9), np.deg2rad(3.4)) - - -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 _vlieg_angles(): - lattice = HKLVlieg.Lattice(np.array([4.0, 4.0, 4.0]), np.array([90.0, 90.0, 90.0])) - return HKLVlieg.VliegAngles(HKLVlieg.UBCalculator(lattice, ENERGY)) - - -def _pixel_grid(det, stride=67): - """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 _q_from_qalpha(det, angles, alpha_i, rows, cols): - """In-plane and out-of-plane momentum transfer via surfaceAngles + QAlpha.""" - gamma, delta = det.surfaceAnglesPoint(rows, cols, alpha_i) - Qxyz = angles.QAlpha(alpha_i, delta, gamma) - q_ip = np.hypot(Qxyz[..., 0], Qxyz[..., 1]) - q_oop = Qxyz[..., 2] - return q_ip, q_oop - - -def _q_from_pyfai(det, alpha_i, orientation, rows, cols): - """Per-pixel q_ip and q_oop from the equations used by FiberIntegrator.""" - param = np.array( - [det.dist, det.poni1, det.poni2, det.rot1, det.rot2, det.rot3], - dtype=np.float64, - ) - # calc_pos_zyx returns (along beam, slow/vertical, fast/horizontal), while - # the pyFAI unit equations expect x=horizontal, y=vertical, z=beam. - t3, t1, t2 = det.calc_pos_zyx(d1=rows, d2=cols, param=param, do_parallax=True) - kwargs = { - "x": t2, - "y": t1, - "z": t3, - "wavelength": det.wavelength, - "incident_angle": alpha_i, - "tilt_angle": 0.0, - "sample_orientation": orientation, - } - # pyFAI returns inverse nanometre, orGUI works in inverse Angstrom. - return pyFAI_units.eq_qip(**kwargs) / 10.0, pyFAI_units.eq_qoop(**kwargs) / 10.0 - - -@requires_fiber -class TestQAlphaMatchesPyFAIPerPixel(unittest.TestCase): - """QAlpha and pyFAI agree pixel by pixel, before any rebinning.""" - - def test_matches_for_every_gui_orientation(self): - angles = _vlieg_angles() - for azimuth_deg, orientation in AZIMUTH_ORIENTATION: - for rotations in (FLAT_DETECTOR, TILTED_DETECTOR): - for alpha_deg in (0.0, 0.08, 2.5): - with self.subTest( - azimuth=azimuth_deg, - orientation=orientation, - 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) - - q_ip, q_oop = _q_from_qalpha(det, angles, alpha_i, rows, cols) - ref_ip, ref_oop = _q_from_pyfai( - det, alpha_i, orientation, rows, cols - ) - - np.testing.assert_allclose(q_oop, ref_oop, atol=1e-9) - np.testing.assert_allclose(q_ip, np.abs(ref_ip), atol=1e-9) - - def test_total_momentum_transfer_is_orientation_independent(self): - """|Q| is a physical invariant, so it must not depend on the labeling.""" - angles = _vlieg_angles() - det = _geometry(90.0) - rows, cols = _pixel_grid(det) - alpha_i = np.deg2rad(0.6) - - q_ip, q_oop = _q_from_qalpha(det, angles, alpha_i, rows, cols) - magnitude = np.hypot(q_ip, q_oop) - - for orientation in range(1, 9): - with self.subTest(orientation=orientation): - ref_ip, ref_oop = _q_from_pyfai(det, alpha_i, orientation, rows, cols) - np.testing.assert_allclose( - magnitude, np.hypot(ref_ip, ref_oop), atol=1e-9 - ) - - -@requires_fiber -class TestSampleOrientationMirrorPairs(unittest.TestCase): - """Document why the comparisons above use |q_ip|.""" - - def test_paired_orientations_share_qoop_and_mirror_qip(self): - det = _geometry(90.0) - rows, cols = _pixel_grid(det, stride=211) - alpha_i = np.deg2rad(0.6) - - for first, second in ((6, 7), (3, 4), (5, 8), (1, 2)): - with self.subTest(pair=(first, second)): - ip_a, oop_a = _q_from_pyfai(det, alpha_i, first, rows, cols) - ip_b, oop_b = _q_from_pyfai(det, alpha_i, second, rows, cols) - - np.testing.assert_allclose(oop_a, oop_b, atol=1e-12) - np.testing.assert_allclose(ip_a, -ip_b, atol=1e-12) - - -@requires_fiber -class TestQAlphaMatchesFiberIntegrator(unittest.TestCase): - """The rebinned FiberIntegrator map places intensity where QAlpha predicts.""" - - def test_single_pixel_peak_lands_at_the_qalpha_position(self): - angles = _vlieg_angles() - alpha_i = np.deg2rad(0.6) - - for azimuth_deg, orientation in AZIMUTH_ORIENTATION: - det = _geometry(azimuth_deg) - integrator = FiberIntegrator( - dist=det.dist, - poni1=det.poni1, - poni2=det.poni2, - wavelength=det.wavelength, - rot1=det.rot1, - rot2=det.rot2, - rot3=det.rot3, - detector=det.detector, - ) - - 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 - - # identical call to the one performed by the Q-plot action - result = integrator.integrate2d_grazing_incidence( - image, - sample_orientation=orientation, - incident_angle=alpha_i, - tilt_angle=0, - unit_oop="qoop_A^-1", - unit_ip="qip_A^-1", - ) - - 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]] - - q_ip, q_oop = _q_from_qalpha( - det, - angles, - alpha_i, - np.array([float(row)]), - np.array([float(col)]), - ) - - # 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) - q_ip[0]), bin_ip) - self.assertLess(abs(q_oop_binned - q_oop[0]), bin_oop) - - -if __name__ == "__main__": - unittest.main() From a18e3bbfcc53c3326214cba38cccb0c3bab39723 Mon Sep 17 00:00:00 2001 From: "T. Fuchs" Date: Wed, 29 Jul 2026 20:30:49 -0400 Subject: [PATCH 10/10] feat(app): show HKL and the selected Q frame in the plot readout --- CHANGELOG.md | 11 ++ doc/source/release_notes.rst | 4 + orgui/app/orGUI.py | 181 ++++++++++++++++++++++-- orgui/app/qconversion.py | 11 ++ orgui/app/test/test_position_readout.py | 181 ++++++++++++++++++++++++ 5 files changed, 373 insertions(+), 15 deletions(-) create mode 100644 orgui/app/test/test_position_readout.py diff --git a/CHANGELOG.md b/CHANGELOG.md index f16b011..d02c4e3 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -66,6 +66,17 @@ 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: - Cancelling or resetting the machine parameter or the crystal parameter diff --git a/doc/source/release_notes.rst b/doc/source/release_notes.rst index 0b3266c..c152202 100644 --- a/doc/source/release_notes.rst +++ b/doc/source/release_notes.rst @@ -28,6 +28,10 @@ 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: - 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. diff --git a/orgui/app/orGUI.py b/orgui/app/orGUI.py index 9c26a16..b4bdea8 100644 --- a/orgui/app/orGUI.py +++ b/orgui/app/orGUI.py @@ -57,6 +57,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 @@ -3135,30 +3136,54 @@ 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 detector pixels to hkl, angles and Q. + + :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. + """ if self.fscan is None: - return np.array([np.nan, np.nan, np.nan, np.nan, np.nan]) + return np.full(8, 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] 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]]), + self._qInSelectedFrame(pos), ) ) - return hkl return xyToHKL + 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. @@ -4117,8 +4142,9 @@ def _selectedQFrame(self): return frame def _onQFrameChanged(self, index): - """GUI hint: re-render the Q-plot after the frame selection changed.""" + """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): @@ -6190,19 +6216,88 @@ 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[") + + +# 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._qFieldName = qFieldName(qconversion.DEFAULT_FRAME) 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 + ("X", lambda x, y: f"{x:.2f}"), + ("Y", lambda x, y: f"{y:.2f}"), + ("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)), @@ -6228,6 +6323,8 @@ def __init__(self, xyHKLConverter, parent=None, backend=None): mask=True, ) + self._tunePositionFieldWidths() + if parent is None: self.setWindowTitle("Plot2D") self.getXAxis().setLabel("Columns") @@ -6257,6 +6354,60 @@ 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 setQFrame(self, frame): + """Relabel the reciprocal-space field for the selected frame. + + silx builds the field titles once, so the existing title widget is + renamed in place. The displayed values follow automatically, because + the converter reads the selected frame when it is called. + + .. note:: + GUI-only. + """ + name = qFieldName(frame) + if name == self._qFieldName: + return + + positionInfo = self.getPositionInfoWidget() + if positionInfo is not None: + previous = f"{self._qFieldName}:" + for widget in positionInfo.findChildren(qt.QLabel): + if widget.text() == previous: + widget.setText(f"{name}:") + break + positionInfo.updateInfo() + self._qFieldName = name + 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 index 15a08b3..091e7f3 100644 --- a/orgui/app/qconversion.py +++ b/orgui/app/qconversion.py @@ -152,6 +152,17 @@ def hasKernel(): 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. diff --git a/orgui/app/test/test_position_readout.py b/orgui/app/test/test_position_readout.py new file mode 100644 index 0000000..27a761a --- /dev/null +++ b/orgui/app/test/test_position_readout.py @@ -0,0 +1,181 @@ +"""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._qFieldName = "Q[alpha]" + self.titles = [_StubTitle("HKL:"), _StubTitle("Q[alpha]:")] + self.positionInfo = _StubPositionInfo(self.titles) + + def getPositionInfoWidget(self): + return self.positionInfo + + +@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[1].text() == f"{name}:" + assert plot.titles[0].text() == "HKL:" # untouched + 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[1].text() == "Q[alpha]:" + + +def test_relabelling_refreshes_the_displayed_values(): + plot = _StubPlot() + + orGUI.Plot2DHKL.setQFrame(plot, "Q_phi") + + assert plot.positionInfo.updates == 1