Unverified Commit ab785083 authored by Coleman Kendrick's avatar Coleman Kendrick Committed by GitHub
Browse files

Merge pull request #30810 from mantidproject/pd188_difc_cal_plot

DIFC calibration diagnostic plot ornl-next
parents 823bac1a 44a6a6d4
......@@ -220,4 +220,76 @@ This can be done using the following code
Here the expected peak positions are vertical lines, the horizontal lines are boundaries between banks.
When run interactively, the zoom/pan tools are available.
DIFC of unwrapped instrument
############################
To check the consistency of pixel-level calibration, the DIFC value of each
pixel can be compared between two different instrument calibrations. The percent
change in DIFC value is plotted over a view of the unwrapped instrument where the
horizontal and vertical axis corresponds to the polar and azimuthal angle, respectively.
The azimuthal angle of 0 corresponds to the direction parallel of the positive Y-axis in
3D space.
Below is an example of the change in DIFC between two different calibrations of the
NOMAD instrument.
.. figure:: /images/NOMAD_difc_calibration.png
:width: 400px
This plot can be generated several different ways: by using calibration files,
calibration workspaces, or resulting workspaces from :ref:`CalculateDIFC <algm-CalculateDIFC>`.
The first input parameter is always required and represents the new calibration.
The second parameter is optional and represents the old calibration. When it is
not specified, the default instrument geometry is used for comparison. Masks can
be included by providing a mask using the ``mask`` parameter. To control the
scale of the plot, a tuple of the minimum and maximum percentage can be specified
for the ``vrange`` parameter.
.. code::
from Calibration.tofpd import diagnostics
# Use filenames to generate the plot
fig, ax = diagnostics.difc_plot2d("NOM_calibrate_d135279_2019_11_28.h5", "NOM_calibrate_d131573_2019_08_18.h5")
When calibration tables are used as inputs, an additional workspace parameter
is needed (``instr_ws``) to hold the instrument definition. This can be the GroupingWorkspace
generated with the calibration tables from :ref:`LoadDiffCal <algm-LoadDiffCal>` as seen below.
.. code::
from mantid.simpleapi import LoadDiffCal
from Calibration.tofpd import diagnostics
# Use calibration tables to generate the plot
LoadDiffCal(Filename="NOM_calibrate_d135279_2019_11_28.h5", WorkspaceName="new")
LoadDiffCal(Filename="NOM_calibrate_d131573_2019_08_18.h5", WorkspaceName="old")
fig, ax = diagnostics.difc_plot2d("new_cal", "old_cal", instr_ws="new_group")
Finally, workspaces with DIFC values can be used directly:
.. code::
from mantid.simpleapi import CalculateDIFC, LoadDiffCal
from Calibration.tofpd import diagnostics
# Use the results from CalculateDIFC directly
LoadDiffCal(Filename="NOM_calibrate_d135279_2019_11_28.h5", WorkspaceName="new")
LoadDiffCal(Filename="NOM_calibrate_d131573_2019_08_18.h5", WorkspaceName="old")
difc_new = CalculateDIFC(InputWorkspace="new_group", CalibrationWorkspace="new_cal")
difc_old = CalculateDIFC(InputWorkspace="old_group", CalibrationWorkspace="old_cal")
fig, ax = diagnostics.difc_plot2d(difc_new, difc_old)
A mask can also be applied with a ``MaskWorkspace`` to hide pixels from the plot:
.. code::
from mantid.simpleapi import LoadDiffCal
from Calibration.tofpd import diagnostics
# Use calibration tables to generate the plot
LoadDiffCal(Filename="NOM_calibrate_d135279_2019_11_28.h5", WorkspaceName="new")
LoadDiffCal(Filename="NOM_calibrate_d131573_2019_08_18.h5", WorkspaceName="old")
fig, ax = diagnostics.difc_plot2d("new_cal", "old_cal", instr_ws="new_group", mask="new_mask")
.. categories:: Calibration
......@@ -17,6 +17,7 @@ New features
- New algorithm :ref:`RebinRagged <algm-RebinRagged>` which can rebin a workspace with different binning parameters for each spectrum
- :ref:`PDCalibration <algm-PDCalibration>` now supports workspaces with grouped detectors (i.e. more than one detector per spectrum).
- New diagnostic plotting tool `Calibration.tofpd..diagnostics.plot2d` which adds markers for expected peak positions
- New diagnostic plotting tool `Calibration.tofpd.diagnostics.difc_plot2d` which plots the change in DIFC between two instrument calibrations.
Improvements
############
......
from mantid.plots.resampling_image.samplingimage import imshow_sampling
from mantid.plots.datafunctions import get_axes_labels
from mantid.simpleapi import mtd
from mantid.simpleapi import CalculateDIFC, LoadDiffCal, mtd
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
from matplotlib.collections import PatchCollection
import numpy as np
......@@ -35,6 +37,43 @@ def _get_xrange(wksp, xmarkers, tolerance, xmin, xmax):
return xmin, xmax
def _get_difc_ws(wksp, instr_ws=None):
if wksp is None:
return None
# Check if given a workspace
ws_str = str(wksp)
difc_ws = None
if not mtd.doesExist(ws_str):
# Check if it was a file instead
if ws_str.endswith(tuple([".h5", ".hd5", ".hdf", ".cal"])):
try:
LoadDiffCal(Filename=ws_str, WorkspaceName="__cal_{}".format(ws_str))
difc_ws = CalculateDIFC(InputWorkspace="__cal_{}_group".format(ws_str),
CalibrationWorkspace="__cal_{}_cal".format(ws_str),
OutputWorkspace="__difc_{}".format(ws_str))
except:
raise RuntimeError("Could not load calibration file {}".format(ws_str))
else:
raise RuntimeError("Could not find workspace {} in ADS and it was not a file".format(ws_str))
else:
# If workspace exists, check if it is a SpecialWorkspace2D (result from CalculateDIFC)
if mtd[ws_str].id() == "SpecialWorkspace2D":
difc_ws = mtd[ws_str]
elif mtd[ws_str].id() == "TableWorkspace":
if not mtd.doesExist(str(instr_ws)):
raise RuntimeError("Expected instrument workspace instr_ws to use with calibration tables")
# Check if the workspace looks like a calibration workspace
col_names = mtd[ws_str].getColumnNames()
# Only need the first two columns for the CalculateDIFC algorithm to work
if len(col_names) >= 2 and col_names[0] == "detid" and col_names[1] == "difc":
# Calculate DIFC on this workspace
difc_ws = CalculateDIFC(InputWorkspace=mtd[str(instr_ws)], CalibrationWorkspace=mtd[ws_str],
OutputWorkspace="__difc_{}".format(ws_str))
else:
raise TypeError("Wrong workspace type. Expects SpecialWorkspace2D, TableWorkspace, or a filename")
return difc_ws
def plot2d(workspace, tolerance: float=0.001, peakpositions: np.ndarray=DIAMOND,
xmin: float=np.nan, xmax: float=np.nan, horiz_markers=[]):
TOLERANCE_COLOR = 'w' # color to mark the area within tolerance
......@@ -78,3 +117,109 @@ def plot2d(workspace, tolerance: float=0.001, peakpositions: np.ndarray=DIAMOND,
# return the figure so others can customize it
return fig, fig.axes
def difc_plot2d(calib_new, calib_old=None, instr_ws=None, mask=None, vrange=(0,1)):
"""
Plots the percent change in DIFC between calib_new and calib_old
:param calib_new: New calibration, can be filename, SpecialWorkspace2D, or calibration table (TableWorkspace)
:param calib_old: Optional old calibration, can be same as types as calib_new; if not specified, default instr used
:param instr_ws: Workspace used for instrument definition, only needed if calib_new and/or calib_old are tables
:param mask: MaskWorkspace used to hide detector pixels from plot
:param vrange: Tuple of (min,max) used for the plot color bar
:return: figure and figure axes
"""
ws_new = _get_difc_ws(calib_new, instr_ws)
if ws_new is None:
raise TypeError("Expected to receive a workspace or filename, got None.")
ws_old = _get_difc_ws(calib_old, instr_ws)
if ws_old is None:
# If no second workspace is given, then load default instrument to compare against
instr_name = ws_new.getInstrument().getName()
ws_old = CalculateDIFC(InputWorkspace=ws_new, OutputWorkspace="__difc_{}".format(instr_name))
delta = ws_new - ws_old
use_mask = False
if mask is not None and mtd.doesExist(str(mask)):
use_mask = True
# Plotting below taken from addie/calibration/CalibrationDiagnostics.py
theta_array = []
phi_array = []
value_array = []
masked_theta_array = []
masked_phi_array = []
info = delta.spectrumInfo()
for det_id in range(info.size()):
pos = info.position(det_id)
phi = np.arctan2(pos[0], pos[1])
theta = np.arccos(pos[2] / pos.norm())
if use_mask and mtd[str(mask)].dataY(det_id):
masked_theta_array.append(theta)
masked_phi_array.append(phi)
else:
theta_array.append(theta)
phi_array.append(phi)
percent = 100.0 * np.sum(np.abs(delta.dataY(det_id))) / np.sum(ws_old.dataY(det_id))
value_array.append(percent)
# Use the largest solid angle for circle radius
sample_position = info.samplePosition()
maximum_solid_angle = 0.0
for det_id in range(info.size()):
maximum_solid_angle = max(maximum_solid_angle, delta.getDetector(det_id).solidAngle(sample_position))
# Convert to degrees for plotting
theta_array = np.rad2deg(theta_array)
phi_array = np.rad2deg(phi_array)
maximum_solid_angle = np.rad2deg(maximum_solid_angle)
if use_mask:
masked_phi_array = np.rad2deg(masked_phi_array)
masked_theta_array = np.rad2deg(masked_theta_array)
# Radius also includes a fudge factor to improve plotting.
# May need to add finer adjustments on a per-instrument basis.
# Small circles seem to alias less than rectangles.
radius = maximum_solid_angle * 8.0
patches = []
for x1, y1 in zip(theta_array, phi_array):
circle = Circle((x1, y1), radius)
patches.append(circle)
masked_patches = []
for x1, y1 in zip(masked_theta_array, masked_phi_array):
circle = Circle((x1, y1), radius)
masked_patches.append(circle)
# Matplotlib requires this to be a Numpy array.
colors = np.array(value_array)
p = PatchCollection(patches)
p.set_array(colors)
p.set_clim(vrange[0], vrange[1])
p.set_edgecolor('face')
fig, ax = plt.subplots()
ax.add_collection(p)
mp = PatchCollection(masked_patches)
mp.set_facecolor('gray')
mp.set_edgecolor('face')
ax.add_collection(mp)
cb = fig.colorbar(p, ax=ax)
cb.set_label('delta DIFC (%)')
ax.set_xlabel(r'in-plane polar (deg)')
ax.set_ylabel(r'azimuthal (deg)')
# Find bounds based on detector pos
xmin = np.min(theta_array)
xmax = np.max(theta_array)
ymin = np.min(phi_array)
ymax = np.max(phi_array)
ax.set_xlim(np.floor(xmin), np.ceil(xmax))
ax.set_ylim(np.floor(ymin), np.ceil(ymax))
fig.show()
return fig, ax
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment