The aim of this document is to explain how to use pyFAI.goniometer for calibrating the position of the detector from the goniometer encoders.
Those data have been acquired at ROBL (ESRF-BM20 German CRG) in winter 2017 by Christoph Henning using a Pilatus 100k detector and LaB6 as calibrant. One hundred and twenty one images have been acquired with the detector moving between 5 and 65 degree with a step size of half a degree. The motor position is registered in the filename.
A prior manual calibration (using pyFAI-calib) has been performed on four images locate at 31.5, 33.5, 35 and 35.5 degrees. Those images were the first with two rings. The control points extrated during this initial calibration has been used as a starting point for this calibration. Then more images have been added to make the model more robust.
In [1]:
# Initialization of the plotting library to be used with the jupyter notebook
%pylab nbagg
In [2]:
#Loading of a few libraries
import time
start_time = time.time()
import os
import glob
import random
import fabio
import pyFAI
from pyFAI.goniometer import GeometryTransformation, GoniometerRefinement, Goniometer
from pyFAI.gui import jupyter
In [3]:
#loading of the list of files, and display of the last one with its headers
image_files = glob.glob("*.cbf")
image_files.sort()
print("List of images: " + ", ".join(image_files) + "." + os.linesep)
fimg = fabio.open(image_files[-1])
print("Image headers:")
for key, value in fimg.header.items():
print("%s: %s"%(key,value))
jupyter.display(fimg.data, label=fimg.filename)
Out[3]:
In [4]:
#Definition of the goniometer translation function:
# The detector rotates vertically, around the horizontal axis, i.e. rot2
goniotrans = GeometryTransformation(param_names = ["dist", "poni1", "poni2", "rot1",
"rot2_offset", "rot2_scale"],
dist_expr="dist",
poni1_expr="poni1",
poni2_expr="poni2",
rot1_expr="rot1",
rot2_expr="rot2_scale * pos + rot2_offset",
rot3_expr="0.0")
#Definition of the function reading the goniometer angle from the filename of the image.
def get_angle(basename):
"""Takes the basename (like del_65.0_0001p ) and returns the angle of the detector"""
return float(basename.split("_")[1])
print('filename', fimg.filename, "angle:",get_angle(fimg.filename))
In [5]:
#Definition of the detector, its mask, the calibrant
mask1 = fabio.open("deviation-mask.edf").data
mask2 = fabio.open("minimum-mask.edf").data
mask = numpy.logical_or(mask1, mask2)
pilatus = pyFAI.detector_factory("Pilatus100k")
pilatus.mask = mask
LaB6 = pyFAI.calibrant.CALIBRANT_FACTORY("LaB6")
wavelength = 7.7490120575e-11
LaB6.wavelength = wavelength
In [6]:
#Definition of the geometry refinement: the parameter order is the same as the param_names
param = {"dist":0.8,
"poni1":0.02,
"poni2":0.04,
"rot1":0,
"rot2_offset":0,
"rot2_scale": numpy.pi/180. # rot2 is in radians, while the motor position is in degrees
}
#Defines the bounds for some variables
bounds = {"dist": (0.79, 0.81),
"rot1": (-0.01, 0.01),
"rot2_offset": (-0.01, 0.01),
"rot2_scale": (numpy.pi/180., numpy.pi/180.) #strict bounds on the scale: we expect the gonio to be precise
}
gonioref = GoniometerRefinement(param, #initial guess
bounds=bounds,
pos_function=get_angle,
trans_function=goniotrans,
detector=pilatus, wavelength=wavelength)
print("Empty refinement object:")
print(gonioref)
In [7]:
#Let's populate the goniometer refinement object with all control point files:
ponis = glob.glob("*.poni")
ponis.sort()
for fn in ponis:
base = os.path.splitext(fn)[0]
fimg = fabio.open(base + ".cbf")
sg =gonioref.new_geometry(base, image=fimg.data, metadata=base, control_points=base+".npt",
geometry=fn, calibrant=LaB6)
print(base, "Angle:", sg.get_position())
print(sg.geometry_refinement)
print()
print("Filled refinement object:")
print(gonioref)
In [8]:
#Display all images with associated calibration:
for sg in gonioref.single_geometries.values():
jupyter.display(sg=sg)
In [9]:
# Initial refinement of the goniometer model with 5 dof
gonioref.refine2()
Out[9]:
In [10]:
# This function adds new images to the pool of data used for the refinement.
# A set of new control points are extractred and a refinement step is performed at each iteration
# The last image of the serie is displayed
def optimize_with_new_images(list_images, pts_per_deg=1):
sg = None
for fname in list_images:
print()
base = os.path.splitext(fname)[0]
fimg = fabio.open(fname)
if base in gonioref.single_geometries:
continue
print(base)
sg = gonioref.new_geometry(base, image=fimg.data, metadata=base,
calibrant=LaB6)
print(sg.extract_cp(pts_per_deg=pts_per_deg))
print("*"*50)
gonioref.refine2()
if sg:
sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param)
jupyter.display(sg=sg)
In [11]:
# Append all other images bewteen 30 and 40 degrees
images_30_40 = glob.glob("del_3?.?_0001p.cbf")
random.shuffle(images_30_40)
optimize_with_new_images(images_30_40, pts_per_deg=3)
In [12]:
# Append all other images
all_images = glob.glob("del_*_0001p.cbf")
random.shuffle(all_images)
optimize_with_new_images(all_images)
In [13]:
# Check the calibration of the first and the last image with rings
print("First & last rings")
print("Total number of images:", len(gonioref.single_geometries) )
fig = plt.figure()
for idx,lbl in enumerate(["del_10.0_0001p", "del_65.0_0001p"]):
sg = gonioref.single_geometries[lbl]
if sg.control_points.get_labels():
sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param)
jupyter.display(sg=sg, ax=fig.add_subplot(2, 1, idx+1))
In [14]:
# Final pass of refinement with all constrains removed, very fine refinement
gonioref.bounds = None
gonioref.refine2("slsqp", eps=1e-13, maxiter=10000, ftol=1e-12)
Out[14]:
In [15]:
#Create a MultiGeometry integrator from the refined geometry:
angles = []
images = []
for sg in gonioref.single_geometries.values():
angles.append(sg.get_position())
images.append(sg.image)
multigeo = gonioref.get_mg(angles)
multigeo.radial_range=(3, 68)
print(multigeo)
In [16]:
# Integrate the whole set of images in a single run:
res = multigeo.integrate1d(images, 10000)
fig, ax = subplots()
ax.plot(*res)
ax.set_xlabel(res.unit.label)
ax.set_ylabel("Intensity")
Out[16]:
In [17]:
# Save the goniometer configuration with 1 angle
gonioref.save("ROBL_v1.json")
In [18]:
#Can the refinement be improved by freeing another degree of freedom ? what about rot1 ?
goniotrans2 = GeometryTransformation(param_names = ["dist", "poni1", "poni2",
"rot1", "rot1_scale",
"rot2_offset", "rot2_scale"],
dist_expr="dist",
poni1_expr="poni1",
poni2_expr="poni2",
rot1_expr="rot1_scale * pos + rot1",
rot2_expr="rot2_scale * pos + rot2_offset",
rot3_expr="0.0")
param2 = (gonioref.nt_param(*gonioref.param))._asdict()
param2["rot1_scale"] = 0
gonioref2 = GoniometerRefinement(param2,
pos_function = get_angle,
trans_function=goniotrans2,
detector=pilatus,
wavelength=wavelength)
gonioref2.single_geometries = gonioref.single_geometries.copy()
print(gonioref2.chi2(), gonioref.chi2())
gonioref2.refine2()
gonioref2.save("ROBL_v2.json")
In [19]:
# Check the calibration of the first and the last image with rings
print("First & last rings")
fig = plt.figure()
for idx,lbl in enumerate(["del_10.0_0001p", "del_65.0_0001p"]):
sg = gonioref.single_geometries[lbl]
if sg.control_points.get_labels():
sg.geometry_refinement.set_param(gonioref2.get_ai(sg.get_position()).param)
jupyter.display(sg=sg, ax=fig.add_subplot(2, 1, idx+1))
In [20]:
#Create a MultiGeometry integrator from the refined geometry and display the integrated image:
multigeo2 = gonioref2.get_mg(angles)
multigeo2.radial_range=(3, 68)
print(multigeo2)
res2 = multigeo2.integrate1d(images, 10000)
#Display the 2 curves with a zoom
fig = figure(figsize=(10,5))
ax = fig.add_subplot(1,2,1)
ax.plot(*res2, label="rot1 & rot2 rotation")
ax.plot(*res, label="rot2 only rotation")
ax.set_xlabel(res.unit.label)
ax.set_ylabel("Intensity")
ax.set_title("Azimuthal integration of 121 images merged")
ax.legend()
ay = fig.add_subplot(1,2,2)
ay.plot(*res2, label="rot1 & rot2 rotation")
ay.plot(*res, label="rot2 only rotation")
ay.set_xlabel(res.unit.label)
ay.set_ylabel("Intensity")
ay.set_xlim(10.5,11)
ay.set_title("Zoom on first peak")
ay.legend()
Out[20]:
In [21]:
print("Total execution time: %.3fs"%(time.time() - start_time))
With the first model, the refinement was not perfect on the very low angles and indicates a miss-fit. Relaxing the constrains on rot1 allowed to spot a non (perfect) orthogonality between the goniometer axis and the incident beam. Releasing the distances is also possible, for example to cope with the sample not perfectly mounted on the center of the goniometer.
This notebook exposes the how to calibrate the goniometer for a small detector moving on a 2theta arm. Once calibrated, the geometry can be saved and restored and stays valid as long as the detector or the goniometer is not unmounted from the beam-line. This configuration can subsequently be used to integrate the data acquired with any sample, in minutes instead of hours. The resolution limit is now the pixel size. Fortunately, pixel detector with small pixel like the MaxiPix or the Lambda detector exists and offer higher resolution.