diff --git a/knife_edge_scan.ipynb b/knife_edge_scan.ipynb new file mode 100644 index 0000000..c117fc2 --- /dev/null +++ b/knife_edge_scan.ipynb @@ -0,0 +1,247 @@ +{ + "cells": [ + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "import sys\n", + "sys.path.append('/home/leonarski_f/aaredaqlib')\n", + "sys.path.append('/home/leonarski_f/aaredaq')\n", + "from aaredaq.devices import BeamlineDevices\n", + "from aaredaqlib.beamline import MXBeamline" + ], + "id": "d2efa08de2dbd46b" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "from epics import PV\n", + "import numpy as np\n", + "import time\n", + "import matplotlib.pyplot as plt\n", + "from scipy.optimize import curve_fit\n", + "from scipy.signal import peak_widths, find_peaks\n", + "from scipy import asarray as ar,exp" + ], + "id": "f8305237c9f98964" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "d = BeamlineDevices(MXBeamline.X06DA)\n", + "b = MXBeamline.X06DA\n", + "beamline = b.value.upper()" + ], + "id": "5407e46ef4c2a26d" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "def gaus(x, *p):\n", + " A, mu, sigma = p\n", + " return A * np.exp(-(x-mu)**2/(2.*sigma**2))\n", + "\n", + "def sigmoid(x, *p):\n", + " L, x0, k, b = p\n", + " return L / (1 + np.exp(-k * (x -x0))) + b\n", + "\n", + "def unpack_data(data):\n", + " x_vals, y_vals = zip(*data)\n", + " x_vals = list(x_vals)\n", + " y_vals = list(y_vals)\n", + " return x_vals, y_vals\n", + "\n", + "def sigmoid_fit(x_vals, y_vals):\n", + " L_guess = np.max(y_vals)\n", + " b_guess = np.min(y_vals)\n", + " x0_guess = x_vals[np.argmin(np.abs(y_vals - (L_guess + b_guess) / 2))]\n", + " k_guess = 1 / (x_vals[-1] - x_vals[0])\n", + "\n", + " s0 = [L_guess, x0_guess, k_guess, b_guess]\n", + " parameters, params_covariance = curve_fit(sigmoid, x_vals, y_vals, p0=s0)\n", + " return parameters\n", + "\n", + "def gauss_fit(x_vals, dydx):\n", + " mean = x_vals[np.argmax(dydx)]\n", + " sigma = (np.max(x_vals) - np.min(x_vals))/4\n", + " A= np.max(dydx)\n", + " dydx[np.isnan(dydx)] = 0\n", + "\n", + " p0 = [A, mean, sigma]\n", + " parameters, covariance = curve_fit( gaus, x_vals, dydx, p0=p0)\n", + " return parameters\n", + "\n", + "def calc_fwhm(gauss_y, x_vals):\n", + " peaks = find_peaks( gauss_y )\n", + " fwhm_peak = peak_widths( gauss_y, peaks[0], rel_height=0.5 )\n", + "\n", + " fwhm_str = x_vals[int( round( fwhm_peak[2][0], 0 ))]\n", + " fwhm_end = x_vals[int( round( fwhm_peak[3][0], 0 ))]\n", + "\n", + " fwhm = fwhm_str - fwhm_end\n", + " fwhm_height = gauss_y[int( round( fwhm_peak[2][0], 0 ))]\n", + "\n", + " # find 1/e2 peak width\n", + " full_peak = peak_widths( gauss_y, peaks[0], rel_height=0.865 )\n", + " full_str = x_vals[int( round( full_peak[2][0], 0 ))]\n", + " full_end = x_vals[int( round( full_peak[3][0], 0 ))]\n", + "\n", + " full = full_str - full_end\n", + " full_height = gauss_y[int( round( full_peak[2][0], 0 ))]\n", + " return fwhm, full" + ], + "id": "d92115ec3d76d57d" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "def knife_edge(motor: str, d: BeamlineDevices, beamline: str, steps: int = 40, start: float = 0.0, step_size: float = 0.002, offset: float = 0.0):\n", + " diode = PV(f\"{beamline}-ES-BS:READOUT\")#X06DA-ES-BS:READOUT\n", + " d.shutter = True\n", + " time.sleep(1.0)\n", + "\n", + " data = []\n", + "\n", + " for i in range(steps): #40\n", + " if motor == 'X':\n", + " d.gmx.move(start-offset+step_size*i, wait=True)\n", + " time.sleep(1.0)\n", + " print(f\"{d.gmx.value:.3f} {diode.get()}\")\n", + " data.append([d.gmx.value, diode.get()])\n", + " else:\n", + " d.gmy.move(start-offset+step_size*i, wait=True)\n", + " time.sleep(1.0)\n", + " print(f\"{d.gmy.value:.3f} {diode.get()}\")\n", + " data.append([d.gmy.value, diode.get()])\n", + "\n", + " if motor == 'X':\n", + " d.gmx.move(start, wait=True)\n", + " else:\n", + " d.gmy.move(start, wait=True)\n", + "\n", + " d.shutter = False\n", + " return data\n" + ], + "id": "512dc280e62e4ed4" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "def plot_knife_edge(data, motor: str):\n", + "\n", + " x_vals, y_vals = unpack_data(data)\n", + "\n", + " plt.plot(x_vals, y_vals, label = 'edge_scan')\n", + " plt.ylabel('Intensity (counts)')\n", + " plt.xlabel('Motor position (mm)')\n", + "\n", + " params=sigmoid_fit(x_vals, y_vals)\n", + "\n", + " plt.plot(x_vals, y_vals)\n", + " plt.plot(x_vals, sigmoid(x_vals, *params))\n", + " plt.show()\n", + "\n", + " y_fit = sigmoid(x_vals, *params)\n", + "\n", + " dydx = np.gradient(y_fit, x_vals)\n", + "\n", + " if motor == 'Y':\n", + " dydx=-dydx\n", + "\n", + " parameters = gauss_fit(x_vals, dydx)\n", + "\n", + " x0_op = parameters[1]\n", + " sigma_op = parameters[2]\n", + " print(parameters)\n", + "\n", + " gauss_y = gaus(x_vals,*parameters)\n", + " fwhm = np.abs(2*np.sqrt(2*np.log(2))*sigma_op)\n", + " print(f\"Manually calculated FWHM to {fwhm} mm\")\n", + " plt.plot(x_vals, dydx, label = 'derivative')\n", + " plt.plot(x_vals, gauss_y,label='Gaussian fit',color ='orange')\n", + " plt.fill_between(x_vals,gauss_y,color='orange',alpha=0.5)\n", + "\n", + " try:\n", + " fwhm, full = calc_fwhm(gauss_y, x_vals)\n", + " print( \"Scipy calculated FWHM = {0} mm\".format( fwhm ) )\n", + " print( \"Scipy calculated 1/e2 = {0} mm\".format( full ) )\n", + " except:\n", + " print('scipy peak finding failed')\n", + "\n", + " plt.axvspan(x0_op+fwhm/2,x0_op-fwhm/2, color='green', alpha=0.75, lw=0, label='FWHM = {0} mm'.format(fwhm))\n", + " plt.legend()\n", + " plt.show()\n", + " #plt.savefig(output_name+filename.split('/')[-1]+'FWHM_{0}.png'.format(fwhm))" + ], + "id": "80460b35a33f1ecf" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "x_data=knife_edge(motor='X', d=d, beamline=beamline, steps=40, start=-17.834, step_size=0.002, offset=0.006)\n", + "plot_knife_edge(x_data, motor='X')" + ], + "id": "5454eb901526eb3" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": [ + "y_data=knife_edge(motor='Y', d=d, beamline=beamline, steps=40, start=0.0, step_size=0.002, offset=0.03)\n", + "plot_knife_edge(y_data, motor='X')" + ], + "id": "f8e0541c9e92e4e6" + }, + { + "metadata": {}, + "cell_type": "code", + "outputs": [], + "execution_count": null, + "source": "", + "id": "d6dd221d89d3b71a" + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 2 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython2", + "version": "2.7.6" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}