Calibration of the 9-Mythen detector at the Cristal beamline at Soleil
Mythen detectors are 1D-strip detector sold by Dectris. On the Cristal beamline at Soleil, 9 of them are mounted on the goniometer.
This notebook explains how to calibrate precisely their position (including the wavelength used) as function of the goniometer position.
All input data are provided in a Nexus file wich contrains both the (approximate) energy, the goniometer positions (500 points have been measured) and the measured signal.
As pyFAI is not made for 1D data, the Mythen detector will be considered as a 1x1280 image.
We start by importing a whole bunch of modules:
[1]:
%matplotlib inline
# use `widget` for better user experience; `inline` is for documentation generation
[2]:
from collections import OrderedDict
from matplotlib import pyplot as plt
import numpy
import os
import h5py
from silx.resources import ExternalResources
from pyFAI import goniometer
from pyFAI.detectors import Detector
from pyFAI.goniometer import ExtendedTransformation, GoniometerRefinement
from pyFAI.control_points import ControlPoints
from pyFAI.geometryRefinement import GeometryRefinement
from pyFAI.gui import jupyter
from pyFAI.units import hc
from pyFAI.calibrant import get_calibrant
from pyFAI.containers import Integrate1dResult
import ipywidgets as widgets
from scipy.signal import find_peaks_cwt
from scipy.interpolate import interp1d
from scipy.optimize import bisect, minimize
from scipy.spatial import distance_matrix
import time
start_time = time.time()
[3]:
#Nota: Useful to configure a proxy if you are behind a firewall
#os.environ["http_proxy"] = "http://proxy.company.fr:3128"
downloader = ExternalResources("detector_calibration", "http://www.silx.org/pub/pyFAI/gonio/")
mythen_ring_file = downloader.getfile("LaB6_17keV_att3_tth2C_24_01_2018_19-43-20_1555.nxs")
The data file can be downoaded from: http://www.silx.org/pub/pyFAI/gonio/LaB6_17keV_att3_tth2C_24_01_2018_19-43-20_1555.nxs
[4]:
#Open the Nexus file and retrieve the actual positions:
h5 = h5py.File(mythen_ring_file, mode="r")
position = h5["/LaB6_17keV_att3_1555/scan_data/actuator_1_1"][:]
print("Positions: ", position[:5], "...")
Positions: [90.00000001 89.79994445 89.5998889 89.39994445 89.19994445] ...
[5]:
#Read all data
data = {}
ds_names = []
for idx in range(1,13):
name = "data_%02i"%idx
ds = h5["/LaB6_17keV_att3_1555/scan_data/"+name][:]
print(name, ds.shape)
if ds.shape[1]<2000:
#Keep only the single modules
data[name] = ds
ds_names.append(name)
print(ds_names)
data_01 (501, 5120)
data_02 (501, 1280)
data_03 (501, 1280)
data_04 (501, 1280)
data_05 (501, 1280)
data_06 (501, 5120)
data_07 (501, 1280)
data_08 (501, 1280)
data_09 (501, 1280)
data_10 (501, 1280)
data_11 (501, 1280)
data_12 (501, 1280)
['data_02', 'data_03', 'data_04', 'data_05', 'data_07', 'data_08', 'data_09', 'data_10', 'data_11', 'data_12']
[6]:
#Define a Mythen-detector mounted vertically:
class MythenV(Detector):
"Verical Mythen dtrip detector from Dectris"
aliases = ["MythenV 1280"]
force_pixel = True
MAX_SHAPE = (1280, 1)
def __init__(self,pixel1=50e-6, pixel2=8e-3):
super(MythenV, self).__init__(pixel1=pixel1, pixel2=pixel2)
[7]:
#Define all modules as single detectors of class MythenV.
# Each one has a mask defined from dummy-values in the dataset
modules = {}
for name, ds in data.items():
one_module = MythenV()
mask = ds[0]<0
#discard the first 20 and last 20 pixels as their intensities are less reliable
mask[:20] = True
mask[-20:] = True
one_module.mask = mask.reshape(-1,1)
modules[name] = one_module
for k,v in modules.items():
print(k, v.name)
data_02 MythenV 1280
data_03 MythenV 1280
data_04 MythenV 1280
data_05 MythenV 1280
data_07 MythenV 1280
data_08 MythenV 1280
data_09 MythenV 1280
data_10 MythenV 1280
data_11 MythenV 1280
data_12 MythenV 1280
[8]:
# Define a peak-picking function based on the dataset-name and the frame_id:
def peak_picking(module_name, frame_id,
threshold=500):
"""Peak-picking base on find_peaks_cwt from scipy plus
second-order tailor exapention refinement for sub-pixel resolution.
The half-pixel offset is accounted here, i.e pixel #0 has its center at 0.5
"""
module = modules[module_name]
msk = module.mask.ravel()
spectrum = data[module_name][frame_id]
guess = find_peaks_cwt(spectrum, [20])
valid = numpy.logical_and(numpy.logical_not(msk[guess]),
spectrum[guess]>threshold)
guess = guess[valid]
#Based on maximum is f'(x) = 0 ~ f'(x0) + (x-x0)*(f''(x0))
df = numpy.gradient(spectrum)
d2f = numpy.gradient(df)
bad = d2f==0
d2f[bad] = 1e-10 #prevent devision by zero. Discared later on
cor = df / d2f
cor[abs(cor)>1] = 0
cor[bad] = 0
ref = guess - cor[guess] + 0.5 #half a pixel offset
x = numpy.zeros_like(ref) + 0.5 #half a pixel offset
return numpy.vstack((ref,x)).T
%timeit peak_picking(ds_names[0], 93)
print(peak_picking(ds_names[0], 93))
15.9 ms ± 202 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
[[287.06072343 0.5 ]]
[9]:
nrj = h5["/LaB6_17keV_att3_1555/CRISTAL/Monochromator/energy"][0]
wl = hc / nrj *1e-10
print("Energy (keV): ",nrj, "\nWavelength (A): ",wl)
LaB6 = get_calibrant("LaB6")
LaB6.wavelength = wl
print(LaB6)
Energy (keV): 17.027082549190933
Wavelength (A): 7.281587910025816e-11
LaB6 Calibrant with 109 reflections at wavelength 7.281587910025816e-11
[10]:
#This cell defines the transformation of coordinates for a simple goniometer mounted vertically.
trans = ExtendedTransformation(dist_expr="dist",
poni1_expr="poni1",
poni2_expr="poni2",
rot1_expr="rot1",
rot2_expr="pi*(offset+scale*angle)/180.",
rot3_expr="0.0",
wavelength_expr="hc/nrj*1e-10",
param_names=["dist", "poni1", "poni2", "rot1", "offset", "scale", "nrj"],
pos_names=["angle"],
constants={"hc": hc})
[11]:
def get_position(idx):
"Returns the postion of the goniometer for the given frame_id"
return position[idx]
#Approximate offset for the module #0 at 0°
print("Approximated offset for the first module: ",get_position(36))
Approximated offset for the first module: 82.79994445106844
[12]:
#This interactive plot lets one visualize any spectra acquired by any module
fig, ax = plt.subplots()
line = ax.plot(data[ds_names[0]][250])[0]
ligne = plt.Line2D(xdata=[640,640], ydata=[-500, 1000], figure=fig, linestyle="--", color='red', axes=ax)
ax.add_line(ligne)
ax.set_title("spectrum")
fig.show()
def update(module_id, frame_id):
spectrum = data[ds_names[module_id]][frame_id]
line.set_data(numpy.arange(spectrum.size), spectrum)
ax.set_title("Module %i, Frame %i"%(module_id, frame_id))
fig.canvas.draw()
interactive_plot = widgets.interactive(update,
module_id=(0, len(data)-1),
frame_id=(0, data[ds_names[0]].shape[0]-1))
display(interactive_plot)
<ipython-input-12-53c6060595b4>:8: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[13]:
#Work with the first module corresponding to:
name = ds_names[0]
print(name)
ds = data[name]
module = modules[name]
#Use the previous widget to select:
## the index where the beam-center is in the middle of the module
zero_pos = 36
## The frame index where the first LaB6 peak enters the right-hand side of the spectrum
peak_zero_start = 74
## The frame index where this first LaB6 leaves the spectrum or the second LaB6 peak appears:
peak_zero_end = 94
# The frames between peak_zero_start and peak_zero_end will be used to calibrate roughly the goniometer
# and used later for finer peak extraction
data_02
[14]:
param0 = {"dist": 0.72,
"poni1": 640*50e-6,
"poni2": 4e-3,
"rot1":0,
"offset": -get_position(zero_pos),
"scale":1,
"nrj": nrj}
#Lock enegy for now and a couple of other parameters
bounds0 = {"nrj": (nrj, nrj),
"dist": (0.71, 0.73),
"poni2": (4e-3, 4e-3),
"rot1": (0,0),
"scale":(1,1),
}
gonioref0 = GoniometerRefinement(param0,
get_position,
trans,
detector=module,
wavelength=wl,
bounds=bounds0
)
goniometers = {name: gonioref0}
print(gonioref0)
GoniometerRefinement with 0 geometries labeled: .
[15]:
# Extract the frames where only the peak zero from LaB6 is present.
for i in range(peak_zero_start, peak_zero_end):
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
peak = peak_picking(name, i)
if len(peak)!=1:
continue
cp.append([peak[0]], ring=0)
img = ds[i].reshape((-1,1)) #Images are vertical ... transpose the spectrum
sg = gonioref0.new_geometry("%s_%04i"%(name,i),
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
print(gonioref0)
print("Residual error before fit:")
print(gonioref0.chi2())
GoniometerRefinement with 20 geometries labeled: data_02_0074, data_02_0075, data_02_0076, data_02_0077, data_02_0078, data_02_0079, data_02_0080, data_02_0081, data_02_0082, data_02_0083, data_02_0084, data_02_0085, data_02_0086, data_02_0087, data_02_0088, data_02_0089, data_02_0090, data_02_0091, data_02_0092, data_02_0093.
Residual error before fit:
6.737408475887519e-07
[16]:
#First refinement:
gonioref0.refine2()
Cost function before refinement: 6.737408475887519e-07
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-8.27999445e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.490988846910722e-11
x: [ 7.200e-01 3.141e-02 4.000e-03 0.000e+00 -8.280e+01
1.000e+00 1.703e+01]
nit: 2
jac: [ 1.584e-07 2.919e-08 nan nan 3.427e-11
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 2.490988846910722e-11
GonioParam(dist=0.719994724358983, poni1=0.031408577292064206, poni2=0.004, rot1=0.0, offset=-82.79995188659902, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.031408577292064206
[16]:
array([ 7.19994724e-01, 3.14085773e-02, 4.00000000e-03, 0.00000000e+00,
-8.27999519e+01, 1.00000000e+00, 1.70270825e+01])
[17]:
#Here we extract all spectra for peaks,
# If there are as many peaks as expected from the theoritical LaB6. perform the assignment.
#Peaks from LaB6:
tths = LaB6.get_2th()
for i in range(peak_zero_end, ds.shape[0]):
peak = peak_picking(name, i)
ai=gonioref0.get_ai(get_position(i))
tth = ai.array_from_unit(unit="2th_rad", scale=False)
tth_low = tth[20]
tth_hi = tth[-20]
ttmin, ttmax = min(tth_low, tth_hi), max(tth_low, tth_hi)
valid_peaks = numpy.logical_and(ttmin<=tths, tths<ttmax)
cnt = valid_peaks.sum()
if (len(peak) == cnt):
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
#revert the order of assignment if needed !!
if tth_hi < tth_low:
peak = peak[-1::-1]
for p, r in zip(peak, numpy.where(valid_peaks)[0]):
#print(p,r)
cp.append([p], ring=r)
img = ds[i].reshape((-1,1))
sg = gonioref0.new_geometry("%s_%04i"%(name,i),
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
#print(sg.label, len(sg.geometry_refinement.data))
#print(gonioref0)
print(" Number of peaks found and used for refinement")
print(sum([len(sg.geometry_refinement.data) for sg in gonioref0.single_geometries.values()]))
print("Residual error before fitting: ", gonioref0.chi2())
Number of peaks found and used for refinement
1203
Residual error before fitting: 3.118662178645908e-06
[18]:
gonioref0.refine2()
Cost function before refinement: 3.118662178645908e-06
[ 7.19994724e-01 3.14085773e-02 4.00000000e-03 0.00000000e+00
-8.27999519e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.7342490266360614e-06
x: [ 7.231e-01 3.185e-02 4.000e-03 0.000e+00 -8.280e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-5.582e-08 -3.614e-10 nan nan 6.362e-10
nan nan]
nfev: 30
njev: 7
Cost function after refinement: 2.7342490266360614e-06
GonioParam(dist=0.7230953424558498, poni1=0.031845370299442156, poni2=0.004, rot1=0.0, offset=-82.79994661585035, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.719994724358983 --> 0.7230953424558498
[18]:
array([ 7.23095342e-01, 3.18453703e-02, 4.00000000e-03, 0.00000000e+00,
-8.27999466e+01, 1.00000000e+00, 1.70270825e+01])
[19]:
gonioref0.set_bounds("poni1", -1, 1)
gonioref0.set_bounds("poni2", -1, 1)
gonioref0.set_bounds("rot1", -1, 1)
gonioref0.set_bounds("scale", 0.9, 1.1)
gonioref0.refine2()
Cost function before refinement: 2.7342490266360614e-06
[ 7.23095342e-01 3.18453703e-02 4.00000000e-03 0.00000000e+00
-8.27999466e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.689172533753603e-06
x: [ 7.231e-01 3.205e-02 3.985e-03 1.095e-05 -8.280e+01
9.994e-01 1.703e+01]
nit: 4
jac: [-5.145e-07 -1.241e-07 5.957e-07 -4.287e-07 -9.292e-10
-8.539e-08 nan]
nfev: 28
njev: 4
Cost function after refinement: 2.689172533753603e-06
GonioParam(dist=0.7230962441248261, poni1=0.032053700795884155, poni2=0.003984854880653816, rot1=1.0947770172731982e-05, offset=-82.79994398859087, scale=0.9994150576472143, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9994150576472143
[19]:
array([ 7.23096244e-01, 3.20537008e-02, 3.98485488e-03, 1.09477702e-05,
-8.27999440e+01, 9.99415058e-01, 1.70270825e+01])
[20]:
# Perform the azimuthal intgration of all data for the first module:
mg = gonioref0.get_mg(position)
mg.radial_range = (0, 95)
images = [i.reshape(-1, 1) for i in ds]
res_mg = mg.integrate1d(images, 50000)
results={name: res_mg}
print(results)
{'data_02': (array([9.50000113e-04, 2.85000034e-03, 4.75000057e-03, ...,
9.49952613e+01, 9.49971613e+01, 9.49990613e+01]), array([0.00000000e+00, 0.00000000e+00, 0.00000000e+00, ...,
1.27158634e+08, 1.31904151e+08, 1.31268034e+08]))}
[21]:
# Plot the integrated pattern vs expected peak positions:
LaB6_new = get_calibrant("LaB6")
LaB6_new.wavelength = hc/gonioref0.param[-1]*1e-10
p = jupyter.plot1d(res_mg, calibrant=LaB6_new)
p.figure.show()
<ipython-input-21-a40be3200594>:6: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
[22]:
#Peak profile function based on a bilinear interpolations:
def calc_fwhm(integrate_result, calibrant, tth_min=None, tth_max=None):
"calculate the tth position and FWHM for each peak"
delta = integrate_result.intensity[1:] - integrate_result.intensity[:-1]
maxima = numpy.where(numpy.logical_and(delta[:-1]>0, delta[1:]<0))[0]
minima = numpy.where(numpy.logical_and(delta[:-1]<0, delta[1:]>0))[0]
maxima += 1
minima += 1
tth = []
FWHM = []
if tth_min is None:
tth_min = integrate_result.radial[0]
if tth_max is None:
tth_max = integrate_result.radial[-1]
for tth_rad in calibrant.get_2th():
tth_deg = tth_rad*integrate_result.unit.scale
if (tth_deg<=tth_min) or (tth_deg>=tth_max):
continue
idx_theo = abs(integrate_result.radial-tth_deg).argmin()
id0_max = abs(maxima-idx_theo).argmin()
id0_min = abs(minima-idx_theo).argmin()
I_max = integrate_result.intensity[maxima[id0_max]]
I_min = integrate_result.intensity[minima[id0_min]]
tth_maxi = integrate_result.radial[maxima[id0_max]]
I_thres = (I_max + I_min)/2.0
if minima[id0_min]>maxima[id0_max]:
if id0_min == 0:
min_lo = integrate_result.radial[0]
else:
min_lo = integrate_result.radial[minima[id0_min-1]]
min_hi = integrate_result.radial[minima[id0_min]]
else:
if id0_min == len(minima) -1:
min_hi = integrate_result.radial[-1]
else:
min_hi = integrate_result.radial[minima[id0_min+1]]
min_lo = integrate_result.radial[minima[id0_min]]
f = interp1d(integrate_result.radial, integrate_result.intensity-I_thres)
try:
tth_lo = bisect(f, min_lo, tth_maxi)
tth_hi = bisect(f, tth_maxi, min_hi)
except:
pass
else:
FWHM.append(tth_hi-tth_lo)
tth.append(tth_deg)
return tth, FWHM
[23]:
# Peak error:
def calc_peak_error(integrate_result, calibrant, tth_min=10, tth_max=95):
"calculate the tth position and FWHM for each peak"
peaks = find_peaks_cwt(integrate_result.intensity, [10])
df = numpy.gradient(integrate_result.intensity)
d2f = numpy.gradient(df)
bad = d2f==0
d2f[bad] = 1e-10
cor = df / d2f
print((abs(cor)>1).sum())
cor[abs(cor)>1] = 0
cor[bad] = 0
got = numpy.interp(peaks-cor[peaks],
numpy.arange(len(integrate_result.radial)),
integrate_result.radial)
mask = numpy.logical_and(got>=tth_min,
got<=tth_max)
got = got[mask]
target = numpy.array(calibrant.get_2th())*integrate_result.unit.scale
mask = numpy.logical_and(target>=tth_min,
target<=tth_max)
target = target[mask]
print(len(got), len(target))
d2 = distance_matrix(target.reshape(-1, 1 ),
got.reshape(-1, 1), p=1)
return target, target-got[d2.argmin(axis=-1)]
[24]:
fig, ax = plt.subplots()
ax.plot(*calc_fwhm(res_mg, LaB6_new), "o", label="FWHM")
ax.plot(*calc_peak_error(res_mg, LaB6_new), "o", label="offset")
ax.set_title("Peak shape & error as function of the angle")
ax.set_xlabel(res_mg.unit.label)
ax.legend()
fig.show()
29218
81 60
<ipython-input-24-b253c38f8aaa>:7: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
Module 1
We can apply the same procdure for the second module … and try to rationalize the procedure.
[25]:
module_id = 1
name = ds_names[module_id]
ds = data[name]
zero_pos = 64
frame_start = 103
frame_stop = 123
[26]:
param1 = {"dist": 0.72,
"poni1": 640*50e-6,
"poni2": 4e-3,
"rot1":0,
"offset": -get_position(zero_pos),
"scale":1,
"nrj": nrj}
#Lock enegy for now and a couple of other parameters
bounds1 = {"nrj": (nrj, nrj),
"dist": (0.7, 0.8),
"poni2": (4e-3, 4e-3),
"rot1": (0,0),
"scale":(1,1), }
gonioref1 = GoniometerRefinement(param1,
get_position,
trans,
detector=modules[name],
wavelength=wl,
bounds=bounds1
)
print(gonioref1)
goniometers[name]=gonioref1
GoniometerRefinement with 0 geometries labeled: .
[27]:
#Exctract frames with peak#0
for i in range(frame_start, frame_stop):
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
peak = peak_picking(name, i)
if len(peak)!=1:
continue
cp.append([peak[0]], ring=0)
img = (ds[i]).reshape((-1,1))
sg = gonioref1.new_geometry("%s_%04i"%(name,i),
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
print(gonioref1)
print(gonioref1.chi2())
GoniometerRefinement with 20 geometries labeled: data_03_0103, data_03_0104, data_03_0105, data_03_0106, data_03_0107, data_03_0108, data_03_0109, data_03_0110, data_03_0111, data_03_0112, data_03_0113, data_03_0114, data_03_0115, data_03_0116, data_03_0117, data_03_0118, data_03_0119, data_03_0120, data_03_0121, data_03_0122.
1.4524664758918898e-06
[28]:
gonioref1.refine2()
Cost function before refinement: 1.4524664758918898e-06
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-7.72000000e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.343186419587841e-11
x: [ 7.200e-01 3.287e-02 4.000e-03 0.000e+00 -7.720e+01
1.000e+00 1.703e+01]
nit: 2
jac: [ 1.374e-07 1.035e-07 nan nan 9.897e-10
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 2.343186419587841e-11
GonioParam(dist=0.7200063780234879, poni1=0.03286839545435715, poni2=0.004, rot1=0.0, offset=-77.19998908851599, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.03286839545435715
[28]:
array([ 7.20006378e-01, 3.28683955e-02, 4.00000000e-03, 0.00000000e+00,
-7.71999891e+01, 1.00000000e+00, 1.70270825e+01])
[29]:
#Exctract all frames with peak>0
tths = LaB6.get_2th()
#print(tths)
for i in range(frame_stop, ds.shape[0]):
frame_name = "%s_%04i"%(name, i)
if frame_name in gonioref1.single_geometries:
continue
peak = peak_picking(name, i)
ai=gonioref1.get_ai(get_position(i))
tth = ai.array_from_unit(unit="2th_rad", scale=False)
tth_low = tth[20]
tth_hi = tth[-20]
ttmin, ttmax = min(tth_low, tth_hi), max(tth_low, tth_hi)
valid_peaks = numpy.logical_and(ttmin<=tths, tths<ttmax)
cnt = valid_peaks.sum()
if (len(peak) == cnt) and cnt>0:
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
#revert the order of assignment if needed !!
if tth_hi < tth_low:
peak = peak[-1::-1]
for p, r in zip(peak, numpy.where(valid_peaks)[0]):
cp.append([p], ring=r)
img = ds[i].reshape((-1,1))
sg = gonioref1.new_geometry(frame_name,
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
#print(frame_name, len(sg.geometry_refinement.data))
print(" Number of peaks found and used for refinement")
print(sum([len(sg.geometry_refinement.data) for sg in gonioref1.single_geometries.values()]))
print("Residual error before fitting: ", gonioref1.chi2())
Number of peaks found and used for refinement
1183
Residual error before fitting: 6.334637851253146e-07
[30]:
gonioref1.refine2()
gonioref1.set_bounds("poni1", -1, 1)
gonioref1.set_bounds("poni2", -1, 1)
gonioref1.set_bounds("rot1", -1, 1)
gonioref1.set_bounds("scale", 0.9, 1.1)
gonioref1.refine2()
Cost function before refinement: 6.334637851253146e-07
[ 7.20006378e-01 3.28683955e-02 4.00000000e-03 0.00000000e+00
-7.71999891e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.3797258269103233e-07
x: [ 7.200e-01 3.338e-02 4.000e-03 0.000e+00 -7.720e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-4.626e-07 2.408e-08 nan nan 1.426e-11
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 1.3797258269103233e-07
GonioParam(dist=0.7200063394923231, poni1=0.033375469747964855, poni2=0.004, rot1=0.0, offset=-77.19998271240225, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.03286839545435715 --> 0.033375469747964855
Cost function before refinement: 1.3797258269103233e-07
[ 7.20006339e-01 3.33754697e-02 4.00000000e-03 0.00000000e+00
-7.71999827e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 8.82850324363492e-10
x: [ 7.207e-01 3.368e-02 4.058e-03 -4.449e-05 -7.720e+01
9.990e-01 1.703e+01]
nit: 13
jac: [ 6.561e-08 6.554e-07 -1.755e-07 1.262e-07 7.935e-09
2.914e-07 nan]
nfev: 91
njev: 13
Cost function after refinement: 8.82850324363492e-10
GonioParam(dist=0.7206820726599381, poni1=0.0336782330159447, poni2=0.004058010681093319, rot1=-4.44855658960725e-05, offset=-77.19997879322948, scale=0.9989672333400683, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9989672333400683
[30]:
array([ 7.20682073e-01, 3.36782330e-02, 4.05801068e-03, -4.44855659e-05,
-7.71999788e+01, 9.98967233e-01, 1.70270825e+01])
[31]:
mg1 = gonioref1.get_mg(position)
mg1.radial_range = (0, 95)
images = [i.reshape(-1, 1) for i in data[name]]
res_mg1 = mg1.integrate1d(images, 50000)
results[name] = res_mg1
[32]:
LaB6_new = get_calibrant("LaB6")
LaB6_new.wavelength = hc/gonioref1.param[-1]*1e-10
p = jupyter.plot1d(res_mg1, calibrant=LaB6_new)
p.figure.show()
<ipython-input-32-c0bc613f7ca9>:4: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
[33]:
fig, ax = plt.subplots()
ax.plot(*calc_fwhm(res_mg1, LaB6_new, 10, 88), "o", label="FWHM")
ax.plot(*calc_peak_error(res_mg1, LaB6_new, 10, 88), "o", label="error")
ax.set_title("Peak shape & error as function of the angle")
ax.set_xlabel(res_mg.unit.label)
ax.legend()
fig.show()
27712
72 53
<ipython-input-33-8212dfc61f3a>:7: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
All other Modules
We define now an automatic procedure for any module. The detection used 3 parameter visually extracted from the Figure1:
zero_pos: the frame where the beam-stop is in the center of the module
frame_start: the frame where the first peak of LaB6 appears (positive)
frame_stop: the frame where the second peak of LaB6 appears (positive)
This is enough for boot-strapping the goniometer configuration.
[34]:
def add_module(name,
zero_pos,
frame_start,
frame_stop,
):
ds = data[name]
param = {"dist": 0.72,
"poni1": 640*50e-6,
"poni2": 4e-3,
"rot1":0,
"offset": -get_position(zero_pos),
"scale":1,
"nrj": nrj}
#Lock enegy for now and a couple of other parameters
bounds = {"nrj": (nrj, nrj),
"dist": (0.7, 0.8),
"poni2": (4e-3, 4e-3),
"rot1": (0,0),
"scale": (1,1)}
gonioref = GoniometerRefinement(param,
get_position,
trans,
detector=modules[name],
wavelength=wl,
bounds=bounds
)
goniometers[name] = gonioref
for i in range(frame_start, frame_stop):
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
peak = peak_picking(name, i)
if len(peak)!=1:
continue
cp.append([peak[0]], ring=0)
img = (ds[i]).reshape((-1,1))
sg = gonioref.new_geometry("%s_%04i"%(name, i),
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
print(gonioref.chi2())
gonioref.refine2()
tths = LaB6.get_2th()
#print(tths)
for i in range(frame_stop, ds.shape[0]):
frame_name = "%s_%04i"%(name, i)
if frame_name in gonioref.single_geometries:
continue
peak = peak_picking(name, i)
ai=gonioref.get_ai(get_position(i))
tth = ai.array_from_unit(unit="2th_rad", scale=False)
tth_low = tth[20]
tth_hi = tth[-20]
ttmin, ttmax = min(tth_low, tth_hi), max(tth_low, tth_hi)
valid_peaks = numpy.logical_and(ttmin<=tths, tths<ttmax)
cnt = valid_peaks.sum()
if (len(peak) == cnt) and cnt>0:
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
#revert the order of assignment if needed !!
if tth_hi < tth_low:
peak = peak[-1::-1]
for p, r in zip(peak, numpy.where(valid_peaks)[0]):
cp.append([p], ring=r)
img = (ds[i]).reshape((-1,1))
sg = gonioref.new_geometry(frame_name,
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
#print(frame_name, len(sg.geometry_refinement.data))
print(" Number of peaks found and used for refinement")
print(sum([len(sg.geometry_refinement.data) for sg in gonioref.single_geometries.values()]))
gonioref.refine2()
gonioref.set_bounds("poni1", -1, 1)
gonioref.set_bounds("poni2", -1, 1)
gonioref.set_bounds("rot1", -1, 1)
gonioref.set_bounds("scale", 0.9, 1.1)
gonioref.refine2()
mg = gonioref.get_mg(position)
mg.radial_range = (0, 95)
images = [i.reshape(-1, 1) for i in ds]
res_mg = mg.integrate1d(images, 50000)
results[name] = res_mg
LaB6_new = get_calibrant("LaB6")
LaB6_new.wavelength = hc/gonioref.param[-1]*1e-10
p = jupyter.plot1d(res_mg, calibrant=LaB6_new)
p.figure.show()
fig, ax = plt.subplots()
ax.plot(*calc_fwhm(res_mg, LaB6_new), "o", label="FWHM")
ax.plot(*calc_peak_error(res_mg, LaB6_new, 10, 89), "o", label="error")
ax.set_title("Peak shape & error as function of the angle")
ax.set_xlabel(res_mg.unit.label)
ax.legend()
fig.show()
[35]:
add_module(ds_names[2],
92,
131,
151)
8.40633024472856e-06
Cost function before refinement: 8.40633024472856e-06
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-7.15999445e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.8085719501136514e-11
x: [ 7.200e-01 3.409e-02 4.000e-03 0.000e+00 -7.160e+01
1.000e+00 1.703e+01]
nit: 2
jac: [ 4.574e-08 8.170e-07 nan nan 9.916e-09
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 1.8085719501136514e-11
GonioParam(dist=0.7200189604924729, poni1=0.034089307266644636, poni2=0.004, rot1=0.0, offset=-71.59991818231154, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.034089307266644636
Number of peaks found and used for refinement
1093
Cost function before refinement: 5.557202432497131e-07
[ 7.20018960e-01 3.40893073e-02 4.00000000e-03 0.00000000e+00
-7.15999182e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.196098457275966e-07
x: [ 7.200e-01 3.457e-02 4.000e-03 0.000e+00 -7.160e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-9.442e-07 2.309e-08 nan nan -1.687e-11
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 1.196098457275966e-07
GonioParam(dist=0.7200183446946208, poni1=0.03456504217138477, poni2=0.004, rot1=0.0, offset=-71.5999122000852, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.034089307266644636 --> 0.03456504217138477
Cost function before refinement: 1.196098457275966e-07
[ 7.20018345e-01 3.45650422e-02 4.00000000e-03 0.00000000e+00
-7.15999122e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.570414630415977e-10
x: [ 7.207e-01 3.482e-02 4.043e-03 -3.377e-05 -7.160e+01
9.990e-01 1.703e+01]
nit: 10
jac: [ 2.182e-08 1.739e-08 -1.608e-07 1.158e-07 -1.824e-10
8.080e-09 nan]
nfev: 70
njev: 10
Cost function after refinement: 6.570414630415977e-10
GonioParam(dist=0.7206864291995188, poni1=0.03482474083000964, poni2=0.004043166292632227, rot1=-3.376839659234428e-05, offset=-71.59990881340703, scale=0.9989828730377052, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9989828730377052
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
25659
71 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[36]:
add_module(ds_names[3],
121,
159,
179)
1.0825377842185958e-06
Cost function before refinement: 1.0825377842185958e-06
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-6.57999444e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.4921937668954493e-11
x: [ 7.200e-01 3.275e-02 4.000e-03 0.000e+00 -6.580e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-8.600e-08 7.994e-08 nan nan 6.189e-10
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 2.4921937668954493e-11
GonioParam(dist=0.7200067449756247, poni1=0.03274966118487238, poni2=0.004, rot1=0.0, offset=-65.7999350248995, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.03274966118487238
Number of peaks found and used for refinement
978
Cost function before refinement: 4.7051020473788555e-07
[ 7.20006745e-01 3.27496612e-02 4.00000000e-03 0.00000000e+00
-6.57999350e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.923942799045285e-08
x: [ 7.205e-01 3.319e-02 4.000e-03 0.000e+00 -6.580e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-1.340e-09 -9.660e-12 nan nan -3.958e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 9.923942799045285e-08
GonioParam(dist=0.7204818160598766, poni1=0.0331886324835486, poni2=0.004, rot1=0.0, offset=-65.79992934734076, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7200067449756247 --> 0.7204818160598766
Cost function before refinement: 9.923942799045285e-08
[ 7.20481816e-01 3.31886325e-02 4.00000000e-03 0.00000000e+00
-6.57999293e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 5.124872429526667e-10
x: [ 7.205e-01 3.341e-02 3.978e-03 1.591e-05 -6.580e+01
9.990e-01 1.703e+01]
nit: 5
jac: [-5.053e-07 -5.371e-09 -1.298e-07 9.555e-08 -4.364e-10
2.566e-09 nan]
nfev: 36
njev: 5
Cost function after refinement: 5.124872429526667e-10
GonioParam(dist=0.7204829252945651, poni1=0.03340994866310485, poni2=0.003977911642230838, rot1=1.5909826714569855e-05, offset=-65.79992655969455, scale=0.9989991199295621, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9989991199295621
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
24025
64 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[37]:
add_module(ds_names[4],
150,
188,
208)
2.319300124385771e-07
Cost function before refinement: 2.319300124385771e-07
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-6.00000556e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.0856624761988845e-10
x: [ 7.200e-01 3.165e-02 4.000e-03 0.000e+00 -6.000e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-2.720e-07 2.346e-08 nan nan -9.144e-11
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 2.0856624761988845e-10
GonioParam(dist=0.719996884226615, poni1=0.03165314623341409, poni2=0.004, rot1=0.0, offset=-60.000059921076755, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.03165314623341409
Number of peaks found and used for refinement
864
Cost function before refinement: 3.88327645072857e-07
[ 7.19996884e-01 3.16531462e-02 4.00000000e-03 0.00000000e+00
-6.00000599e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 8.355695006101576e-08
x: [ 7.206e-01 3.205e-02 4.000e-03 0.000e+00 -6.000e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-2.494e-09 -1.973e-11 nan nan -2.941e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 8.355695006101576e-08
GonioParam(dist=0.7206442263327708, poni1=0.03204960860935797, poni2=0.004, rot1=0.0, offset=-60.00005480029715, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.719996884226615 --> 0.7206442263327708
Cost function before refinement: 8.355695006101576e-08
[ 7.20644226e-01 3.20496086e-02 4.00000000e-03 0.00000000e+00
-6.00000548e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 3.3388350811582614e-10
x: [ 7.206e-01 3.224e-02 3.976e-03 1.754e-05 -6.000e+01
9.990e-01 1.703e+01]
nit: 5
jac: [-3.349e-07 -1.351e-09 -1.135e-07 8.314e-08 -3.434e-10
3.682e-09 nan]
nfev: 36
njev: 5
Cost function after refinement: 3.3388350811582614e-10
GonioParam(dist=0.7206453076627379, poni1=0.0322371168044828, poni2=0.003975651575863446, rot1=1.7542244320799353e-05, offset=-60.000052437901985, scale=0.9990118655645687, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9990118655645687
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
22251
56 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[38]:
add_module(ds_names[5],
178,
216,
236)
1.4965819410817444e-06
Cost function before refinement: 1.4965819410817444e-06
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-5.43998333e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.0171535051300557e-10
x: [ 7.200e-01 3.288e-02 4.000e-03 0.000e+00 -5.440e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-2.790e-07 1.127e-07 nan nan 1.017e-09
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 2.0171535051300557e-10
GonioParam(dist=0.7200081377941879, poni1=0.032881397945235105, poni2=0.004, rot1=0.0, offset=-54.399822256645166, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.032881397945235105
Number of peaks found and used for refinement
754
Cost function before refinement: 3.084364540025083e-07
[ 7.20008138e-01 3.28813979e-02 4.00000000e-03 0.00000000e+00
-5.43998223e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.550962572082017e-08
x: [ 7.205e-01 3.324e-02 4.000e-03 0.000e+00 -5.440e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-1.597e-09 -1.653e-11 nan nan -3.335e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 6.550962572082017e-08
GonioParam(dist=0.7205268321085441, poni1=0.033236366248085354, poni2=0.004, rot1=0.0, offset=-54.3998176624037, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7200081377941879 --> 0.7205268321085441
Cost function before refinement: 6.550962572082017e-08
[ 7.20526832e-01 3.32363662e-02 4.00000000e-03 0.00000000e+00
-5.43998177e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.9801327203467305e-10
x: [ 7.205e-01 3.339e-02 3.973e-03 1.948e-05 -5.440e+01
9.990e-01 1.703e+01]
nit: 5
jac: [-6.288e-07 -6.165e-09 -8.418e-08 6.319e-08 -4.255e-10
1.871e-09 nan]
nfev: 36
njev: 5
Cost function after refinement: 2.9801327203467305e-10
GonioParam(dist=0.7205282661866322, poni1=0.03339481352066476, poni2=0.003972951169578437, rot1=1.948366989428078e-05, offset=-54.399815665343084, scale=0.9990259510565376, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9990259510565376
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
20549
54 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[39]:
add_module(ds_names[6],
207,
245,
266)
7.41484379835199e-08
Cost function before refinement: 7.41484379835199e-08
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-4.85998889e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 5.651340369404608e-11
x: [ 7.200e-01 3.180e-02 4.000e-03 0.000e+00 -4.860e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-2.813e-07 1.905e-08 nan nan -1.234e-10
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 5.651340369404608e-11
GonioParam(dist=0.7199985717445607, poni1=0.03180385768846998, poni2=0.004, rot1=0.0, offset=-48.599891358723, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.03180385768846998
Number of peaks found and used for refinement
626
Cost function before refinement: 2.383303034849477e-07
[ 7.19998572e-01 3.18038577e-02 4.00000000e-03 0.00000000e+00
-4.85998914e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 5.198330810552424e-08
x: [ 7.206e-01 3.211e-02 4.000e-03 0.000e+00 -4.860e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-2.199e-09 -2.359e-11 nan nan -2.599e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 5.198330810552424e-08
GonioParam(dist=0.7206068541828302, poni1=0.032113848914074164, poni2=0.004, rot1=0.0, offset=-48.59988734247567, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7199985717445607 --> 0.7206068541828302
Cost function before refinement: 5.198330810552424e-08
[ 7.20606854e-01 3.21138489e-02 4.00000000e-03 0.00000000e+00
-4.85998873e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.556084755886736e-10
x: [ 7.206e-01 3.225e-02 3.969e-03 2.238e-05 -4.860e+01
9.990e-01 1.703e+01]
nit: 6
jac: [-6.637e-07 8.690e-08 -6.746e-08 5.130e-08 7.874e-10
2.881e-09 nan]
nfev: 42
njev: 6
Cost function after refinement: 2.556084755886736e-10
GonioParam(dist=0.7206092933284303, poni1=0.0322472572293433, poni2=0.003968935506720823, rot1=2.2375489038523132e-05, offset=-48.59988566020599, scale=0.9990409081102396, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9990409081102396
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
18650
54 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[40]:
add_module(ds_names[7],
236,
273,
293)
4.377487894658635e-06
Cost function before refinement: 4.377487894658635e-06
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-4.27998889e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.109058146293743e-11
x: [ 7.200e-01 3.049e-02 4.000e-03 0.000e+00 -4.280e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-3.538e-07 3.502e-09 nan nan -3.652e-10
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 9.109058146293743e-11
GonioParam(dist=0.7199860337927957, poni1=0.030492505461750096, poni2=0.004, rot1=0.0, offset=-42.79990784449265, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.030492505461750096
Number of peaks found and used for refinement
543
Cost function before refinement: 1.9115344534723355e-07
[ 7.19986034e-01 3.04925055e-02 4.00000000e-03 0.00000000e+00
-4.27999078e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 4.1213003335648066e-08
x: [ 7.207e-01 3.077e-02 4.000e-03 0.000e+00 -4.280e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-3.060e-09 -2.445e-11 nan nan -3.621e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 4.1213003335648066e-08
GonioParam(dist=0.7207055546218871, poni1=0.030768430906894705, poni2=0.004, rot1=0.0, offset=-42.79990417397939, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7199860337927957 --> 0.7207055546218871
Cost function before refinement: 4.1213003335648066e-08
[ 7.20705555e-01 3.07684309e-02 4.00000000e-03 0.00000000e+00
-4.27999042e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.9055738692529888e-10
x: [ 7.207e-01 3.086e-02 3.965e-03 2.507e-05 -4.280e+01
9.991e-01 1.703e+01]
nit: 5
jac: [-5.207e-07 -3.451e-07 -5.785e-08 4.380e-08 -4.683e-09
9.352e-08 nan]
nfev: 35
njev: 5
Cost function after refinement: 1.9055738692529888e-10
GonioParam(dist=0.7207071581327096, poni1=0.030863392389826882, poni2=0.003965206218806199, rot1=2.50696591730669e-05, offset=-42.799902973938835, scale=0.9990506059525829, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9990506059525829
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
16803
45 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[41]:
add_module(ds_names[8],
264,
302,
322)
9.734136493812359e-08
Cost function before refinement: 9.734136493812359e-08
[ 7.20000000e-01 3.20000000e-02 4.00000000e-03 0.00000000e+00
-3.72000000e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.3838346593948452e-10
x: [ 7.200e-01 3.178e-02 4.000e-03 0.000e+00 -3.720e+01
1.000e+00 1.703e+01]
nit: 2
jac: [-5.987e-07 2.236e-08 nan nan -2.112e-10
nan nan]
nfev: 9
njev: 2
Cost function after refinement: 1.3838346593948452e-10
GonioParam(dist=0.7199980533133429, poni1=0.031775349005215406, poni2=0.004, rot1=0.0, offset=-37.20000282728476, scale=1.0, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032 --> 0.031775349005215406
Number of peaks found and used for refinement
454
Cost function before refinement: 1.4308139308482607e-07
[ 7.19998053e-01 3.17753490e-02 4.00000000e-03 0.00000000e+00
-3.72000028e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 3.226788268488173e-08
x: [ 7.208e-01 3.201e-02 4.000e-03 0.000e+00 -3.720e+01
1.000e+00 1.703e+01]
nit: 7
jac: [-3.925e-09 -5.278e-11 nan nan -3.042e-10
nan nan]
nfev: 29
njev: 7
Cost function after refinement: 3.226788268488173e-08
GonioParam(dist=0.7208093283946152, poni1=0.03201295264882314, poni2=0.004, rot1=0.0, offset=-37.19999969774526, scale=1.0, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7199980533133429 --> 0.7208093283946152
Cost function before refinement: 3.226788268488173e-08
[ 7.20809328e-01 3.20129526e-02 4.00000000e-03 0.00000000e+00
-3.71999997e+01 1.00000000e+00 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.851171758694135e-10
x: [ 7.208e-01 3.208e-02 3.960e-03 2.909e-05 -3.720e+01
9.991e-01 1.703e+01]
nit: 5
jac: [-5.574e-07 1.084e-09 -4.859e-08 3.728e-08 -3.498e-10
3.440e-09 nan]
nfev: 35
njev: 5
Cost function after refinement: 1.851171758694135e-10
GonioParam(dist=0.7208107860401025, poni1=0.03208175644866908, poni2=0.003959633637265279, rot1=2.90905939606598e-05, offset=-37.199998826575715, scale=0.9990588539456343, nrj=17.027082549190933)
maxdelta on: scale (5) 1.0 --> 0.9990588539456343
<ipython-input-34-ced5c09b5da0>:99: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
p.figure.show()
16875
45 54
<ipython-input-34-ced5c09b5da0>:107: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[42]:
len(goniometers)
[42]:
9
[43]:
# print all the parameters to be able to compare them visually
goniometers["data_12"] = goniometers["data_11"]
for name in ds_names:
print(name, *["%8.4e"%i for i in goniometers[name].param])
data_02 7.2310e-01 3.2054e-02 3.9849e-03 1.0948e-05 -8.2800e+01 9.9942e-01 1.7027e+01
data_03 7.2068e-01 3.3678e-02 4.0580e-03 -4.4486e-05 -7.7200e+01 9.9897e-01 1.7027e+01
data_04 7.2069e-01 3.4825e-02 4.0432e-03 -3.3768e-05 -7.1600e+01 9.9898e-01 1.7027e+01
data_05 7.2048e-01 3.3410e-02 3.9779e-03 1.5910e-05 -6.5800e+01 9.9900e-01 1.7027e+01
data_07 7.2065e-01 3.2237e-02 3.9757e-03 1.7542e-05 -6.0000e+01 9.9901e-01 1.7027e+01
data_08 7.2053e-01 3.3395e-02 3.9730e-03 1.9484e-05 -5.4400e+01 9.9903e-01 1.7027e+01
data_09 7.2061e-01 3.2247e-02 3.9689e-03 2.2375e-05 -4.8600e+01 9.9904e-01 1.7027e+01
data_10 7.2071e-01 3.0863e-02 3.9652e-03 2.5070e-05 -4.2800e+01 9.9905e-01 1.7027e+01
data_11 7.2081e-01 3.2082e-02 3.9596e-03 2.9091e-05 -3.7200e+01 9.9906e-01 1.7027e+01
data_12 7.2081e-01 3.2082e-02 3.9596e-03 2.9091e-05 -3.7200e+01 9.9906e-01 1.7027e+01
Use the negative part of the spectum …
Until now, we used only the data where 2th >0 For the last modules, this thows away half of the data.
We setup here a way to assign the peaks for the negative part of the spectrum.
[44]:
def complete_gonio(module_id=None, name=None):
"Scan missing frames for un-indexed peaks"
if name is None:
name = ds_names[module_id]
gonioref = goniometers[name]
ds = data[name]
print("Number of peaks previously found:",
sum([len(sg.geometry_refinement.data) for sg in gonioref.single_geometries.values()]))
tths = LaB6.get_2th()
for i in range(ds.shape[0]):
frame_name = "%s_%04i"%(name, i)
if frame_name in gonioref.single_geometries:
continue
peak = peak_picking(name, i)
ai=gonioref.get_ai(get_position(i))
tth = ai.array_from_unit(unit="2th_rad", scale=False)
tth_low = tth[20]
tth_hi = tth[-20]
ttmin, ttmax = min(tth_low, tth_hi), max(tth_low, tth_hi)
valid_peaks = numpy.logical_and(ttmin<=tths, tths<ttmax)
cnt = valid_peaks.sum()
if (len(peak) == cnt) and cnt>0:
cp = ControlPoints(calibrant=LaB6, wavelength=wl)
#revert the order of assignment if needed !!
if tth_hi < tth_low:
peak = peak[-1::-1]
for p, r in zip(peak, numpy.where(valid_peaks)[0]):
cp.append([p], ring=r)
img = ds[i].reshape((-1,1))
sg = gonioref.new_geometry(frame_name,
image=img,
metadata=i,
control_points=cp,
calibrant=LaB6)
sg.geometry_refinement.data = numpy.array(cp.getList())
#print(frame_name, len(sg.geometry_refinement.data))
print("Number of peaks found after re-scan:",
sum([len(sg.geometry_refinement.data) for sg in gonioref.single_geometries.values()]))
return gonioref
[45]:
gonio8 = complete_gonio(module_id=8)
gonio8.refine2()
Number of peaks previously found: 454
Number of peaks found after re-scan: 1006
Cost function before refinement: 4.70708266491699e-09
[ 7.20810786e-01 3.20817564e-02 3.95963364e-03 2.90905940e-05
-3.71999988e+01 9.99058854e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.80271531341448e-10
x: [ 7.208e-01 3.208e-02 3.962e-03 2.771e-05 -3.720e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-8.293e-07 1.798e-08 -2.226e-07 1.638e-07 -1.374e-10
-3.004e-08 nan]
nfev: 29
njev: 4
Cost function after refinement: 9.80271531341448e-10
GonioParam(dist=0.7208134356283892, poni1=0.032084918686756565, poni2=0.003961539926573026, rot1=2.770579470925171e-05, offset=-37.19999878568301, scale=0.9989913574687115, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9990588539456343 --> 0.9989913574687115
[45]:
array([ 7.20813436e-01, 3.20849187e-02, 3.96153993e-03, 2.77057947e-05,
-3.71999988e+01, 9.98991357e-01, 1.70270825e+01])
[46]:
gonio7 = complete_gonio(module_id=7)
Number of peaks previously found: 543
Number of peaks found after re-scan: 1004
[47]:
gonio7.refine2()
Cost function before refinement: 2.049771586732014e-09
[ 7.20707158e-01 3.08633924e-02 3.96520622e-03 2.50696592e-05
-4.27999030e+01 9.99050606e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.98216618419778e-10
x: [ 7.207e-01 3.087e-02 3.967e-03 2.399e-05 -4.280e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-6.616e-07 1.821e-08 -2.311e-07 1.692e-07 -1.059e-10
-2.989e-08 nan]
nfev: 29
njev: 4
Cost function after refinement: 6.98216618419778e-10
GonioParam(dist=0.720709214442203, poni1=0.03086624488929287, poni2=0.0039666907721100705, rot1=2.3991415326512732e-05, offset=-42.79990293705691, scale=0.9990066875069397, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9990506059525829 --> 0.9990066875069397
[47]:
array([ 7.20709214e-01, 3.08662449e-02, 3.96669077e-03, 2.39914153e-05,
-4.27999029e+01, 9.99006688e-01, 1.70270825e+01])
[48]:
gonio6 = complete_gonio(module_id=6)
gonio6.refine2()
Number of peaks previously found: 626
Number of peaks found after re-scan: 990
Cost function before refinement: 9.054150799194566e-10
[ 7.20609293e-01 3.22472572e-02 3.96893551e-03 2.23754890e-05
-4.85998857e+01 9.99040908e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.093707809212381e-10
x: [ 7.206e-01 3.225e-02 3.970e-03 2.155e-05 -4.860e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-8.022e-07 3.637e-08 -2.332e-07 1.713e-07 1.531e-10
-6.906e-08 nan]
nfev: 29
njev: 4
Cost function after refinement: 6.093707809212381e-10
GonioParam(dist=0.7206117521177711, poni1=0.03224903534860273, poni2=0.003970067520600401, rot1=2.1549818979054536e-05, offset=-48.59988563691411, scale=0.9990180421528526, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9990409081102396 --> 0.9990180421528526
[48]:
array([ 7.20611752e-01, 3.22490353e-02, 3.97006752e-03, 2.15498190e-05,
-4.85998856e+01, 9.99018042e-01, 1.70270825e+01])
[49]:
gonio5 = complete_gonio(module_id=5)
gonio5.refine2()
Number of peaks previously found: 754
Number of peaks found after re-scan: 1038
Cost function before refinement: 5.322293237710593e-10
[ 7.20528266e-01 3.33948135e-02 3.97295117e-03 1.94836699e-05
-5.43998157e+01 9.99025951e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 5.247776760667795e-10
x: [ 7.205e-01 3.340e-02 3.974e-03 1.898e-05 -5.440e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-6.548e-07 1.292e-07 -1.627e-07 1.198e-07 1.289e-09
-2.774e-07 nan]
nfev: 29
njev: 4
Cost function after refinement: 5.247776760667795e-10
GonioParam(dist=0.7205306819340235, poni1=0.033395338635169546, poni2=0.003973636918904946, rot1=1.8979819424984276e-05, offset=-54.39981565750616, scale=0.9990218992359475, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9990259510565376 --> 0.9990218992359475
[49]:
array([ 7.20530682e-01, 3.33953386e-02, 3.97363692e-03, 1.89798194e-05,
-5.43998157e+01, 9.99021899e-01, 1.70270825e+01])
[50]:
gonio4 = complete_gonio(module_id=4)
gonio4.refine2()
Number of peaks previously found: 864
Number of peaks found after re-scan: 1081
Cost function before refinement: 4.522892051482172e-10
[ 7.20645308e-01 3.22371168e-02 3.97565158e-03 1.75422443e-05
-6.00000524e+01 9.99011866e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 4.3299107269375743e-10
x: [ 7.206e-01 3.224e-02 3.976e-03 1.742e-05 -6.000e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-3.236e-07 -1.338e-08 -8.079e-08 5.953e-08 -4.927e-10
3.166e-08 nan]
nfev: 29
njev: 4
Cost function after refinement: 4.3299107269375743e-10
GonioParam(dist=0.7206464170800312, poni1=0.032235988252634566, poni2=0.003975819586239342, rot1=1.7416686730645905e-05, offset=-60.00005245097538, scale=0.9990197113106509, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9990118655645687 --> 0.9990197113106509
[50]:
array([ 7.20646417e-01, 3.22359883e-02, 3.97581959e-03, 1.74166867e-05,
-6.00000525e+01, 9.99019711e-01, 1.70270825e+01])
[51]:
gonio3 = complete_gonio(module_id=3)
gonio3.refine2()
Number of peaks previously found: 978
Number of peaks found after re-scan: 1156
Cost function before refinement: 7.053066687141914e-10
[ 7.20482925e-01 3.34099487e-02 3.97791164e-03 1.59098267e-05
-6.57999266e+01 9.98999120e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.316703679729696e-10
x: [ 7.205e-01 3.341e-02 3.978e-03 1.605e-05 -6.580e+01
9.990e-01 1.703e+01]
nit: 4
jac: [-4.998e-07 -1.632e-08 9.666e-09 -4.950e-09 -5.754e-10
3.901e-08 nan]
nfev: 29
njev: 4
Cost function after refinement: 6.316703679729696e-10
GonioParam(dist=0.7204847938379576, poni1=0.03340748825027609, poni2=0.003977700789074132, rot1=1.6054211747525078e-05, offset=-65.79992658926287, scale=0.9990156265501668, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9989991199295621 --> 0.9990156265501668
[51]:
array([ 7.20484794e-01, 3.34074883e-02, 3.97770079e-03, 1.60542117e-05,
-6.57999266e+01, 9.99015627e-01, 1.70270825e+01])
[52]:
gonio2 = complete_gonio(module_id=2)
gonio2.refine2()
Number of peaks previously found: 1093
Number of peaks found after re-scan: 1229
Cost function before refinement: 8.306329045500625e-10
[ 7.20686429e-01 3.48247408e-02 4.04316629e-03 -3.37683966e-05
-7.15999088e+01 9.98982873e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 7.77171881108816e-10
x: [ 7.207e-01 3.482e-02 4.043e-03 -3.357e-05 -7.160e+01
9.990e-01 1.703e+01]
nit: 4
jac: [ 9.152e-08 -1.052e-09 5.616e-08 -4.083e-08 -4.257e-10
2.851e-09 nan]
nfev: 29
njev: 4
Cost function after refinement: 7.77171881108816e-10
GonioParam(dist=0.7206860334778521, poni1=0.03482262187955647, poni2=0.004042896907016237, rot1=-3.3572702282314364e-05, offset=-71.5999088383209, scale=0.9989982448524656, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9989828730377052 --> 0.9989982448524656
[52]:
array([ 7.20686033e-01, 3.48226219e-02, 4.04289691e-03, -3.35727023e-05,
-7.15999088e+01, 9.98998245e-01, 1.70270825e+01])
[53]:
gonio1 = complete_gonio(module_id=1)
gonio1.refine2()
Number of peaks previously found: 1183
Number of peaks found after re-scan: 1285
Cost function before refinement: 9.822430104344508e-10
[ 7.20682073e-01 3.36782330e-02 4.05801068e-03 -4.44855659e-05
-7.71999788e+01 9.98967233e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.683583557485058e-10
x: [ 7.207e-01 3.368e-02 4.058e-03 -4.454e-05 -7.720e+01
9.990e-01 1.703e+01]
nit: 4
jac: [ 1.052e-07 -1.409e-09 -2.647e-09 1.496e-09 -3.432e-10
4.363e-09 nan]
nfev: 29
njev: 4
Cost function after refinement: 9.683583557485058e-10
GonioParam(dist=0.7206815883828414, poni1=0.03367666939062012, poni2=0.004058082310428429, rot1=-4.453529297359146e-05, offset=-77.19997881135518, scale=0.998975905641853, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9989672333400683 --> 0.998975905641853
[53]:
array([ 7.20681588e-01, 3.36766694e-02, 4.05808231e-03, -4.45352930e-05,
-7.71999788e+01, 9.98975906e-01, 1.70270825e+01])
[54]:
gonio0 = complete_gonio(module_id=0)
gonio0.refine2()
Number of peaks previously found: 1203
Number of peaks found after re-scan: 1255
Cost function before refinement: 2.580092915553823e-06
[ 7.23096244e-01 3.20537008e-02 3.98485488e-03 1.09477702e-05
-8.27999440e+01 9.99415058e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 2.579987978563334e-06
x: [ 7.231e-01 3.207e-02 3.981e-03 1.351e-05 -8.280e+01
9.994e-01 1.703e+01]
nit: 4
jac: [-1.439e-07 3.613e-08 5.528e-07 -3.992e-07 1.021e-09
-1.146e-07 nan]
nfev: 29
njev: 4
Cost function after refinement: 2.579987978563334e-06
GonioParam(dist=0.7230969947734238, poni1=0.03206584436046777, poni2=0.0039813020995294445, rot1=1.3513751444584455e-05, offset=-82.79994383845344, scale=0.9993962397688153, nrj=17.027082549190933)
maxdelta on: scale (5) 0.9994150576472143 --> 0.9993962397688153
[54]:
array([ 7.23096995e-01, 3.20658444e-02, 3.98130210e-03, 1.35137514e-05,
-8.27999438e+01, 9.99396240e-01, 1.70270825e+01])
[55]:
#Rescan module0 which looks much different:
gonio0.single_geometries.clear()
gonio0 = complete_gonio(module_id=0)
gonio0.refine2()
Number of peaks previously found: 0
Number of peaks found after re-scan: 1250
Cost function before refinement: 9.454139688783061e-07
[ 7.23096995e-01 3.20658444e-02 3.98130210e-03 1.35137514e-05
-8.27999438e+01 9.99396240e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.315197978919561e-07
x: [ 7.220e-01 3.220e-02 3.944e-03 4.514e-05 -8.280e+01
9.991e-01 1.703e+01]
nit: 9
jac: [-1.206e-08 -1.145e-08 1.037e-07 -7.484e-08 1.332e-10
-9.529e-09 nan]
nfev: 64
njev: 9
Cost function after refinement: 9.315197978919561e-07
GonioParam(dist=0.7219737516574156, poni1=0.03220031516651301, poni2=0.003943811306578325, rot1=4.514209993152478e-05, offset=-82.79994226182613, scale=0.9991407871099522, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7230969947734238 --> 0.7219737516574156
[55]:
array([ 7.21973752e-01, 3.22003152e-02, 3.94381131e-03, 4.51420999e-05,
-8.27999423e+01, 9.99140787e-01, 1.70270825e+01])
Discard wronly assigned peaks
We have seen previously that some modules have a much higher residual error, while all have almost the same number of peaks recorded and fitted.
Some frames are contributing much more than all the other in those badly-fitted data. Let’s spot them and re-assign them
[56]:
#search for mis-assigned peaks in module #0
labels = []
errors = []
for lbl,sg in gonio0.single_geometries.items():
labels.append(lbl)
errors.append(sg.geometry_refinement.chi2())
s = numpy.argsort(errors)
for i in s[-10:]:
print(labels[i], errors[i])
data_02_0455 6.313037376229225e-07
data_02_0464 6.523610354376743e-07
data_02_0465 6.891635784533209e-07
data_02_0457 7.370412438549099e-07
data_02_0456 7.376223648967419e-07
data_02_0460 7.458554081403769e-07
data_02_0466 7.671413725143998e-07
data_02_0461 7.935694034310535e-07
data_02_0462 7.971541383592315e-07
data_02_0480 0.0011307660991106366
[57]:
#remove wrongly assigned peak for frame 480
print(gonio0.single_geometries.pop("data_02_0480").control_points)
gonio0.refine2()
gonio0 = complete_gonio(module_id=0)
gonio0.refine2()
ControlPoints instance containing 4 group of point:
LaB6 Calibrant with 109 reflections at wavelength 7.281587910025816e-11
Containing 4 groups of points:
#psg ring 52: 1 points
#psh ring 53: 1 points
#psi ring 54: 1 points
#psj ring 55: 1 points
Cost function before refinement: 8.041212761797189e-09
[ 7.21973752e-01 3.22003152e-02 3.94381131e-03 4.51420999e-05
-8.27999423e+01 9.99140787e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.214244938591865e-10
x: [ 7.206e-01 3.228e-02 4.043e-03 -2.053e-05 -8.280e+01
9.990e-01 1.703e+01]
nit: 9
jac: [-3.087e-08 -3.829e-08 -1.829e-07 1.319e-07 -8.443e-10
-1.615e-08 nan]
nfev: 64
njev: 9
Cost function after refinement: 9.214244938591865e-10
GonioParam(dist=0.7205957893874516, poni1=0.03228440830812512, poni2=0.00404252549764979, rot1=-2.053391818186055e-05, offset=-82.79994103787172, scale=0.998975174801073, nrj=17.027082549190933)
maxdelta on: dist (0) 0.7219737516574156 --> 0.7205957893874516
Number of peaks previously found: 1246
Number of peaks found after re-scan: 1263
Cost function before refinement: 9.385076296783189e-10
[ 7.20595789e-01 3.22844083e-02 4.04252550e-03 -2.05339182e-05
-8.27999410e+01 9.98975175e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.382098917207145e-10
x: [ 7.206e-01 3.228e-02 4.043e-03 -2.057e-05 -8.280e+01
9.990e-01 1.703e+01]
nit: 1
jac: [-2.666e-08 -1.466e-06 -1.817e-07 1.310e-07 -1.881e-08
-4.389e-08 nan]
nfev: 9
njev: 1
Cost function after refinement: 9.382098917207145e-10
GonioParam(dist=0.7205957963976899, poni1=0.03228479377438343, poni2=0.004042573274141478, rot1=-2.0568374189891214e-05, offset=-82.79994103292513, scale=0.9989751863414432, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.03228440830812512 --> 0.03228479377438343
[57]:
array([ 7.20595796e-01, 3.22847938e-02, 4.04257327e-03, -2.05683742e-05,
-8.27999410e+01, 9.98975186e-01, 1.70270825e+01])
[58]:
def search_outliers(module_id=None, name=None, threshold=1.2):
"Search for wrongly assigned peaks"
if name is None:
name = ds_names[module_id]
gonioref = goniometers[name]
labels = []
errors = []
for lbl,sg in gonioref.single_geometries.items():
labels.append(lbl)
errors.append(sg.geometry_refinement.chi2())
s = numpy.argsort(errors)
last = errors[s[-1]]
to_remove = []
for i in s[-1::-1]:
lbl = labels[i]
current = errors[i]
print(lbl , current, last, last/current)
if threshold*current<last:
break
last=current
to_remove.append(lbl)
return to_remove
for lbl in search_outliers(8):
gonio8.single_geometries.pop(lbl)
gonio8.refine2()
gonio8 = complete_gonio(module_id=8)
gonio8.refine2()
data_11_0495 1.465326424860063e-06 1.465326424860063e-06 1.0
data_11_0496 1.4360706723304346e-06 1.465326424860063e-06 1.0203720841134876
data_11_0499 1.4166274578619404e-06 1.4360706723304346e-06 1.0137250018419373
data_11_0497 1.4084298974827966e-06 1.4166274578619404e-06 1.0058203538520412
data_11_0498 1.3832745249726834e-06 1.4084298974827966e-06 1.0181853797319154
data_11_0500 1.3708837484278925e-06 1.3832745249726834e-06 1.0090385319389776
data_11_0492 1.329563509601561e-06 1.3708837484278925e-06 1.0310780481924584
data_11_0489 1.328744438946604e-06 1.329563509601561e-06 1.0006164245214877
data_11_0491 1.3215889751287176e-06 1.328744438946604e-06 1.0054142883699446
data_11_0490 1.301797579259692e-06 1.3215889751287176e-06 1.015203128492742
data_11_0494 1.2865330381985967e-06 1.301797579259692e-06 1.0118648651903015
data_11_0493 1.2861198600451933e-06 1.2865330381985967e-06 1.0003212594457478
data_11_0483 1.224424745072622e-06 1.2861198600451933e-06 1.0503870206975539
data_11_0486 1.207349392777562e-06 1.224424745072622e-06 1.0141428424921617
data_11_0485 1.2017864447065466e-06 1.207349392777562e-06 1.0046288989990844
data_11_0484 1.199457920840169e-06 1.2017864447065466e-06 1.0019413135099784
data_11_0487 1.1764421082118776e-06 1.199457920840169e-06 1.0195639143376753
data_11_0477 1.1389376207700437e-06 1.1764421082118776e-06 1.032929360447745
data_11_0478 1.1156529806074908e-06 1.1389376207700437e-06 1.020870862685164
data_11_0481 1.1104608438150525e-06 1.1156529806074908e-06 1.0046756594988082
data_11_0479 1.1057542147152647e-06 1.1104608438150525e-06 1.0042564875965676
data_11_0480 1.1002197562639265e-06 1.1057542147152647e-06 1.0050303209152796
data_11_0482 8.673122808629696e-07 1.1002197562639265e-06 1.2685393491364099
Cost function before refinement: 1.0594739815701638e-09
[ 7.20813436e-01 3.20849187e-02 3.96153993e-03 2.77057947e-05
-3.71999988e+01 9.98991357e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 1.0589425492798256e-09
x: [ 7.208e-01 3.209e-02 3.962e-03 2.762e-05 -3.720e+01
9.990e-01 1.703e+01]
nit: 1
jac: [-8.451e-07 -1.253e-06 -2.518e-07 1.849e-07 -1.614e-08
9.174e-08 nan]
nfev: 9
njev: 1
Cost function after refinement: 1.0589425492798256e-09
GonioParam(dist=0.7208138053467276, poni1=0.032085466901404895, poni2=0.003961650094174231, rot1=2.7624888491723393e-05, offset=-37.19999877862184, scale=0.9989913173307039, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032084918686756565 --> 0.032085466901404895
Number of peaks previously found: 918
Number of peaks found after re-scan: 1006
Cost function before refinement: 9.804577380710388e-10
[ 7.20813805e-01 3.20854669e-02 3.96165009e-03 2.76248885e-05
-3.71999988e+01 9.98991317e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 9.797992160755848e-10
x: [ 7.208e-01 3.209e-02 3.962e-03 2.759e-05 -3.720e+01
9.990e-01 1.703e+01]
nit: 1
jac: [-8.288e-07 2.044e-06 -2.226e-07 1.638e-07 2.536e-08
1.004e-06 nan]
nfev: 9
njev: 1
Cost function after refinement: 9.797992160755848e-10
GonioParam(dist=0.7208139923453484, poni1=0.03208500579330537, poni2=0.003961700327630729, rot1=2.7587922995196375e-05, offset=-37.19999878434437, scale=0.9989910907851859, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.032085466901404895 --> 0.03208500579330537
[58]:
array([ 7.20813992e-01, 3.20850058e-02, 3.96170033e-03, 2.75879230e-05,
-3.71999988e+01, 9.98991091e-01, 1.70270825e+01])
[59]:
print(gonio7.chi2())
for lbl in search_outliers(7):
gonio7.single_geometries.pop(lbl)
gonio7.refine2()
gonio7 = complete_gonio(module_id=7)
gonio7.refine2()
6.98216618419778e-10
data_10_0292 4.478313610295236e-06 4.478313610295236e-06 1.0
data_10_0291 4.441618075732054e-06 4.478313610295236e-06 1.008261749195339
data_10_0288 4.421114039991817e-06 4.441618075732054e-06 1.0046377531895274
data_10_0285 4.4151969079424084e-06 4.421114039991817e-06 1.0013401739883365
data_10_0290 4.4114083219467785e-06 4.4151969079424084e-06 1.0008588155344365
data_10_0289 4.40450577270408e-06 4.4114083219467785e-06 1.0015671563618955
data_10_0284 4.386350135873239e-06 4.40450577270408e-06 1.0041391216543243
data_10_0275 4.376873169865872e-06 4.386350135873239e-06 1.0021652366060354
data_10_0286 4.374944729160871e-06 4.376873169865872e-06 1.0004407920155303
data_10_0283 4.360029737297362e-06 4.374944729160871e-06 1.0034208463616476
data_10_0281 4.3581187741366675e-06 4.360029737297362e-06 1.0004384834970619
data_10_0282 4.3552714424913155e-06 4.3581187741366675e-06 1.0006537667474804
data_10_0278 4.352333646868408e-06 4.3552714424913155e-06 1.0006749932016406
data_10_0280 4.3515501506852325e-06 4.352333646868408e-06 1.0001800499031481
data_10_0277 4.349866613398178e-06 4.3515501506852325e-06 1.0003870319337769
data_10_0276 4.3396700589153925e-06 4.349866613398178e-06 1.0023496151422475
data_10_0279 4.338090818701843e-06 4.3396700589153925e-06 1.0003640403761813
data_10_0274 4.337324373540883e-06 4.338090818701843e-06 1.0001767092094
data_10_0287 4.319689618966426e-06 4.337324373540883e-06 1.004082412425427
data_10_0496 1.950808600710568e-06 4.319689618966426e-06 2.2143072454124972
Cost function before refinement: 6.590244869729802e-10
[ 7.20709214e-01 3.08662449e-02 3.96669077e-03 2.39914153e-05
-4.27999029e+01 9.99006688e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.580238959302392e-10
x: [ 7.207e-01 3.087e-02 3.967e-03 2.397e-05 -4.280e+01
9.990e-01 1.703e+01]
nit: 2
jac: [-6.411e-07 8.001e-09 -1.445e-07 1.068e-07 -2.377e-10
-3.751e-07 nan]
nfev: 15
njev: 2
Cost function after refinement: 6.580238959302392e-10
GonioParam(dist=0.7207093528388456, poni1=0.03086564460420738, poni2=0.003966722984588224, rot1=2.396764015678069e-05, offset=-42.799902944538964, scale=0.9990064423877438, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.03086624488929287 --> 0.03086564460420738
Number of peaks previously found: 985
Number of peaks found after re-scan: 1004
Cost function before refinement: 6.991408418442968e-10
[ 7.20709353e-01 3.08656446e-02 3.96672298e-03 2.39676402e-05
-4.27999029e+01 9.99006442e-01 1.70270825e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.980019649241363e-10
x: [ 7.207e-01 3.087e-02 3.967e-03 2.393e-05 -4.280e+01
9.990e-01 1.703e+01]
nit: 2
jac: [-6.609e-07 1.522e-07 -2.321e-07 1.700e-07 1.580e-09
9.613e-08 nan]
nfev: 15
njev: 2
Cost function after refinement: 6.980019649241363e-10
GonioParam(dist=0.7207094971365863, poni1=0.030866238233237652, poni2=0.003966772614681739, rot1=2.393128798437873e-05, offset=-42.7999029369946, scale=0.9990067735020023, nrj=17.027082549190933)
maxdelta on: poni1 (1) 0.03086564460420738 --> 0.030866238233237652
[59]:
array([ 7.20709497e-01, 3.08662382e-02, 3.96677261e-03, 2.39312880e-05,
-4.27999029e+01, 9.99006774e-01, 1.70270825e+01])
[60]:
print(gonio0.chi2())
print(len(search_outliers(0)))
# for lbl in search_outliers(7):
# gonio7.single_geometries.pop(lbl)
# gonio7.refine2()
# gonio7 = complete_gonio(module_id=7)
# gonio7.refine2()
9.382098917207145e-10
data_02_0462 7.971541383592315e-07 7.971541383592315e-07 1.0
data_02_0461 7.935694034310535e-07 7.971541383592315e-07 1.0045172292589397
data_02_0466 7.671413725143998e-07 7.935694034310535e-07 1.0344500138612425
data_02_0460 7.458554081403769e-07 7.671413725143998e-07 1.0285389958183646
data_02_0456 7.376223648967419e-07 7.458554081403769e-07 1.0111615965505432
data_02_0457 7.370412438549099e-07 7.376223648967419e-07 1.000788451184621
data_02_0465 6.891635784533209e-07 7.370412438549099e-07 1.0694721353514358
data_02_0464 6.523610354376743e-07 6.891635784533209e-07 1.0564143794868979
data_02_0455 6.313037376229225e-07 6.523610354376743e-07 1.0333552560516113
data_02_0463 6.240793018859902e-07 6.313037376229225e-07 1.0115761502025462
data_02_0459 6.163884625324197e-07 6.240793018859902e-07 1.012477260398374
data_02_0448 6.082013535587733e-07 6.163884625324197e-07 1.0134611817710386
data_02_0454 5.868915926138719e-07 6.082013535587733e-07 1.0363095352073335
data_02_0439 5.868104673105387e-07 5.868915926138719e-07 1.0001382478804528
data_02_0435 5.859904555714726e-07 5.868104673105387e-07 1.0013993602306481
data_02_0426 5.755880410766981e-07 5.859904555714726e-07 1.0180726730793705
data_02_0458 5.724204045086455e-07 5.755880410766981e-07 1.005533759004925
data_02_0453 5.590326810656076e-07 5.724204045086455e-07 1.0239480157358936
data_02_0434 5.570140355534997e-07 5.590326810656076e-07 1.0036240478394802
data_02_0444 5.524794596805079e-07 5.570140355534997e-07 1.0082076822830919
data_02_0430 5.489352907178469e-07 5.524794596805079e-07 1.0064564421755908
data_02_0443 5.459467871390915e-07 5.489352907178469e-07 1.0054739832693511
data_02_0447 5.424205085181835e-07 5.459467871390915e-07 1.0065010053372454
data_02_0468 5.360712035478987e-07 5.424205085181835e-07 1.011844144823044
data_02_0425 5.340449951225031e-07 5.360712035478987e-07 1.0037940781093375
data_02_0467 5.066723540804425e-07 5.340449951225031e-07 1.054024342993292
data_02_0438 5.014370092167499e-07 5.066723540804425e-07 1.0104406830119506
data_02_0429 4.994038918066035e-07 5.014370092167499e-07 1.0040710884386415
data_02_0442 4.988910402085017e-07 4.994038918066035e-07 1.00102798318023
data_02_0421 4.887205487635657e-07 4.988910402085017e-07 1.020810443658788
data_02_0452 4.818395002871005e-07 4.887205487635657e-07 1.014280789500166
data_02_0441 4.7080825196106415e-07 4.818395002871005e-07 1.0234304481284848
data_02_0433 4.690750465319513e-07 4.7080825196106415e-07 1.0036949427216968
data_02_0437 4.674410331405129e-07 4.690750465319513e-07 1.0034956567258553
data_02_0420 4.5997632178011953e-07 4.674410331405129e-07 1.0162284687427057
data_02_0451 4.5329039220439244e-07 4.5997632178011953e-07 1.0147497712078406
data_02_0446 4.453788918120156e-07 4.5329039220439244e-07 1.0177635279485049
data_02_0428 4.424565896011663e-07 4.453788918120156e-07 1.006604720733131
data_02_0424 4.3749908369213206e-07 4.424565896011663e-07 1.0113314658106183
data_02_0432 4.3734558209974527e-07 4.3749908369213206e-07 1.000350984664461
data_02_0436 4.3318674322151904e-07 4.3734558209974527e-07 1.0096005682152178
data_02_0450 4.262100641221257e-07 4.3318674322151904e-07 1.0163691092413862
data_02_0445 4.1944403134349654e-07 4.262100641221257e-07 1.0161309549618749
data_02_0427 4.184056122213615e-07 4.1944403134349654e-07 1.0024818479767084
data_02_0423 4.1311974462174437e-07 4.184056122213615e-07 1.012795001130863
data_02_0419 3.8768794142524627e-07 4.1311974462174437e-07 1.0655986438551683
data_02_0408 3.542344482717018e-07 3.8768794142524627e-07 1.0944388478217264
data_02_0418 3.436568182824559e-07 3.542344482717018e-07 1.0307796308017727
data_02_0412 3.436508350200033e-07 3.436568182824559e-07 1.0000174108770963
data_02_0398 3.357912171790021e-07 3.436508350200033e-07 1.0234062638893007
data_02_0389 3.2715501976740725e-07 3.357912171790021e-07 1.0263978752877911
data_02_0394 3.2138367696552696e-07 3.2715501976740725e-07 1.017957796912316
data_02_0416 3.194680972472967e-07 3.2138367696552696e-07 1.005996153402283
data_02_0403 3.071335990321196e-07 3.194680972472967e-07 1.0401600419297896
data_02_0411 3.0688039545695554e-07 3.071335990321196e-07 1.0008250887932644
data_02_0407 3.062881042731559e-07 3.0688039545695554e-07 1.0019337714248
data_02_0402 3.049021910877446e-07 3.062881042731559e-07 1.0045454353098187
data_02_0393 2.9809089008092075e-07 3.049021910877446e-07 1.022849745609384
data_02_0388 2.945098238482041e-07 2.9809089008092075e-07 1.0121594118183386
data_02_0397 2.9208260134706256e-07 2.945098238482041e-07 1.0083100550664346
data_02_0384 2.824644866830826e-07 2.9208260134706256e-07 1.0340507041324851
data_02_0406 2.8172360974302977e-07 2.824644866830826e-07 1.002629800678504
data_02_0401 2.8153111478043374e-07 2.8172360974302977e-07 1.0006837431192859
data_02_0417 2.810858764680368e-07 2.8153111478043374e-07 1.0015839938953588
data_02_0392 2.739812566500345e-07 2.810858764680368e-07 1.025931043257741
data_02_0410 2.7025876923260697e-07 2.739812566500345e-07 1.0137737895721106
data_02_0415 2.650179609765315e-07 2.7025876923260697e-07 1.0197752946131058
data_02_0405 2.6433964094778335e-07 2.650179609765315e-07 1.002566092721909
data_02_0414 2.5775849616492943e-07 2.6433964094778335e-07 1.0255322128301172
data_02_0396 2.5700836086138224e-07 2.5775849616492943e-07 1.002918719457348
data_02_0391 2.485179358769336e-07 2.5700836086138224e-07 1.034164234281477
data_02_0387 2.4757391546634e-07 2.485179358769336e-07 1.0038130851096143
data_02_0383 2.436324550470703e-07 2.4757391546634e-07 1.016177895586646
data_02_0400 2.411211184995836e-07 2.436324550470703e-07 1.0104152492453333
data_02_0395 2.386285802215298e-07 2.411211184995836e-07 1.0104452629929737
data_02_0399 2.307242924224641e-07 2.386285802215298e-07 1.0342585850673784
data_02_0386 2.2832964938605536e-07 2.307242924224641e-07 1.0104876569593462
data_02_0409 2.209925804911588e-07 2.2832964938605536e-07 1.0332005213866902
data_02_0404 2.200632471515531e-07 2.209925804911588e-07 1.0042230283867695
data_02_0390 2.1855319254474957e-07 2.200632471515531e-07 1.006909323031254
data_02_0382 2.173831522586647e-07 2.1855319254474957e-07 1.0053823871534102
data_02_0381 2.0686221829839521e-07 2.173831522586647e-07 1.0508596207021874
data_02_0413 2.0388155188199395e-07 2.0686221829839521e-07 1.0146195984329494
data_02_0385 2.014030584110105e-07 2.0388155188199395e-07 1.0123061362152976
data_02_0360 1.9864899804127048e-07 2.014030584110105e-07 1.0138639529869053
data_02_0375 1.9263548760005817e-07 1.9864899804127048e-07 1.03121704373442
data_02_0370 1.9203022481550336e-07 1.9263548760005817e-07 1.003151914159015
data_02_0365 1.911260702481893e-07 1.9203022481550336e-07 1.004730671049428
data_02_0380 1.896072463025197e-07 1.911260702481893e-07 1.0080103686714923
data_02_0369 1.7565562213219791e-07 1.896072463025197e-07 1.0794260041379253
data_02_0355 1.7104081410688557e-07 1.7565562213219791e-07 1.0269807416984609
data_02_0374 1.6999505684652706e-07 1.7104081410688557e-07 1.0061516921712768
data_02_0364 1.694072422991218e-07 1.6999505684652706e-07 1.0034698312741988
data_02_0350 1.6598945975060975e-07 1.694072422991218e-07 1.0205903588917458
data_02_0363 1.5748137975093087e-07 1.6598945975060975e-07 1.0540259427059573
data_02_0373 1.5391538880364847e-07 1.5748137975093087e-07 1.0231685146949898
data_02_0354 1.5253345660434522e-07 1.5391538880364847e-07 1.0090598628659404
data_02_0358 1.5111666453375018e-07 1.5253345660434522e-07 1.0093754853243113
data_02_0349 1.452482237806647e-07 1.5111666453375018e-07 1.0404028400508858
data_02_0368 1.4209152242578323e-07 1.452482237806647e-07 1.0222159725013171
data_02_0359 1.372053520546384e-07 1.4209152242578323e-07 1.0356120974726921
data_02_0379 1.370664689286464e-07 1.372053520546384e-07 1.0010132538401082
data_02_0367 1.3261333165400434e-07 1.370664689286464e-07 1.033579861233413
data_02_0378 1.325185368645295e-07 1.3261333165400434e-07 1.0007153322977882
data_02_0362 1.2827295093725083e-07 1.325185368645295e-07 1.0330980607856721
data_02_0372 1.2748598803605932e-07 1.2827295093725083e-07 1.006172936440426
data_02_0344 1.2738010149606113e-07 1.2748598803605932e-07 1.0008312643713937
data_02_0353 1.2704938988446473e-07 1.2738010149606113e-07 1.0026030161333095
data_02_0345 1.2601085957679468e-07 1.2704938988446473e-07 1.0082415937099227
data_02_0343 1.216743623001515e-07 1.2601085957679468e-07 1.0356401890641984
data_02_0357 1.1836471180962112e-07 1.216743623001515e-07 1.0279614628374516
data_02_0348 1.1627477364692665e-07 1.1836471180962112e-07 1.017974132282903
data_02_0371 1.1505417473274633e-07 1.1627477364692665e-07 1.0106089059090257
data_02_0366 1.1444872514963532e-07 1.1505417473274633e-07 1.005290138289609
data_02_0352 1.1166505000103068e-07 1.1444872514963532e-07 1.0249287950758
data_02_0356 1.1066799539378786e-07 1.1166505000103068e-07 1.0090094214112673
data_02_0377 1.0535397743098524e-07 1.1066799539378786e-07 1.050439652041459
data_02_0351 1.024389242356505e-07 1.0535397743098524e-07 1.028456499490652
data_02_0361 1.0153767963564204e-07 1.024389242356505e-07 1.008875962137824
data_02_0347 1.0134012736906253e-07 1.0153767963564204e-07 1.0019493982463636
data_02_0376 9.832303420881312e-08 1.0134012736906253e-07 1.030685517229278
data_02_0342 9.289929706956034e-08 9.832303420881312e-08 1.0583829728570673
data_02_0329 9.160053376660972e-08 9.289929706956034e-08 1.0141785560579784
data_02_0346 9.11694229508457e-08 9.160053376660972e-08 1.0047286776839255
data_02_0111 9.091592748881475e-08 9.11694229508457e-08 1.0027882404000348
data_02_0074 9.069444008367644e-08 9.091592748881475e-08 1.0024421277085338
data_02_0137 9.032979636101275e-08 9.069444008367644e-08 1.004036804436117
data_02_0095 8.745287729949787e-08 9.032979636101275e-08 1.0328967913961522
data_02_0148 8.648924272107816e-08 8.745287729949787e-08 1.011141669739523
data_02_0125 8.424422939458488e-08 8.648924272107816e-08 1.0266488677340502
data_02_0334 8.326863601499629e-08 8.424422939458488e-08 1.0117162166486418
data_02_0112 8.320700382561938e-08 8.326863601499629e-08 1.0007407091535956
data_02_0075 8.138272522950375e-08 8.320700382561938e-08 1.0224160421142332
data_02_0126 8.103538161033434e-08 8.138272522950375e-08 1.0042863205215673
data_02_0096 8.04529530328609e-08 8.103538161033434e-08 1.007239368544192
data_02_0076 8.03161118419718e-08 8.04529530328609e-08 1.0017037825630595
data_02_0341 7.826190779854336e-08 8.03161118419718e-08 1.0262478145653724
data_02_0138 7.73441165303848e-08 7.826190779854336e-08 1.011866335919656
data_02_0313 7.682170527795377e-08 7.73441165303848e-08 1.0068003079408463
data_02_0097 7.680681819199071e-08 7.682170527795377e-08 1.000193825057638
data_02_0176 7.66773188177075e-08 7.680681819199071e-08 1.001688887617355
data_02_0328 7.592391064574658e-08 7.66773188177075e-08 1.0099232002876701
data_02_0113 7.555670244097317e-08 7.592391064574658e-08 1.004860034820343
data_02_0327 7.544625046424394e-08 7.555670244097317e-08 1.0014639823191953
data_02_0127 7.131655318571188e-08 7.544625046424394e-08 1.057906574197131
data_02_0312 6.9114916843674e-08 7.131655318571188e-08 1.0318547202627415
data_02_0139 6.909313831848126e-08 6.9114916843674e-08 1.0003152053260682
data_02_0340 6.902051779125408e-08 6.909313831848126e-08 1.0010521585399694
data_02_0333 6.774037272934449e-08 6.902051779125408e-08 1.018897815443449
data_02_0307 6.766458522400373e-08 6.774037272934449e-08 1.0011200468471042
data_02_0185 6.763927671203422e-08 6.766458522400373e-08 1.0003741688734678
data_02_0150 6.745632682542282e-08 6.763927671203422e-08 1.0027121234615233
data_02_0177 6.735833860900495e-08 6.745632682542282e-08 1.0014547303042414
data_02_0193 6.662207488596927e-08 6.735833860900495e-08 1.011051347834721
data_02_0208 6.626149627555743e-08 6.662207488596927e-08 1.0054417517061842
data_02_0098 6.592270735973779e-08 6.626149627555743e-08 1.0051391838925985
data_02_0302 6.582978203072585e-08 6.592270735973779e-08 1.0014116001321192
data_02_0140 6.57620486951575e-08 6.582978203072585e-08 1.0010299760562864
data_02_0077 6.529940874175905e-08 6.57620486951575e-08 1.007084902640207
data_02_0201 6.30318615248224e-08 6.529940874175905e-08 1.0359746192176742
data_02_0114 6.245921319803077e-08 6.30318615248224e-08 1.0091683563956533
data_02_0318 6.239185556929697e-08 6.245921319803077e-08 1.001079590086225
data_02_0178 6.168574049096606e-08 6.239185556929697e-08 1.0114469741744336
data_02_0326 6.156214131519181e-08 6.168574049096606e-08 1.0020077140452512
data_02_0215 6.154482629620753e-08 6.156214131519181e-08 1.0002813399602584
data_02_0311 6.138472410942559e-08 6.154482629620753e-08 1.0026081763680577
data_02_0099 6.0730378062034e-08 6.138472410942559e-08 1.0107746084953266
data_02_0301 6.069713498226437e-08 6.0730378062034e-08 1.0005476877908546
data_02_0186 6.061854749380371e-08 6.069713498226437e-08 1.0012964264521298
data_02_0078 6.055514982039741e-08 6.061854749380371e-08 1.0010469410709797
data_02_0079 5.949129213300473e-08 6.055514982039741e-08 1.0178825782606002
data_02_0306 5.91244260488755e-08 5.949129213300473e-08 1.0062049834331745
data_02_0332 5.727316434403084e-08 5.91244260488755e-08 1.0323233703960273
data_02_0141 5.719419670667194e-08 5.727316434403084e-08 1.001380693180532
data_02_0115 5.716266762862783e-08 5.719419670667194e-08 1.0005515676463692
data_02_0128 5.69555024196704e-08 5.716266762862783e-08 1.0036373168553752
data_02_0168 5.6797071300329324e-08 5.69555024196704e-08 1.0027894240972977
data_02_0100 5.659626681287352e-08 5.6797071300329324e-08 1.0035480164817183
data_02_0216 5.52523403110799e-08 5.659626681287352e-08 1.024323431265121
data_02_0217 5.47644941180883e-08 5.52523403110799e-08 1.008908074489644
data_02_0278 5.440684246618147e-08 5.47644941180883e-08 1.00657365205726
data_02_0194 5.336290930438235e-08 5.440684246618147e-08 1.019562898189162
data_02_0169 5.33316427844242e-08 5.336290930438235e-08 1.0005862658325477
data_02_0325 5.3061671911245655e-08 5.33316427844242e-08 1.0050878697081034
data_02_0284 5.276553591771626e-08 5.3061671911245655e-08 1.0056122995508128
data_02_0151 5.236707142103079e-08 5.276553591771626e-08 1.0076090658857324
data_02_0116 5.235342965295829e-08 5.236707142103079e-08 1.0002605706667724
data_02_0129 5.1837413048961746e-08 5.235342965295829e-08 1.0099545207533631
data_02_0209 5.183708040770919e-08 5.1837413048961746e-08 1.0000064170522325
data_02_0101 5.1669753255601975e-08 5.183708040770919e-08 1.0032383965775775
data_02_0331 5.107873307696695e-08 5.1669753255601975e-08 1.0115707681657737
data_02_0117 5.10137341413657e-08 5.107873307696695e-08 1.0012741458098544
data_02_0170 4.9929659581962395e-08 5.10137341413657e-08 1.0217120358616452
data_02_0202 4.97999691571996e-08 4.9929659581962395e-08 1.0026042270097277
data_02_0187 4.912206238348923e-08 4.97999691571996e-08 1.0138004542321135
data_02_0195 4.909545697928861e-08 4.912206238348923e-08 1.0005419117335408
data_02_0179 4.884858914365607e-08 4.909545697928861e-08 1.005053735224707
data_02_0152 4.8736593849818535e-08 4.884858914365607e-08 1.0022979712981717
data_02_0338 4.863241454752601e-08 4.8736593849818535e-08 1.0021421782829787
data_02_0080 4.776368500103665e-08 4.863241454752601e-08 1.0181880762857913
data_02_0300 4.730771729074619e-08 4.776368500103665e-08 1.0096383367535606
data_02_0289 4.716165846814214e-08 4.730771729074619e-08 1.0030969823231028
data_02_0305 4.6948542791565873e-08 4.716165846814214e-08 1.0045393459286356
data_02_0130 4.6783812999108644e-08 4.6948542791565873e-08 1.0035210852193763
data_02_0317 4.676572950679156e-08 4.6783812999108644e-08 1.000386682566653
data_02_0142 4.6327681479980905e-08 4.676572950679156e-08 1.0094554273561034
data_02_0203 4.572493317313839e-08 4.6327681479980905e-08 1.0131820489394747
data_02_0180 4.515687271376087e-08 4.572493317313839e-08 1.0125797121288342
data_02_0236 4.4905309913886836e-08 4.515687271376087e-08 1.0056020724577215
data_02_0210 4.400392577456567e-08 4.4905309913886836e-08 1.0204841755242249
data_02_0081 4.364525646736216e-08 4.400392577456567e-08 1.0082178302118976
data_02_0131 4.354492864855354e-08 4.364525646736216e-08 1.0023040069629774
data_02_0310 4.337831336396784e-08 4.354492864855354e-08 1.0038409811646596
data_02_0153 4.317816058687787e-08 4.337831336396784e-08 1.0046355095809893
data_02_0339 4.314213555749596e-08 4.317816058687787e-08 1.0008350312036338
data_02_0196 4.23775629408846e-08 4.314213555749596e-08 1.01804192038315
data_02_0188 4.2228768481487626e-08 4.23775629408846e-08 1.0035235330024412
data_02_0229 4.208906135417072e-08 4.2228768481487626e-08 1.0033193215249279
data_02_0290 4.1538868891419455e-08 4.208906135417072e-08 1.0132452442118594
data_02_0172 4.115915270056277e-08 4.1538868891419455e-08 1.0092255589812347
data_02_0102 4.110164335007903e-08 4.115915270056277e-08 1.0013991983238701
data_02_0330 4.090697832851984e-08 4.110164335007903e-08 1.004758724049375
data_02_0337 4.0521742966527407e-08 4.090697832851984e-08 1.0095068803508933
data_02_0309 3.974405981438906e-08 4.0521742966527407e-08 1.0195672801362077
data_02_0232 3.9588237640243336e-08 3.974405981438906e-08 1.0039360725163307
data_02_0171 3.952384984565909e-08 3.9588237640243336e-08 1.0016290871166569
data_02_0230 3.945230415583387e-08 3.952384984565909e-08 1.0018134730367743
data_02_0242 3.921650266905879e-08 3.945230415583387e-08 1.0060128127376622
data_02_0316 3.909765867391486e-08 3.921650266905879e-08 1.003039670383721
data_02_0239 3.8972549190214436e-08 3.909765867391486e-08 1.0032101950295782
data_02_0237 3.875843721078654e-08 3.8972549190214436e-08 1.0055242676133576
data_02_0299 3.8590618949182734e-08 3.875843721078654e-08 1.004348680227824
data_02_0324 3.8269657665395815e-08 3.8590618949182734e-08 1.0083868344627795
data_02_0249 3.8109111998342565e-08 3.8269657665395815e-08 1.0042127895045214
data_02_0143 3.81075391790093e-08 3.8109111998342565e-08 1.000041273180246
data_02_0118 3.8101098363225644e-08 3.81075391790093e-08 1.0001690454097216
data_02_0189 3.800245193459355e-08 3.8101098363225644e-08 1.0025957911558412
data_02_0082 3.7742316519203115e-08 3.800245193459355e-08 1.0068924072336174
data_02_0238 3.756009448605654e-08 3.7742316519203115e-08 1.0048514796259158
data_02_0323 3.751518306258622e-08 3.756009448605654e-08 1.001197153253801
data_02_0181 3.706270956353247e-08 3.751518306258622e-08 1.0122083221756393
data_02_0231 3.682098674305962e-08 3.706270956353247e-08 1.0065648110453862
data_02_0103 3.677846288057728e-08 3.682098674305962e-08 1.0011562164144927
data_02_0255 3.6284794666359044e-08 3.677846288057728e-08 1.0136053743381366
data_02_0204 3.616483795118732e-08 3.6284794666359044e-08 1.0033169432511666
data_02_0272 3.5992364955312874e-08 3.616483795118732e-08 1.0047919328471075
data_02_0154 3.5938992759256385e-08 3.5992364955312874e-08 1.001485077681893
data_02_0083 3.589503422205851e-08 3.5938992759256385e-08 1.0012246411836783
data_02_0315 3.540837869950097e-08 3.589503422205851e-08 1.013744078108959
data_02_0271 3.540526187431156e-08 3.540837869950097e-08 1.0000880328240607
data_02_0132 3.5391466146159544e-08 3.540526187431156e-08 1.00038980380454
data_02_0211 3.531823528704799e-08 3.5391466146159544e-08 1.0020734574792984
data_02_0283 3.5082251142443155e-08 3.531823528704799e-08 1.0067265964104377
data_02_0304 3.4693944862226126e-08 3.5082251142443155e-08 1.0111923357738373
data_02_0266 3.438015840623516e-08 3.4693944862226126e-08 1.0091269636481388
data_02_0282 3.424046100777206e-08 3.438015840623516e-08 1.0040798924532992
data_02_0308 3.4026211474780814e-08 3.424046100777206e-08 1.0062966026397044
data_02_0259 3.394382336929748e-08 3.4026211474780814e-08 1.0024271899068935
data_02_0119 3.3927979328715124e-08 3.394382336929748e-08 1.0004669903983627
data_02_0253 3.3906648726021295e-08 3.3927979328715124e-08 1.0006290979349268
data_02_0197 3.380955191293169e-08 3.3906648726021295e-08 1.002871875183074
data_02_0251 3.359312442011387e-08 3.380955191293169e-08 1.0064426127832347
data_02_0144 3.324770433476301e-08 3.359312442011387e-08 1.010389291298819
data_02_0182 3.304672707119537e-08 3.324770433476301e-08 1.0060816087213313
data_02_0336 3.291690626091583e-08 3.304672707119537e-08 1.0039438946433943
data_02_0264 3.273172069285836e-08 3.291690626091583e-08 1.0056576789773803
data_02_0104 3.252849058747102e-08 3.273172069285836e-08 1.0062477570190613
data_02_0298 3.249594062884612e-08 3.252849058747102e-08 1.0010016622998137
data_02_0250 3.242477208092828e-08 3.249594062884612e-08 1.0021948819791302
data_02_0256 3.23632897138752e-08 3.242477208092828e-08 1.0018997564090872
data_02_0173 3.2333302089077816e-08 3.23632897138752e-08 1.0009274532095351
data_02_0258 3.22693739170363e-08 3.2333302089077816e-08 1.0019810787840469
data_02_0243 3.21705518617872e-08 3.22693739170363e-08 1.0030718172219633
data_02_0190 3.214098309363574e-08 3.21705518617872e-08 1.0009199708691334
data_02_0262 3.193217453326343e-08 3.214098309363574e-08 1.0065391274920785
data_02_0269 3.151817992901874e-08 3.193217453326343e-08 1.0131351050465807
data_02_0120 3.150506414970302e-08 3.151817992901874e-08 1.0004163070182432
data_02_0191 3.141980375375175e-08 3.150506414970302e-08 1.0027135877938476
data_02_0261 3.124496701944196e-08 3.141980375375175e-08 1.0055956767117404
data_02_0218 3.055804546327236e-08 3.124496701944196e-08 1.022479237325411
data_02_0252 3.0186530033837795e-08 3.055804546327236e-08 1.0123073247908294
data_02_0322 3.0117442642990743e-08 3.0186530033837795e-08 1.0022939328437015
data_02_0314 3.009793464370751e-08 3.0117442642990743e-08 1.0006481507623086
data_02_0265 3.005155814323877e-08 3.009793464370751e-08 1.0015432311445445
data_02_0268 2.9905636453243234e-08 3.005155814323877e-08 1.0048794042629283
data_02_0155 2.9791696455670336e-08 2.9905636453243234e-08 1.0038245555348766
data_02_0267 2.959674789256824e-08 2.9791696455670336e-08 1.0065868237892126
data_02_0198 2.8804755286080425e-08 2.959674789256824e-08 1.027495203435057
data_02_0084 2.8541810283806125e-08 2.8804755286080425e-08 1.0092126252560611
data_02_0303 2.846204854188475e-08 2.8541810283806125e-08 1.0028023893573226
data_02_0257 2.8461991094541266e-08 2.846204854188475e-08 1.0000020183880773
data_02_0281 2.8445644744530024e-08 2.8461991094541266e-08 1.000574652118384
data_02_0219 2.8276179870352716e-08 2.8445644744530024e-08 1.0059932025809113
data_02_0174 2.7830212773979406e-08 2.8276179870352716e-08 1.0160245665383585
data_02_0212 2.7693025488395505e-08 2.7830212773979406e-08 1.0049538569067287
data_02_0133 2.7647622286096888e-08 2.7693025488395505e-08 1.0016422100182354
data_02_0297 2.7424003583534e-08 2.7647622286096888e-08 1.0081541231527973
data_02_0205 2.719201753826118e-08 2.7424003583534e-08 1.008531402458328
data_02_0321 2.6960865849891603e-08 2.719201753826118e-08 1.0085736003308108
data_02_0277 2.6917632473839484e-08 2.6960865849891603e-08 1.0016061359071655
data_02_0105 2.682414247631298e-08 2.6917632473839484e-08 1.0034852930567701
data_02_0199 2.6596847855458792e-08 2.682414247631298e-08 1.008545923264645
data_02_0145 2.6322022042612023e-08 2.6596847855458792e-08 1.0104409080883627
data_02_0134 2.5886623017233462e-08 2.6322022042612023e-08 1.0168194601933478
data_02_0244 2.586085837913866e-08 2.5886623017233462e-08 1.000996279308176
data_02_0288 2.5814320024567325e-08 2.586085837913866e-08 1.001802811560678
data_02_0280 2.5567505217879908e-08 2.5814320024567325e-08 1.009653456783684
data_02_0200 2.5161855134809113e-08 2.5567505217879908e-08 1.0161216285880932
data_02_0286 2.5066605124049343e-08 2.5161855134809113e-08 1.003799876779819
data_02_0085 2.5060032212241507e-08 2.5066605124049343e-08 1.0002622866464084
data_02_0270 2.501478009837129e-08 2.5060032212241507e-08 1.0018090150579881
data_02_0263 2.4817553536481463e-08 2.501478009837129e-08 1.0079470589879016
data_02_0156 2.461345700260151e-08 2.4817553536481463e-08 1.008292071034897
data_02_0276 2.4560455607394448e-08 2.461345700260151e-08 1.0021579972315784
data_02_0285 2.4538554136531483e-08 2.4560455607394448e-08 1.0008925330620992
data_02_0121 2.405860616838822e-08 2.4538554136531483e-08 1.0199491177828035
data_02_0296 2.405840293667572e-08 2.405860616838822e-08 1.0000084474315705
data_02_0207 2.3903482341567128e-08 2.405840293667572e-08 1.0064810889432285
data_02_0449 2.3732450803006916e-08 2.3903482341567128e-08 1.0072066530330084
data_02_0275 2.3517630583804178e-08 2.3732450803006916e-08 1.0091344329284038
data_02_0183 2.350059371293589e-08 2.3517630583804178e-08 1.0007249549129011
data_02_0206 2.3061833868340542e-08 2.350059371293589e-08 1.0190253666339033
data_02_0220 2.2994991724143692e-08 2.3061833868340542e-08 1.0029068131443017
data_02_0146 2.274100321319186e-08 2.2994991724143692e-08 1.0111687469796626
data_02_0221 2.25730278406735e-08 2.274100321319186e-08 1.0074414196316053
data_02_0295 2.2222329238385043e-08 2.25730278406735e-08 1.0157813611042485
data_02_0106 2.220331189048818e-08 2.2222329238385043e-08 1.0008565095149167
data_02_0287 2.1985352902847122e-08 2.220331189048818e-08 1.00991382711045
data_02_0335 2.1740529587537427e-08 2.1985352902847122e-08 1.0112611477252162
data_02_0157 2.1600877912403635e-08 2.1740529587537427e-08 1.0064650925624463
data_02_0086 2.128119346174364e-08 2.1600877912403635e-08 1.0150219230530788
data_02_0135 2.0605336844278892e-08 2.128119346174364e-08 1.0328000761439822
data_02_0245 2.0214721336813185e-08 2.0605336844278892e-08 1.019323318929673
data_02_0147 2.012244450124676e-08 2.0214721336813185e-08 1.004585766682607
data_02_0122 2.00493566983083e-08 2.012244450124676e-08 1.0036453939165355
data_02_0107 1.938732662337907e-08 2.00493566983083e-08 1.0341475690686972
data_02_0184 1.9329828054146482e-08 1.938732662337907e-08 1.0029746032438325
data_02_0233 1.924222296865813e-08 1.9329828054146482e-08 1.004552752851427
data_02_0175 1.9132002251777512e-08 1.924222296865813e-08 1.0057610654352906
data_02_0192 1.906093847893854e-08 1.9132002251777512e-08 1.0037282410264057
data_02_0214 1.8668537337049456e-08 1.906093847893854e-08 1.0210193832973902
data_02_0246 1.8085101445965826e-08 1.8668537337049456e-08 1.0322605816078392
data_02_0224 1.7865150745960694e-08 1.8085101445965826e-08 1.012311718111579
data_02_0123 1.7784121263186588e-08 1.7865150745960694e-08 1.0045562826284726
data_02_0234 1.7743078014120953e-08 1.7784121263186588e-08 1.0023131978021498
data_02_0087 1.768952147617265e-08 1.7743078014120953e-08 1.0030275854561947
data_02_0136 1.6856231726262295e-08 1.768952147617265e-08 1.049435114766017
data_02_0222 1.685224351916152e-08 1.6856231726262295e-08 1.0002366573386055
data_02_0223 1.6841140219094164e-08 1.685224351916152e-08 1.0006592962188372
data_02_0108 1.6692853350424656e-08 1.6841140219094164e-08 1.0088832547412114
data_02_0124 1.6497853114554667e-08 1.6692853350424656e-08 1.0118197340293906
data_02_0274 1.6468220100291053e-08 1.6497853114554667e-08 1.0017994060124986
data_02_0158 1.6092249696010172e-08 1.6468220100291053e-08 1.0233634458440013
data_02_0240 1.5068606562892854e-08 1.6092249696010172e-08 1.067932169364491
data_02_0320 1.503695553233086e-08 1.5068606562892854e-08 1.0021048828995964
data_02_0247 1.417393337054734e-08 1.503695553233086e-08 1.060887979308329
data_02_0248 1.4169020848167241e-08 1.417393337054734e-08 1.0003467086704678
data_02_0088 1.3875326621797155e-08 1.4169020848167241e-08 1.021166653180525
data_02_0235 1.346834787193108e-08 1.3875326621797155e-08 1.0302174218943547
data_02_0160 1.3120044854721955e-08 1.346834787193108e-08 1.0265473953074002
data_02_0161 1.2298610270804686e-08 1.3120044854721955e-08 1.0667908459435655
data_02_0260 1.2167262041760474e-08 1.2298610270804686e-08 1.0107952165896812
data_02_0159 1.1999987829094322e-08 1.2167262041760474e-08 1.01393953186024
data_02_0294 1.1975931948430658e-08 1.1999987829094322e-08 1.0020086854841235
data_02_0273 1.1970562880034427e-08 1.1975931948430658e-08 1.0004485226342352
data_02_0241 1.1800058078610208e-08 1.1970562880034427e-08 1.0144494883235609
data_02_0109 1.1756312711164799e-08 1.1800058078610208e-08 1.0037210108747672
data_02_0254 1.1681364901973924e-08 1.1756312711164799e-08 1.0064160147225785
data_02_0089 1.0944847533250584e-08 1.1681364901973924e-08 1.0672935247828526
data_02_0319 1.0710358622772656e-08 1.0944847533250584e-08 1.0218936562945102
data_02_0291 9.647637171764193e-09 1.0710358622772656e-08 1.1101535466236994
data_02_0110 9.297541194942495e-09 9.647637171764193e-09 1.0376546841235978
data_02_0090 8.946145284495091e-09 9.297541194942495e-09 1.0392790301602213
data_02_0440 8.921394664214538e-09 8.946145284495091e-09 1.002774299446681
data_02_0292 8.845596524401773e-09 8.921394664214538e-09 1.0085690252323476
data_02_0293 8.013189234840733e-09 8.845596524401773e-09 1.10387964955848
data_02_0091 5.400196979058821e-09 8.013189234840733e-09 1.4838698043635659
376
Overlay of the differents results
We are getting to an end. Here are the first actually integrated data
[61]:
fig, ax = plt.subplots()
summed, counted, radial = None, None, None
for i in range(9):
name = ds_names[i]
ds = data[name]
gonioref = goniometers[name]
mg = gonioref.get_mg(position)
mg.radial_range = (0, 95)
images = [i.reshape(-1, 1) for i in ds]
res_mg = mg.integrate1d(images, 50000)
results[name] = res_mg
if summed is None:
summed = res_mg.sum
counted = res_mg.count
else:
summed += res_mg.sum
counted += res_mg.count
radial = res_mg.radial
jupyter.plot1d(res_mg, label="%i %s"%(i, name), calibrant=LaB6, ax=ax )
ax.plot(radial, summed/counted, label="Merged")
ax.legend()
fig.show()
<ipython-input-61-e8eac1ef1cce>:22: RuntimeWarning: invalid value encountered in divide
ax.plot(radial, summed/counted, label="Merged")
<ipython-input-61-e8eac1ef1cce>:24: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
Multi-Gonio fit
Can we fit everything togeather ? Just assume energy and scale parameter of the goniometer are the same for all modules and fit everything.
[62]:
from scipy.optimize import minimize
class MultiGoniometer:
def __init__(self, list_of_goniometers,
param_name_split,
param_name_common):
self.nb_gonio = len(list_of_goniometers)
self.goniometers = list_of_goniometers
self.names_split = param_name_split
self.names_common = param_name_common
self.param = None
def init_param(self):
param = []
for gonio in self.goniometers:
param += list(gonio.param[:len(self.names_split)])
param += list(gonio.param[len(self.names_split):])
self.param = numpy.array(param)
def residu2(self, param):
"Actually performs the calulation of the average of the error squared"
sumsquare = 0.0
npt = 0
for idx, gonio in enumerate(self.goniometers):
gonio_param = numpy.concatenate((param[len(self.names_split)*idx:len(self.names_split)*(1+idx)],
param[len(self.names_split)*len(self.goniometers):]))
sumsquare += gonio.residu2(gonio_param)
return sumsquare
def chi2(self, param=None):
"""Calculate the average of the square of the error for a given parameter set
"""
if param is not None:
return self.residu2(param)
else:
if self.param is None:
self.init_param()
return self.residu2(self.param)
def refine2(self, method="slsqp", **options):
"""Geometry refinement tool
See https://docs.scipy.org/doc/scipy-0.18.1/reference/generated/scipy.optimize.minimize.html
:param method: name of the minimizer
:param options: options for the minimizer
"""
if method.lower() in ["simplex", "nelder-mead"]:
method = "Nelder-Mead"
former_error = self.chi2()
print("Cost function before refinement: %s" % former_error)
param = numpy.asarray(self.param, dtype=numpy.float64)
print(param)
res = minimize(self.residu2, param, method=method,
tol=1e-12,
options=options)
print(res)
newparam = res.x
new_error = res.fun
print("Cost function after refinement: %s" % new_error)
if new_error < former_error:
self.param = newparam
return self.param
def integrate(self, list_of_dataset, npt=50000, radial_range=(0,100)):
summed = None
counted = None
param = self.param
for idx, ds in enumerate(list_of_dataset):
gonio = self.goniometers[idx]
gonio_param = numpy.concatenate((param[len(self.names_split)*idx:len(self.names_split)*(1+idx)],
param[len(self.names_split)*len(self.goniometers):]))
print(gonio_param)
gonio.param = gonio_param
mg = gonio.get_mg(position)
mg.radial_range = radial_range
images = [i.reshape(-1, 1) for i in ds]
res_mg = mg.integrate1d(images, 50000)
if summed is None:
summed = res_mg.sum
counted = res_mg.count
else:
summed += res_mg.sum
counted += res_mg.count
radial = res_mg.radial
res = Integrate1dResult(radial, summed/numpy.maximum(counted, 1e-10))
res._set_unit(res_mg.unit)
res._set_count(counted)
res._set_sum(summed)
return res
[63]:
multigonio = MultiGoniometer([goniometers[ds_names[i]] for i in range(9)],
["dist", "poni1", "poni2", "rot1", "offset"],
["scale", "nrj"])
[64]:
%time print(multigonio.chi2())
multigonio.param = numpy.array([ 7.20594053e-01, 3.22408604e-02, 4.05228023e-03, -2.75578440e-05,
-8.27999414e+01, 7.20612302e-01, 3.36369797e-02, 4.02094516e-03,
-1.74996556e-05, -7.71999791e+01, 7.20636130e-01, 3.47920978e-02,
4.01341931e-03, -1.21330600e-05, -7.15999090e+01, 7.20757808e-01,
3.33850817e-02, 3.95036100e-03, 3.46517345e-05, -6.57999267e+01,
7.20813915e-01, 3.22167822e-02, 3.97128822e-03, 2.00055269e-05,
-6.00000525e+01, 7.20881596e-01, 3.33801850e-02, 3.97760147e-03,
1.47074593e-05, -5.43998157e+01, 7.21048510e-01, 3.22346939e-02,
4.02104962e-03, -1.69519259e-05, -4.85998856e+01, 7.21074630e-01,
3.08484557e-02, 4.09385968e-03, -6.91378973e-05, -4.27999030e+01,
7.21154891e-01, 3.20619921e-02, 4.24950906e-03, -1.81328256e-04,
-3.71999987e+01, 9.99038595e-01, 1.70266104e+01])
%time print(multigonio.chi2())
8.603025165471922e-09
CPU times: user 387 ms, sys: 20 ms, total: 407 ms
Wall time: 406 ms
6.1395059557657925e-09
CPU times: user 382 ms, sys: 24 ms, total: 406 ms
Wall time: 405 ms
[65]:
%time multigonio.refine2()
Cost function before refinement: 6.1395059557657925e-09
[ 7.20594053e-01 3.22408604e-02 4.05228023e-03 -2.75578440e-05
-8.27999414e+01 7.20612302e-01 3.36369797e-02 4.02094516e-03
-1.74996556e-05 -7.71999791e+01 7.20636130e-01 3.47920978e-02
4.01341931e-03 -1.21330600e-05 -7.15999090e+01 7.20757808e-01
3.33850817e-02 3.95036100e-03 3.46517345e-05 -6.57999267e+01
7.20813915e-01 3.22167822e-02 3.97128822e-03 2.00055269e-05
-6.00000525e+01 7.20881596e-01 3.33801850e-02 3.97760147e-03
1.47074593e-05 -5.43998157e+01 7.21048510e-01 3.22346939e-02
4.02104962e-03 -1.69519259e-05 -4.85998856e+01 7.21074630e-01
3.08484557e-02 4.09385968e-03 -6.91378973e-05 -4.27999030e+01
7.21154891e-01 3.20619921e-02 4.24950906e-03 -1.81328256e-04
-3.71999987e+01 9.99038595e-01 1.70266104e+01]
message: Optimization terminated successfully
success: True
status: 0
fun: 6.1395059557657925e-09
x: [ 7.206e-01 3.224e-02 ... 9.990e-01 1.703e+01]
nit: 1
jac: [ 5.910e-08 -8.339e-08 ... -4.210e-08 6.458e-07]
nfev: 48
njev: 1
Cost function after refinement: 6.1395059557657925e-09
CPU times: user 19.2 s, sys: 464 ms, total: 19.7 s
Wall time: 19.7 s
[65]:
array([ 7.20594053e-01, 3.22408604e-02, 4.05228023e-03, -2.75578440e-05,
-8.27999414e+01, 7.20612302e-01, 3.36369797e-02, 4.02094516e-03,
-1.74996556e-05, -7.71999791e+01, 7.20636130e-01, 3.47920978e-02,
4.01341931e-03, -1.21330600e-05, -7.15999090e+01, 7.20757808e-01,
3.33850817e-02, 3.95036100e-03, 3.46517345e-05, -6.57999267e+01,
7.20813915e-01, 3.22167822e-02, 3.97128822e-03, 2.00055269e-05,
-6.00000525e+01, 7.20881596e-01, 3.33801850e-02, 3.97760147e-03,
1.47074593e-05, -5.43998157e+01, 7.21048510e-01, 3.22346939e-02,
4.02104962e-03, -1.69519259e-05, -4.85998856e+01, 7.21074630e-01,
3.08484557e-02, 4.09385968e-03, -6.91378973e-05, -4.27999030e+01,
7.21154891e-01, 3.20619921e-02, 4.24950906e-03, -1.81328256e-04,
-3.71999987e+01, 9.99038595e-01, 1.70266104e+01])
[66]:
LaB6_new = get_calibrant("LaB6")
LaB6_new.wavelength = 1e-10*hc/multigonio.param[-1]
print(LaB6,"\n", LaB6_new)
LaB6 Calibrant with 109 reflections at wavelength 7.281789822635113e-11
LaB6 Calibrant with 109 reflections at wavelength 7.281789829007909e-11
[67]:
%time res = multigonio.integrate([data[ds_names[i]] for i in range(9)])
[ 7.20594053e-01 3.22408604e-02 4.05228023e-03 -2.75578440e-05
-8.27999414e+01 9.99038595e-01 1.70266104e+01]
[ 7.20612302e-01 3.36369797e-02 4.02094516e-03 -1.74996556e-05
-7.71999791e+01 9.99038595e-01 1.70266104e+01]
[ 7.20636130e-01 3.47920978e-02 4.01341931e-03 -1.21330600e-05
-7.15999090e+01 9.99038595e-01 1.70266104e+01]
[ 7.20757808e-01 3.33850817e-02 3.95036100e-03 3.46517345e-05
-6.57999267e+01 9.99038595e-01 1.70266104e+01]
[ 7.20813915e-01 3.22167822e-02 3.97128822e-03 2.00055269e-05
-6.00000525e+01 9.99038595e-01 1.70266104e+01]
[ 7.20881596e-01 3.33801850e-02 3.97760147e-03 1.47074593e-05
-5.43998157e+01 9.99038595e-01 1.70266104e+01]
[ 7.21048510e-01 3.22346939e-02 4.02104962e-03 -1.69519259e-05
-4.85998856e+01 9.99038595e-01 1.70266104e+01]
[ 7.21074630e-01 3.08484557e-02 4.09385968e-03 -6.91378973e-05
-4.27999030e+01 9.99038595e-01 1.70266104e+01]
[ 7.21154891e-01 3.20619921e-02 4.24950906e-03 -1.81328256e-04
-3.71999987e+01 9.99038595e-01 1.70266104e+01]
CPU times: user 6min 5s, sys: 913 ms, total: 6min 6s
Wall time: 35.7 s
[68]:
ax = jupyter.plot1d(res, calibrant=LaB6_new)
ax.figure.show()
<ipython-input-68-5d6e268af85e>:2: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
ax.figure.show()
[69]:
fig, ax = plt.subplots()
ax.plot(*calc_fwhm(res, LaB6_new, 10, 95), "o", label="FWHM")
ax.plot(*calc_peak_error(res, LaB6_new, 10, 95), "o", label="error")
ax.set_title("Peak shape & error as function of the angle")
ax.set_xlabel(res.unit.label)
ax.legend()
fig.show()
28454
68 60
<ipython-input-69-c84371e16330>:7: UserWarning: Matplotlib is currently using module://matplotlib_inline.backend_inline, which is a non-GUI backend, so cannot show the figure.
fig.show()
[70]:
print("total run time: ", time.time()-start_time)
total run time: 956.6056206226349
Conclusion
The calibration works and the FWHM of every single peak is pretty small: 0.02°. The geometry has been refined with the wavelength: The goniometer scale parameter refines to 0.999 instead of 1 and the wavelength is fitted with a change at the 5th digit which is pretty precise.