diff --git a/CHANGES.rst b/CHANGES.rst index a4705d059..c37d00a94 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -30,6 +30,11 @@ describe future plans. release expected by 2024-03-31 +New Features +------------ + +* Add new plan for edge alignment called edge_align + Fixes ----- diff --git a/apstools/plans/__init__.py b/apstools/plans/__init__.py index 1337677a5..7e6770c52 100644 --- a/apstools/plans/__init__.py +++ b/apstools/plans/__init__.py @@ -2,6 +2,7 @@ from .alignment import TuneResults from .alignment import lineup from .alignment import lineup2 +from .alignment import edge_align from .alignment import tune_axes from .command_list import CommandFileReadError from .command_list import command_list_as_table diff --git a/apstools/plans/alignment.py b/apstools/plans/alignment.py index ce5b19f71..aa599d227 100644 --- a/apstools/plans/alignment.py +++ b/apstools/plans/alignment.py @@ -16,6 +16,9 @@ import numpy as np import pyRestTable +from scipy.optimize import curve_fit +from scipy.special import erf + from bluesky import plan_stubs as bps from bluesky import plans as bp from bluesky import preprocessors as bpp @@ -25,6 +28,7 @@ from ophyd import Signal from ophyd.scaler import ScalerCH from ophyd.scaler import ScalerChannel +import warnings from .. import utils from .doc_run import write_stream @@ -34,8 +38,14 @@ def lineup( # fmt: off - detectors, axis, minus, plus, npts, - time_s=0.1, peak_factor=4, width_factor=0.8, + detectors, + axis, + minus, + plus, + npts, + time_s=0.1, + peak_factor=4, + width_factor=0.8, feature="cen", rescan=True, bec=None, @@ -135,7 +145,8 @@ def peak_analysis(): if det0.name not in bec.peaks[feature]: logger.error( - "No statistical analysis of scan peak for feature '%s'!" " (bec.peaks=%s, bec=%s)", + "No statistical analysis of scan peak for feature '%s'!" + " (bec.peaks=%s, bec=%s)", feature, bec.peaks, bec, @@ -213,10 +224,128 @@ def peak_analysis(): scaler.stage_sigs = old_sigs +def edge_align(detectors, mover, start, end, points, cat=None, md={}): + + """ + Align to the edge given mover & detector data, relative to absolute position. + + This plan can be used in the queueserver, Jupyter notebooks, and IPython + consoles. + + PARAMETERS + ---------- + detectors *Readable* or [*Readable*]: + Detector object or list of detector objects (each is a Device or + Signal). + + mover *Movable*: + Mover object, such as motor or other positioner. + + start *float*: + Starting point for the scan. This is an absolute motor location. + + end *float*: + Ending point for the scan. This is an absolute motor location. + + points *int*: + Number of points in the scan. + + cat *databroker.temp().v2*: + Catalog where bluesky data is saved and can be retrieved from. + + md *dict*: + User-supplied metadata for this scan. + """ + + def guess_erf_params(x_data, y_data): + """ + Provide an initial guess for the parameters of an error function. + + Parameters + ---------- + x_data : A numpy array of the values on the x_axis + y_data : A numpy array of the values on the y_axis + + Returns + ------- + guess : dict + A dictionary containing the guessed parameters 'low_y_data', 'high_y_data', 'width', and 'midpoint'. + """ + + # Sort data to make finding the mid-point easier and to assist in other estimations + y_data_sorted = np.sort(y_data) + x_data_sorted = np.sort(x_data) + + # Estimate low and high as the first and last elements (assuming sorted data) + low_y_data = np.min(y_data_sorted) + high_y_data = np.max(y_data_sorted) + + low_x_data = np.min(x_data_sorted) + high_x_data = np.max(x_data_sorted) + + # Estimate wid as a fraction of the range. This is very arbitrary and might need tuning! + width = ( + (high_x_data - low_x_data) / 10 + ) # This is a guess and might need adjustment based on your data's characteristics + + # Estimate the midpoint of the x values + midpoint = x_data[int(len(x_data) / 2)] + + return [low_y_data, high_y_data, width, midpoint] + + def erf_model(x, low, high, width, midpoint): + """ + Create error function for fitting and simulation + + Parameters + ---------- + x : input upon which error function is evaluated + low : min value of error function + high : max value of error function + width : "spread" of error function transition region + midpoint: location of error function's "center" + """ + return (high - low) * 0.5 * (1 - erf((x - midpoint) / width)) + low + + if not isinstance(detectors, (tuple, list)): + detectors = [detectors] + + _md = dict(purpose="edge_align") + _md.update(md or {}) + + uid = yield from bp.scan(detectors, mover, start, end, points, md=_md) + cat = cat or utils.getCatalog() + run = cat[uid] # return uids + ds = run.primary.read() + + x = ds["mover"] + y = ds["noisy"] + + initial_guess = guess_erf_params(x, y) + + try: + with warnings.catch_warnings(): + warnings.filterwarnings( + "error", "Covariance of the parameters could not be estimated" + ) + popt, pcov = curve_fit(erf_model, x, y, p0=initial_guess) + + yield from bps.mv(mover, popt[3]) + + except UserWarning as e: + print(f"Warning caught as exception: {e}") + print("\n Check to see if signal is linear.\n") + + def lineup2( # fmt: off - detectors, mover, rel_start, rel_end, points, - peak_factor=2.5, width_factor=0.8, + detectors, + mover, + rel_start, + rel_end, + points, + peak_factor=2.5, + width_factor=0.8, feature="centroid", nscans=2, signal_stats=None, @@ -621,7 +750,9 @@ def _scan(width=1, step_factor=10, num=10, snake=True): } _md.update(md or {}) - yield from self.tune(width=width, num=num, peak_factor=peak_factor, md=_md) + yield from self.tune( + width=width, num=num, peak_factor=peak_factor, md=_md + ) if not self.tune_ok: return @@ -633,7 +764,9 @@ def _scan(width=1, step_factor=10, num=10, snake=True): if snake: width *= -1 - return (yield from _scan(width=width, step_factor=step_factor, num=num, snake=snake)) + return ( + yield from _scan(width=width, step_factor=step_factor, num=num, snake=snake) + ) def multi_pass_tune_summary(self): t = pyRestTable.Table() @@ -724,7 +857,10 @@ class TuneResults(Device): peakstats_attrs = "x y cen com fwhm min max crossings".split() def report(self, title=None): - keys = self.peakstats_attrs + "tune_ok center initial_position final_position".split() + keys = ( + self.peakstats_attrs + + "tune_ok center initial_position final_position".split() + ) t = pyRestTable.Table() t.addLabel("key") t.addLabel("result") diff --git a/apstools/plans/edge_align.py b/apstools/plans/edge_align.py deleted file mode 100644 index cfa42096d..000000000 --- a/apstools/plans/edge_align.py +++ /dev/null @@ -1,67 +0,0 @@ -import numpy as np -from scipy.optimize import curve_fit -import matplotlib.pyplot as plt - -# Eample data -xpos = np.array([ - -0.89500, - -0.78900, - -0.68400, - -0.57900, - -0.47400, - -0.36800, - -0.26300, - -0.15800, - -0.05300, - 0.05300, - 0.15800, - 0.26300, - 0.36800, - 0.47400, - 0.57900, - 0.68400, - 0.78900, - 0.89500, - 1.00000, -]) - -sig = np.array([ - 5225.44680, - 10406.54446, - 29364.11991, - 91847.11330, - 39397.99591, - 12719.50183, - 5950.02145, - 3483.21205, - 2210.18858, - 1552.42419, - 1088.73029, - 882.58759, - 680.85617, - 562.82269, - 466.21100, - 385.87965, - 315.28561, - 284.60984, - 246.65506, -]) - - -dy_dx = np.gradient(sig, xpos) - -# Plotting -plt.figure(figsize=(10, 5)) - -# Plot the original curve -plt.subplot(1, 2, 1) -plt.plot(xpos, sig, label='Original curve') -plt.legend() - -# Plot the derivative -plt.subplot(1, 2, 2) -plt.plot(xpos, dy_dx, label='Derivative', color='red') -plt.legend() - -plt.tight_layout() -plt.show() \ No newline at end of file diff --git a/docs/source/api/_plans.rst b/docs/source/api/_plans.rst index 49fdd7cfc..47744806a 100644 --- a/docs/source/api/_plans.rst +++ b/docs/source/api/_plans.rst @@ -29,14 +29,15 @@ Custom Scans .. autosummary:: - ~apstools.plans.doc_run.documentation_run - ~apstools.plans.labels_to_streams.label_stream_decorator + ~apstools.plans.alignment.edge_align ~apstools.plans.alignment.lineup ~apstools.plans.alignment.lineup2 + ~apstools.plans.alignment.tune_axes + ~apstools.plans.alignment.TuneAxis + ~apstools.plans.doc_run.documentation_run + ~apstools.plans.labels_to_streams.label_stream_decorator ~apstools.plans.nscan_support.nscan ~apstools.plans.sscan_support.sscan_1D - ~apstools.plans.alignment.TuneAxis - ~apstools.plans.alignment.tune_axes .. _plans.overall: