Calibration of a detector on a translation table
The aim of this document is to explain how to use pyFAI.goniometer for calibrating the position detector from the translation table encoders.
Those data have been acquired at ESRF-ID29 in summer 2013 on a Pilatus 6M using Ceria (CeO2) as calibrant. Seven images have been acquired with the detector moved between 15 cm and 45 cm from the sample position. A prior calibration has been performed using the MX-calibrate script from the pyFAI suite. The control points extracted during this initial calibration have been used as a starting point for this calibration.
The raw data files are available at: http://www.silx.org/pub/pyFAI/gonio/MX-ceria/
[1]:
# Initialization of the plotting library for use in the Jupyter notebook
%matplotlib inline
import numpy
from matplotlib.pyplot import subplots
[2]:
# Loading of a few libraries
import time
start_time =time.perf_counter()
import os
import fabio
import pyFAI
from pyFAI.goniometer import GeometryTransformation, GoniometerRefinement, Goniometer
from pyFAI.gui import jupyter
print("PyFAI version:", pyFAI.version)
PyFAI version: 0.21.0-dev1
[3]:
#Download all images
#Nota: comment-out to configure a proxy if you are behind a firewall
#os.environ["http_proxy"] = "http://proxy.company.com:3128"
from silx.resources import ExternalResources
downloader = ExternalResources("pyFAI", "http://www.silx.org/pub/pyFAI/testimages", "PYFAI_DATA")
all_files = downloader.getdir("MX_ceria.tar.bz2")
print(os.linesep.join(os.path.basename(i) for i in all_files))
MX_ceria
ceria_300_1_0001.npt
ceria_150_1_0001.cbf
ceria_200_1_0001.npt
ceria_250_1_0001.cbf
ceria_250_1_0001.npt
ceria_400_1_0001.cbf
ceria_400_1_0001.npt
ceria_450_1_0001.cbf
ceria_150_1_0001.poni
ceria_150_1_0001.npt
ceria_200_1_0001.poni
ceria_350_1_0001.poni
ceria_300_1_0001.cbf
ceria_450_1_0001.npt
ceria_250_1_0001.poni
ceria_350_1_0001.cbf
ceria_350_1_0001.npt
ceria_400_1_0001.poni
ceria_450_1_0001.poni
ceria_200_1_0001.cbf
ceria_300_1_0001.poni
[4]:
# Loading of the list of files, and display of the first one with its headers
image_files = [i for i in all_files if i.endswith(".cbf")]
image_files.sort()
print("List of images: " + ", ".join(image_files) + "." + os.linesep)
fimg = fabio.open(image_files[0])
print("Image headers:")
for key, value in fimg.header.items():
print("%s: %s"%(key,value))
jupyter.display(fimg.data, label=os.path.basename(fimg.filename))
List of images: /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_150_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_200_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_250_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_300_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_350_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_400_1_0001.cbf, /tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_450_1_0001.cbf.
Image headers:
_array_data.header_contents: # Detector: PILATUS 6M, S/N 60-0104, ESRF ID29
# 2013/Aug/29 17:26:59.699
# Pixel_size 172e-6 m x 172e-6 m
# Silicon sensor, thickness 0.000320 m
# Start_angle 0.000000 deg.
# Exposure_time 0.037000 s
# Exposure_period 0.040000 s
# Tau = 0 s
# Count_cutoff 1048500
# Threshold_setting 7612 eV
# N_excluded_pixels = 321
# Excluded_pixels: badpix_mask.tif
# Flat_field: (nil)
# Trim_directory: (nil)
# Wavelength 0.972386 A
# Detector_distance 0.150000 m
# Energy_range (0, 0) eV
# Detector_Voffset 0.0000 m
# Beam_xy (1230.90, 1254.09) pixels
# Flux 2.823146e+11 ph/s
# Transmission 20.1173
# Angle_increment 1.0000 deg.
# Detector_2theta 0.0000 deg.
# Polarization 0.99
# Alpha 0.0000 deg.
# Kappa 0.0020 deg.
# Phi 0.0000 deg.
# Chi 0.0000 deg.
# Oscillation_axis omega
# N_oscillations 1
# file_comments
Content-Type: application/octet-stream;
conversions: x-CBF_BYTE_OFFSET
Content-Transfer-Encoding: BINARY
X-Binary-Size: 6262451
X-Binary-ID: 0
X-Binary-Element-Type: signed 32-bit integer
X-Binary-Element-Byte-Order: LITTLE_ENDIAN
Content-MD5: BIfsFrKJBFklJn97/hjO/A==
X-Binary-Number-of-Elements: 6224001
X-Binary-Size-Fastest-Dimension: 2463
X-Binary-Size-Second-Dimension: 2527
X-Binary-Size-Padding: 128
[4]:
<AxesSubplot:title={'center':'ceria_150_1_0001.cbf'}>
[5]:
# Definition of the geometry translation function:
geotrans = GeometryTransformation(param_names = ["dist_offset", "dist_scale",
"poni1", "poni2", "rot1","rot2"],
dist_expr="pos * dist_scale + dist_offset",
poni1_expr="poni1",
poni2_expr="poni2",
rot1_expr="rot1",
rot2_expr="rot2",
rot3_expr="0.0")
# Definition of the function reading the detector position from the header of the image.
def get_distance(header):
"""Takes the header of the CBF-file and returns the distance of the detector"""
dist = 0
for line in header.get("_array_data.header_contents","").split("\n"):
words = line.split()
if words[1] == "Detector_distance":
dist = float(words[2])
break
return dist
print("Distance:",get_distance(fimg.header))
Distance: 0.15
[6]:
# Definition of the detector, the calibrant and extraction of the wavelength used from the headers
pilatus = pyFAI.detector_factory("Pilatus6M")
CeO2 = pyFAI.calibrant.get_calibrant("CeO2")
for line in fimg.header.get("_array_data.header_contents","").split("\n"):
words = line.split()
if words[1] == "Wavelength":
wavelength = float(words[2])*1e-10
break
print("Wavelength:", wavelength)
CeO2.wavelength = wavelength
print(CeO2)
Wavelength: 9.72386e-11
CeO2 Calibrant with 42 reflections at wavelength 9.72386e-11
[7]:
# Definition of the geometry refinement: the parameter order is the same as the param_names
param = {"dist_offset":0,
"dist_scale":1,
"poni1":0.2,
"poni2":0.2,
"rot1":0,
"rot2":0}
gonioref = GoniometerRefinement(param, #initial guess
pos_function=get_distance,
trans_function=geotrans,
detector=pilatus,
wavelength=wavelength)
print("Empty refinement object:")
print(gonioref)
Empty refinement object:
GoniometerRefinement with 0 geometries labeled: .
[8]:
# Let's populate the goniometer refinement object with all control point files:
ponis = [i for i in all_files if i.endswith(".poni")]
ponis.sort()
print(ponis)
for fn in ponis:
base = os.path.splitext(fn)[0]
fimg = fabio.open(base + ".cbf")
gonioref.new_geometry(os.path.basename(base), image=fimg.data, metadata=fimg.header, control_points=base+".npt",
geometry=fn, calibrant=CeO2)
print("Filled refinement object:")
print(gonioref)
print(os.linesep+"\tLabel \t Distance")
for k, v in gonioref.single_geometries.items():
print(k,v.get_position())
['/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_150_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_200_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_250_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_300_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_350_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_400_1_0001.poni', '/tmp/pyFAI_testdata_kieffer/MX_ceria.tar.bz2__content/MX_ceria/ceria_450_1_0001.poni']
Filled refinement object:
GoniometerRefinement with 7 geometries labeled: ceria_150_1_0001, ceria_200_1_0001, ceria_250_1_0001, ceria_300_1_0001, ceria_350_1_0001, ceria_400_1_0001, ceria_450_1_0001.
Label Distance
ceria_150_1_0001 0.15
ceria_200_1_0001 0.2
ceria_250_1_0001 0.25
ceria_300_1_0001 0.3
ceria_350_1_0001 0.35
ceria_400_1_0001 0.4
ceria_450_1_0001 0.45
[9]:
# Display all images with associated calibration:
fig, ax = subplots(4,2, figsize=(15, 20))
for idx, sg in enumerate(gonioref.single_geometries.values()):
jupyter.display(sg=sg, ax=ax[idx//2,idx%2])
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
[10]:
# Initial refinement of the translation table model
gonioref.refine3()
Free parameters: ['dist_offset', 'dist_scale', 'poni1', 'poni2', 'rot1', 'rot2']
Fixed: {}
fun: 5.11938033229523e-07
jac: array([ 1.15391254e-06, 3.42231246e-07, 8.12186798e-08, 5.02960944e-08,
-7.17341919e-08, -1.11380359e-07])
message: 'Optimization terminated successfully'
nfev: 130
nit: 18
njev: 18
status: 0
success: True
x: array([-0.00118793, 1.00190471, 0.21548514, 0.21309905, 0.00661241,
0.00280329])
Constrained Least square 0.0016696847686478021 --> 5.11938033229523e-07
maxdelta on poni1: 0.2 --> 0.2154851363228725
[10]:
5.11938033229523e-07
[11]:
# Save the result of the fitting to a file and display the content of the JSON file:
gonioref.save("ID29.json")
with open("ID29.json") as fd:
print(fd.read())
{
"content": "Goniometer calibration v2",
"detector": "Pilatus 6M",
"detector_config": {},
"wavelength": 9.72386e-11,
"param": [
-0.001187934613201132,
1.0019047068951217,
0.2154851363228725,
0.2130990512962314,
0.006612408551460147,
0.0028032884874256087
],
"param_names": [
"dist_offset",
"dist_scale",
"poni1",
"poni2",
"rot1",
"rot2"
],
"pos_names": [
"pos"
],
"trans_function": {
"content": "GeometryTransformation",
"param_names": [
"dist_offset",
"dist_scale",
"poni1",
"poni2",
"rot1",
"rot2"
],
"pos_names": [
"pos"
],
"dist_expr": "pos * dist_scale + dist_offset",
"poni1_expr": "poni1",
"poni2_expr": "poni2",
"rot1_expr": "rot1",
"rot2_expr": "rot2",
"rot3_expr": "0.0",
"constants": {
"pi": 3.141592653589793
}
}
}
[12]:
# Restore the translation table setting from the file
transtable = Goniometer.sload("ID29.json")
print("Translation table: \n",transtable)
Translation table:
Goniometer with param GonioParam(dist_offset=-0.001187934613201132, dist_scale=1.0019047068951217, poni1=0.2154851363228725, poni2=0.2130990512962314, rot1=0.006612408551460147, rot2=0.0028032884874256087)
with Detector Pilatus 6M PixelSize= 1.720e-04, 1.720e-04 m
[13]:
# Create a multi-geometry object for all images in this set:
distances = [get_distance(fabio.open(fn).header) for fn in image_files]
print("Distances: ", distances)
multigeo = transtable.get_mg(distances)
multigeo.radial_range=(0, 65)
print(multigeo)
Distances: [0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45]
MultiGeometry integrator with 7 geometries on (0, 65) radial range (2th_deg) and (-180, 180) azimuthal range (deg)
[14]:
# Integrate the set of images in a single run:
res = multigeo.integrate1d([fabio.open(fn).data for fn in image_files], 10000)
# Display the result using matplotlib
fig, ax = subplots()
ax.plot(*res)
ax.set_xlabel(res.unit.label)
ax.set_ylabel("Intensity")
ax.set_xlim(17, 22)
ax.set_title("Zoom on the two first rings")
pass
Accoring to the provious image, peaks look double which indicates a bad modeling of the setup or a bad fitting. As the fitting ended successfully, the bug is likely in the model: let’s allow the PONI to move with the distance
[15]:
# Let's refine poni1 and poni2 also as function of the distance:
geotrans2 = GeometryTransformation(param_names = ["dist_offset", "dist_scale",
"poni1_offset", "poni1_scale",
"poni2_offset", "poni2_scale",
"rot1","rot2"],
dist_expr="pos * dist_scale + dist_offset",
poni1_expr="pos * poni1_scale + poni1_offset",
poni2_expr="pos * poni2_scale + poni2_offset",
rot1_expr="rot1",
rot2_expr="rot2",
rot3_expr="0.0")
#initial guess from former parameter set
param2 = (gonioref.nt_param(*gonioref.param))._asdict()
param2["poni1_offset"] = 0
param2["poni2_offset"] = 0
param2["poni1_scale"] = 1
param2["poni2_scale"] = 1
gonioref2 = GoniometerRefinement(param2,
pos_function = get_distance,
trans_function=geotrans2,
detector=pilatus,
wavelength=wavelength)
gonioref2.single_geometries = gonioref.single_geometries.copy()
print(gonioref2)
GoniometerRefinement with 7 geometries labeled: ceria_150_1_0001, ceria_200_1_0001, ceria_250_1_0001, ceria_300_1_0001, ceria_350_1_0001, ceria_400_1_0001, ceria_450_1_0001.
[16]:
# Refinement of the second model with all distances free
gonioref2.refine2()
Cost function before refinement: 0.04364487420511431
[-0.00118793 1.00190471 0. 1. 0. 1.
0.00661241 0.00280329]
fun: 1.6219985734022523e-07
jac: array([ 5.32202590e-07, 8.74261445e-08, 1.65539717e-07, 1.80200939e-08,
-2.92968990e-07, -7.38344390e-08, 8.45039771e-08, -1.73615469e-08])
message: 'Optimization terminated successfully'
nfev: 310
nit: 34
njev: 34
status: 0
success: True
x: array([-0.00118686, 1.00184287, 0.21574533, -0.00429667, 0.21300993,
0.00138094, 0.00735187, 0.00492121])
Cost function after refinement: 1.6219985734022523e-07
GonioParam(dist_offset=-0.0011868649118234032, dist_scale=1.0018428737070253, poni1_offset=0.21574533062901727, poni1_scale=-0.004296673901268426, poni2_offset=0.21300993413666963, poni2_scale=0.0013809412146645605, rot1=0.007351867557285051, rot2=0.004921206654305685)
maxdelta on: poni1_scale (3) 1 --> -0.004296673901268426
[16]:
array([-0.00118686, 1.00184287, 0.21574533, -0.00429667, 0.21300993,
0.00138094, 0.00735187, 0.00492121])
[17]:
# Integration of all images with the second model
multigeo2 = gonioref2.get_mg(distances)
multigeo2.radial_range=(0, 65)
print(multigeo2)
res2 = multigeo2.integrate1d([fabio.open(fn).data for fn in image_files], 10000)
# Display the result, zooming on the two first rings
fig, ax = subplots(figsize=(15,8))
ax.plot(*res)
ax.plot(*res, label="only distance free")
ax.plot(*res2, label="distance and PONI free")
ax.set_ylabel("Intensity")
ax.set_xlim(17, 22)
ax.set_title("Zoom on the two first rings")
ax.set_xlabel(res2.unit.label)
ax.legend()
pass
MultiGeometry integrator with 7 geometries on (0, 65) radial range (2th_deg) and (-180, 180) azimuthal range (deg)
[18]:
# Re-extract many more control points from images for a better fit
fig, ax = subplots(4,2, figsize=(15, 20))
for idx, sg in enumerate(gonioref2.single_geometries.values()):
sg.extract_cp(pts_per_deg=3)
jupyter.display(sg=sg, ax=ax[idx//2, idx%2])
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
WARNING:matplotlib.legend:No handles with labels found to put in legend.
[19]:
# Refine again the model
gonioref2.refine2()
# Build the MultiGeometry integrator object
multigeo3 = gonioref2.get_mg(distances)
multigeo3.radial_range=(0, 65)
print(multigeo3)
# Perform the azimuthal integration
res3 = multigeo3.integrate1d([fabio.open(fn).data for fn in image_files], 10000)
# Display the result
fig, ax = subplots(figsize=(15,8))
ax.plot(*res, label="only distance free")
ax.plot(*res2, label="distance and PONI free")
ax.plot(*res2, linestyle="--", label="distance and PONI free, more points")
ax.set_xlabel(res2.unit.label)
ax.set_xlim(17, 22)
ax.legend()
pass
Cost function before refinement: 5.378166443047018e-08
[-0.00118686 1.00184287 0.21574533 -0.00429667 0.21300993 0.00138094
0.00735187 0.00492121]
fun: 4.993441406328794e-08
jac: array([ 1.02208007e-07, -3.59430263e-10, 1.31174470e-07, 3.72928654e-07,
9.23606140e-07, 3.12824059e-07, 1.70757263e-07, -1.73712272e-07])
message: 'Optimization terminated successfully'
nfev: 103
nit: 11
njev: 11
status: 0
success: True
x: array([-0.00120984, 1.00183186, 0.21574401, -0.00410805, 0.21301368,
0.00124719, 0.00726244, 0.00480959])
Cost function after refinement: 4.993441406328794e-08
GonioParam(dist_offset=-0.0012098389850379633, dist_scale=1.0018318622865574, poni1_offset=0.2157440074986014, poni1_scale=-0.0041080541102543324, poni2_offset=0.21301368148936745, poni2_scale=0.0012471921732800323, rot1=0.00726243612929008, rot2=0.00480959379409951)
maxdelta on: poni1_scale (3) -0.004296673901268426 --> -0.0041080541102543324
MultiGeometry integrator with 7 geometries on (0, 65) radial range (2th_deg) and (-180, 180) azimuthal range (deg)
This re-extraction of control point did not help to get a sharper diffraction profile. This step was not needed.
Conclusion
This notebook exposes the how to calibrate a translation table for a moving detector. It allows to: * Check the proper alignement of the table regarding the actual beam * Check the encoder’s precision (usually good) and offsets (arbitrary) * Perform azimuthal integration to retrieve powder diffraction patterns at any position of the detector.
[20]:
print(f"Total execution time: {time.perf_counter() - start_time:.3f} s")
Total execution time: 100.145 s