diff options
Diffstat (limited to 'contrib/python')
-rw-r--r-- | contrib/python/20170063.py | 254 | ||||
-rw-r--r-- | contrib/python/cirpad/blender_pyfai.py | 358 | ||||
-rw-r--r-- | contrib/python/common.py | 68 | ||||
-rw-r--r-- | contrib/python/mars.py | 368 | ||||
-rw-r--r-- | contrib/python/mythen.py | 306 | ||||
-rw-r--r-- | contrib/python/swing/pinhole1.smv | bin | 0 -> 33554944 bytes | |||
-rw-r--r-- | contrib/python/swing/plot.py | 12 |
7 files changed, 1366 insertions, 0 deletions
diff --git a/contrib/python/20170063.py b/contrib/python/20170063.py new file mode 100644 index 0000000..e5a5357 --- /dev/null +++ b/contrib/python/20170063.py @@ -0,0 +1,254 @@ +#!/usr/bin/env python3 +""" coding: utf-8 +""" + +from typing import Iterator, List, NamedTuple, Text, Tuple, Union + +import os + +from functools import partial +from math import pi, radians + +import h5py +import pylab +import pyFAI + +from fabio.edfimage import edfimage +from numpy import ndarray + +from pyFAI.goniometer import GeometryTransformation, GoniometerRefinement +from pyFAI.gui import jupyter + +from common import * + +ROOT = "/home/experiences/diffabs/reguer/20170063/" +PUBLISHED = os.path.join(ROOT, "published-data") +CALIB = os.path.join(PUBLISHED, "xrd", "calibrationCeO2") + +MetaDataSource = NamedTuple("MetaDataSource", [("image", H5Path), + ("delta", H5Path)]) + +MetaData = NamedTuple("MetaData", [("image", ndarray), + ("delta", float)]) + +_Diffabs = NamedTuple("_Diffabs", [("filename", Text), + ("basedir", Text), + ("metasources", MetaDataSource), + ("idxs", List[int]), + ("calibrant", Text), + ("detector", Text), + ("wavelength", float)]) + + +class Diffabs(_Diffabs): + def __len__(self) -> int: + with h5py.File(self.filename, mode='r') as f: + return get_shape(f, self.metasources.image)[0] + + def __item(self, f: h5py.File, index: int) -> MetaData: + return MetaData(get_item_at_index(f, + self.metasources.image, index), + get_item_at_index(f, + self.metasources.delta, index)) + + def __item__(self, index: int) -> MetaData: + with h5py.File(self.filename, mode='r') as f: + return self.__item(f, index) + + def frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in self.idxs: + yield self.__item(f, index) + + def all_frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in range(len(self)): + yield self.__item(f, index) + + +def save_as_edf(calib: Diffabs, + basedir: Text) -> None: + """Save the multi calib images into edf files in order to do the first + calibration""" + for idx, metadata in zip(calib.idxs, calib.frames()): + base = os.path.splitext(os.path.basename(calib.filename))[0] + output = os.path.join(basedir, base + '_%d.edf' % (idx,)) + edfimage(metadata.spectrum).write(output) + + +def get_wavelength(multicalib: Diffabs) -> float: + """Return the wavelength""" + return multicalib.wavelength + + +def get_calibrant(multicalib: Diffabs) -> pyFAI.calibrant.Calibrant: + """Return the calibrant with the right wavelength""" + calibrant = pyFAI.calibrant.get_calibrant(multicalib.calibrant) + calibrant.wavelength = get_wavelength(multicalib) + return calibrant + + +def get_detector(multicalib: Diffabs) -> pyFAI.Detector: + return pyFAI.detector_factory(multicalib.detector) + + +def optimize_with_new_images(multicalib: Diffabs, + gonioref: GoniometerRefinement, + calibrant: pyFAI.calibrant.Calibrant, + pts_per_deg: float=1) -> None: + """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 + + """ + sg = None + for idx, metadata in enumerate(multicalib.all_frames()): + print() + label = "CeO2_img%d_poni" % (idx,) + if label in gonioref.single_geometries: + continue + print(label) + + sg = gonioref.new_geometry(label, image=metadata.image, + metadata=metadata, calibrant=calibrant) + print(sg.extract_cp(pts_per_deg=pts_per_deg))# , "watershed")) + print("*"*50) + gonioref.refine2() + if sg: + sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param) # noqa + jupyter.display(sg=sg) + + +def calibration(json: str, multicalib:Diffabs) -> None: + """Do a calibration with a bunch of images""" + + # Definition of the geometry refinement: the parameter order is + # the same as the param_names + calibrant = get_calibrant(multicalib) + detector = get_detector(multicalib) + + distance = 0.39 + poni1 = -0.01 + poni2 = 0.033 + rot1 = 0 + rot2_scale = -1 * pi / 180. + rot2_offset = radians(16) + rot3 = radians(-90) + + parameters = [Parameter("dist", distance, (distance, distance)), + Parameter("poni1", poni1, (poni1-1, poni1+1)), + Parameter("poni2", poni2, (poni2-1, poni2+1)), + Parameter("rot1", rot1, (rot1, rot1)), + Parameter("rot2_scale", rot2_scale, (rot2_scale, rot2_scale)), + Parameter("rot2_offset", rot2_offset, (rot2_offset-0.1, rot2_offset+0.1)), + Parameter("rot3", rot3, (rot3, rot3))] + + params = {p.name: p.value for p in parameters} + bounds = {p.name: p.bounds for p in parameters} + param_names = [p.name for p in parameters] + + # Let's refine poni1 and poni2 also as function of the distance: + + trans_function = GeometryTransformation(param_names=param_names, + pos_names=["delta"], + dist_expr="dist", + poni1_expr="poni1", + poni2_expr="poni2", + rot1_expr="rot1", + rot2_expr="delta * rot2_scale + rot2_offset", + rot3_expr="rot3") + + def pos_function(metadata: MetaData) -> Tuple[float]: + """Definition of the function reading the detector position from the + header of the image.""" + return metadata.delta, + + gonioref = GoniometerRefinement(params, # initial guess + bounds=bounds, + pos_function=pos_function, + trans_function=trans_function, + detector=detector, + wavelength=multicalib.wavelength) + + print("Empty refinement object:") + print(gonioref) + + # Let's populate the goniometer refinement object with the known poni + for idx, metadata in zip(multicalib.idxs, multicalib.frames()): + label = "CeO2_img%d_poni" % (idx,) + control_points = os.path.join(multicalib.basedir, "img%d_pts.npt" % (idx,)) + ai = pyFAI.load(os.path.join(multicalib.basedir, label)) # noqa + print(ai) + + gonioref.new_geometry(label, metadata.image, metadata, + control_points, calibrant, ai) + + print("Filled refinement object:") + print(gonioref) + print(os.linesep + "\tlabel \t delta") + for k, v in gonioref.single_geometries.items(): + print(k, v.get_position()) + + for g in gonioref.single_geometries.values(): + ai = gonioref.get_ai(g.get_position()) + print(ai) + + for sg in gonioref.single_geometries.values(): + jupyter.display(sg=sg) + + gonioref.refine2() + + for multi in [multicalib]: + optimize_with_new_images(multi, gonioref, calibrant) + + for idx, sg in enumerate(gonioref.single_geometries.values()): + sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param) + jupyter.display(sg=sg) + + gonioref.save(json) + pylab.show() + + +def integrate(json: Text) -> None: + """Integrate a file with a json calibration file""" + gonio = pyFAI.goniometer.Goniometer.sload(json) + filename = os.path.join(ROOT, "tth2C-16_11_15-55-22_008.nxs") + wavelength = 0.67186e-10 + multicalib = Diffabs(filename, + MetaDataSource(H5PathWithAttribute("interpretation", b"spectrum"), + H5PathContains("scan_data/actuator_1_1")), + [], 0, "LaB6", "mythen", wavelength) + + # print(len(multicalib)) + # metadata = multicalib.__item__(21) + + # ai = gonio.get_ai(metadata.tth) + # res = ai.integrate1d(metadata.spectrum, 1280, unit="2th_deg") + # jupyter.plot1d(res) + + images = [] + positions = [] + for metadata in multicalib.all_frames(): + images.append(metadata.spectrum) + positions.append((metadata.tth,)) + mai = gonio.get_mg(positions) + res = mai.integrate1d(images, 30000) + jupyter.plot1d(res) + pylab.show() + + +def main() -> None: + json = os.path.join(CALIB, "diffabs.json") + filename = os.path.join(ROOT, "2017", "Run5", "2017-11-15", "scan_81.nxs") + wavelength = 6.81198566418e-11 + calib = Diffabs(filename, + CALIB, + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), + H5PathContains("scan_data/actuator_1_1")), + [3, 4], "CeO2", "imxpad_s140", wavelength) + + calibration(json, calib) + +if __name__ == "__main__": + main() diff --git a/contrib/python/cirpad/blender_pyfai.py b/contrib/python/cirpad/blender_pyfai.py new file mode 100644 index 0000000..b2d1faa --- /dev/null +++ b/contrib/python/cirpad/blender_pyfai.py @@ -0,0 +1,358 @@ +import bpy +import numpy +import pyFAI +import math +import functools +import random +import bmesh + + +SCALE = 100 + +# Definition du détecteur Cirpad + +class Cirpad(pyFAI.detectors.Detector): # (Detector): + MAX_SHAPE = (11200, 120) + IS_FLAT = False + IS_CONTIGUOUS = False + force_pixel = True + uniform_pixel = False + aliases = ["XCirpad"] + MEDIUM_MODULE_SIZE = (560, 120) + MODULE_SIZE = (80, 120) # number of pixels per module (y, x) + PIXEL_SIZE = (130e-6, 130e-6) + DIFFERENT_PIXEL_SIZE = 2.5 + rot = [0,0,-6.74] + + def _calc_pixels_size(self,length, module_size, pixel_size): + size = numpy.ones(length) + n = (length // module_size) + for i in range(1, n): + size[i * module_size - 1] = self.DIFFERENT_PIXEL_SIZE + size[i * module_size] = self.DIFFERENT_PIXEL_SIZE + return pixel_size * size + + def passage(self,corners,rot): #,u) # A définir plus précisément comment on obtient rot et u + shape = corners.shape + origine = corners[0][0][0,:] + nmd = rotation(corners,rot) + u = corners[shape[0]-1][0][1,:] - corners[0][0][0,:] + u = u / numpy.linalg.norm(u) + s = rotation(u,rot) + s = s / numpy.linalg.norm(s) + v = numpy.array([-u[1],u[0],u[2]]) + r = origine - nmd[0][0][0,:] + w = (0.1e-3+0.24e-3+75.14e-3)*u + (0.8e-3)*v + (0.55e-3)*s + r + nmd = translation(nmd,w) + return(nmd) + + def get_pixel_corners(self): + pixel_size1 = self._calc_pixels_size(self.MEDIUM_MODULE_SIZE[0], self.MODULE_SIZE[0], self.PIXEL_SIZE[0]) + pixel_size2 = (numpy.ones(self.MEDIUM_MODULE_SIZE[1]) * self.PIXEL_SIZE[1]).astype(numpy.float32) + # half pixel offset + pixel_center1 = pixel_size1 / 2.0 # half pixel offset + pixel_center2 = pixel_size2 / 2.0 + # size of all preceeding pixels + pixel_center1[1:] += numpy.cumsum(pixel_size1[:-1]) + pixel_center2[1:] += numpy.cumsum(pixel_size2[:-1]) + + pixel_center1.shape = -1, 1 + pixel_center1.strides = pixel_center1.strides[0], 0 + + pixel_center2.shape = 1, -1 + pixel_center2.strides = 0, pixel_center2.strides[1] + + pixel_size1.shape = -1, 1 + pixel_size1.strides = pixel_size1.strides[0], 0 + + pixel_size2.shape = 1, -1 + pixel_size2.strides = 0, pixel_size2.strides[1] + + # On calcule la position du premier module + corners = numpy.zeros((self.MEDIUM_MODULE_SIZE[0], self.MEDIUM_MODULE_SIZE[1], 4, 3), dtype=numpy.float32) + corners[:, :, 0, 1] = pixel_center1 - pixel_size1 / 2.0 + corners[:, :, 0, 2] = pixel_center2 - pixel_size2 / 2.0 + corners[:, :, 1, 1] = pixel_center1 + pixel_size1 / 2.0 + corners[:, :, 1, 2] = pixel_center2 - pixel_size2 / 2.0 + corners[:, :, 2, 1] = pixel_center1 + pixel_size1 / 2.0 + corners[:, :, 2, 2] = pixel_center2 + pixel_size2 / 2.0 + corners[:, :, 3, 1] = pixel_center1 - pixel_size1 / 2.0 + corners[:, :, 3, 2] = pixel_center2 + pixel_size2 / 2.0 + + n_corners = corners + # Puis on calcule les coins pour les 19 modules restant + for i in range(1, 20): + n_corners = self.passage(n_corners, self.rot) + corners = numpy.concatenate((corners, n_corners), axis=0) # A voir la disposition finale souhaité + return corners + + # Pas fait encore + def calc_cartesian_positions(self,d1 = None, d2 = None, center = True, use_cython = True): + if (d1 is None) or d2 is None: + d1 = pyFAI.utils.expand2d(numpy.arange(self.MAX_SHAPE[0]).astype(numpy.float32), self.MAX_SHAPE[1], False) + d2 = pyFAI.utils.expand2d(numpy.arange(self.MAX_SHAPE[1]).astype(numpy.float32), self.MAX_SHAPE[0], True) + corners = self.get_pixel_corners() + if center: + # avoid += It modifies in place and segfaults + d1 = d1 + 0.5 + d2 = d2 + 0.5 + if False and use_cython: + p1, p2, p3 = bilinear.calc_cartesian_positions(d1.ravel(), d2.ravel(), corners, is_flat=False) + p1.shape = d1.shape + p2.shape = d2.shape + p3.shape = d2.shape + else: # A verifier + i1 = d1.astype(int).clip(0, corners.shape[0] - 1) + i2 = d2.astype(int).clip(0, corners.shape[1] - 1) + delta1 = d1 - i1 + delta2 = d2 - i2 + pixels = corners[i1, i2] + if pixels.ndim == 3: + A0 = pixels[:, 0, 0] + A1 = pixels[:, 0, 1] + A2 = pixels[:, 0, 2] + B0 = pixels[:, 1, 0] + B1 = pixels[:, 1, 1] + B2 = pixels[:, 1, 2] + C0 = pixels[:, 2, 0] + C1 = pixels[:, 2, 1] + C2 = pixels[:, 2, 2] + D0 = pixels[:, 3, 0] + D1 = pixels[:, 3, 1] + D2 = pixels[:, 3, 2] + else: + A0 = pixels[:, :, 0, 0] + A1 = pixels[:, :, 0, 1] + A2 = pixels[:, :, 0, 2] + B0 = pixels[:, :, 1, 0] + B1 = pixels[:, :, 1, 1] + B2 = pixels[:, :, 1, 2] + C0 = pixels[:, :, 2, 0] + C1 = pixels[:, :, 2, 1] + C2 = pixels[:, :, 2, 2] + D0 = pixels[:, :, 3, 0] + D1 = pixels[:, :, 3, 1] + D2 = pixels[:, :, 3, 2] + + # points A and D are on the same dim1 (Y), they differ in dim2 (X) + # points B and C are on the same dim1 (Y), they differ in dim2 (X) + # points A and B are on the same dim2 (X), they differ in dim1 (Y) + # points C and D are on the same dim2 (X), they differ in dim1 ( + p1 = A1 * (1.0 - delta1) * (1.0 - delta2) \ + + B1 * delta1 * (1.0 - delta2) \ + + C1 * delta1 * delta2 \ + + D1 * (1.0 - delta1) * delta2 + p2 = A2 * (1.0 - delta1) * (1.0 - delta2) \ + + B2 * delta1 * (1.0 - delta2) \ + + C2 * delta1 * delta2 \ + + D2 * (1.0 - delta1) * delta2 + p3 = A0 * (1.0 - delta1) * (1.0 - delta2) \ + + B0 * delta1 * (1.0 - delta2) \ + + C0 * delta1 * delta2 \ + + D0 * (1.0 - delta1) * delta2 + # To ensure numerical consitency with cython procedure. + p1 = p1.astype(numpy.float32) + p2 = p2.astype(numpy.float32) + p3 = p3.astype(numpy.float32) + return p1, p2, p3 + +# Fonction utilisée pour définir le Cirpad mais je ne sais pas si je les met dans la classe ou pas car elles sont banals + +def M(theta, u): + """ + :param theta: the axis value in radian + :type theta: float + :param u: the axis vector [x, y, z] + :type u: [float, float, float] + :return: the rotation matrix + :rtype: numpy.ndarray (3, 3) + """ + c = math.cos(theta) + one_minus_c = 1 - c + s = math.sin(theta) + return [[c + u[0]**2 * one_minus_c, + u[0] * u[1] * one_minus_c - u[2] * s, + u[0] * u[2] * one_minus_c + u[1] * s], + [u[0] * u[1] * one_minus_c + u[2] * s, + c + u[1]**2 * one_minus_c, + u[1] * u[2] * one_minus_c - u[0] * s], + [u[0] * u[2] * one_minus_c - u[1] * s, + u[1] * u[2] * one_minus_c + u[0] * s, + c + u[2]**2 * one_minus_c]] + + +def rotation(md,rot): + shape = md.shape + axe = numpy.array([[1,0,0],[0,1,0],[0,0,1]]) # A voir si ce n'est pas une entrée + P = functools.reduce(numpy.dot, [M(numpy.radians(rot[i]),axe[i]) for i in range(len(rot))]) + try : + nmd = numpy.transpose(numpy.reshape(numpy.tensordot(P,numpy.reshape(numpy.transpose(md),(3,shape[0]*shape[1]*4)), axes=1),(3,4,shape[1],shape[0]))) + except IndexError : + nmd = numpy.transpose(numpy.tensordot(P,numpy.transpose(md), axes=1)) + return(nmd) + + +def translation(md, u): + nmd = md + u + return(nmd) + +# Change coordonée pour passer de celle de pyFAI à celle de l'affichage. + +def inter2(array): + ar = numpy.empty(numpy.shape(array)) + ar[0] = array[2] + ar[1] = array[1] + ar[2] = array[0] + return(ar) + +# Affiche points coins pixels module +def print_corners(corners, name): + shape = corners.shape + vertices = [] + vertices = numpy.reshape(corners, (shape[0] * shape [1] * shape[2], 3)) + contour_mesh = bpy.data.meshes.new("contour_mesh") + contour_mesh .from_pydata(vertices, [], []) + contour_mesh .update() + contour_objet = bpy.data.objects.new(name, contour_mesh ) + scene = bpy.context.scene + contour_objet.active_material = mat1 + scene.objects.link(contour_objet) + + bpy.context.scene.objects.active = bpy.context.scene.objects[name] + bpy.ops.object.mode_set(mode = 'EDIT') + bpy.ops.mesh.remove_doubles() + bpy.ops.object.mode_set(mode = 'OBJECT') + + + +# Materiau classique utilisé +mat1 = bpy.data.materials.new('Mat1') +mat1.diffuse_color = (1,1,1) + +mat2 = bpy.data.materials.new('Mat2') +mat2.diffuse_color = (0,1,0) + + +# Affiche points centres pixels +def print_centers(centers, name): + shape = centers[0].shape + m = numpy.zeros((3, shape[0], shape[1])) + for i in range(3): + if centers[i] is not None: + m[i] = centers[i] * SCALE + + m = inter2(m) + + # Creation des vertices + vertices = numpy.reshape(numpy.transpose(m), (shape [0] * shape[1], 3)) + + # Creation du mesh à partir des vertices + contour_mesh = bpy.data.meshes.new("contour_mesh") + contour_mesh.from_pydata(vertices, [], []) + contour_mesh.update() + + # Creation d'un object avec ce mesh + contour_objet = bpy.data.objects.new(name, contour_mesh) + + # si pas de materiaux, le menu n'apparait pas ??? + contour_objet.active_material = mat2 + + # ajout de l'objet à la scène pour qu'il devienne visible + scene = bpy.context.scene + scene.objects.link(contour_objet) + + +# Affiche le detecteur avec ses couleurs, pas tourner dans le bon sens pour l'instant +def print_pixel(corners, image): + + shape = corners.shape + + # Creer des vertices du detecteur + vertices = numpy.reshape(corners * SCALE, (shape[0] * shape [1] * shape[2], 3)) + + # creation de la list contenant les indices des vertices formant une face (un pixel) + faces = [[4*i,4*i+1,4*i+2,4*i+3] for i in range(shape[0]*shape[1])] + + # creation du mesh avec ces vertices et ces faces. + contour_mesh = bpy.data.meshes.new("contour_mesh") + contour_mesh.from_pydata(vertices, [], faces) + contour_mesh.update() + + # creation de l'object a partir du mesh + contour_objet = bpy.data.objects.new('detector', contour_mesh ) + + contour_objet.active_material = mat1 + + # ajout de l'object dans la scene. + scene = bpy.context.scene + scene.objects.link(contour_objet) + + # Ajoute les couleurs + # on cree un groupe de vertex par couleurs de la LUT + # ici la lut est simple 0-255 (sur le rouge) + max = int(image.max()) + 1 + MyVertices = [[] for i in range(max)] + + # On partour toute l'image et pour chaque intensité, on rajoute les index des + # vertex de la list vertices formant une face dans le bon groupe de vertex. + # TODO improve speed :) + width = shape[1] + + # TODO voir comment utiliser plutot les faces pour ces couleurs. + for i in range(shape[0]): + i_width = i * width + for j in range(shape[1]): + i_width_j = 4 * (i_width + j) + l = range(i_width_j, i_width_j + 4) + MyVertices[int(image[i,j])].extend(l) + bpy.context.scene.objects.active = bpy.context.scene.objects['detector'] + + # Met les couleurs par groupe + for k in range(max): + + name = "MyVertexGroup"+str(k) + + #Create a Vertex Group + Group1 = bpy.context.object.vertex_groups.new(name) + + #Add the vertices to the vertex group + #with Weight = 1.0 The Weight isn't + #relevant in this case + Group1.add(MyVertices[k], 1.0, 'ADD') + + #Select the vertex group + bpy.ops.object.vertex_group_set_active(group= name) + + bpy.ops.object.material_slot_add() + bpy.context.object.material_slots[bpy.context.object.material_slots.__len__() - 1].material = MakeMaterial_1(k,max) + bpy.ops.object.mode_set(mode = 'EDIT') + #Deselect all the vertices + bpy.ops.mesh.select_all(action='DESELECT') + #Select the vertices of the vertex group + bpy.ops.object.vertex_group_select() + #Assign the material on the selected vertices + bpy.ops.object.material_slot_assign() + bpy.ops.object.mode_set(mode = 'OBJECT') + + bpy.ops.object.mode_set(mode = 'EDIT') + bpy.ops.mesh.remove_doubles() + bpy.ops.object.mode_set(mode = 'OBJECT') + + +# Fonction servant a créé un materiau par couleur, dépend de l'intensité +def MakeMaterial_1(k, max): + mat = bpy.data.materials.new("Mat2") + mat.diffuse_shader = 'MINNAERT' + mat.diffuse_color = (1, round(1-k/max,4), round(1-k/max,4)) + mat.darkness = 0.8 + return mat + +detector = pyFAI.detectors.ImXPadS70() +corners = detector.get_pixel_corners() +centers = detector.calc_cartesian_positions() + +shape = corners.shape[0:2] + +image = numpy.random.randint(256, size=shape) + +print_pixel(corners,image) diff --git a/contrib/python/common.py b/contrib/python/common.py new file mode 100644 index 0000000..ef40fee --- /dev/null +++ b/contrib/python/common.py @@ -0,0 +1,68 @@ +#!/usr/bin/env python3 +""" coding: utf-8 +""" + +from typing import Iterator, List, NamedTuple, Text, Tuple, Union + +import os + +from functools import partial + +import h5py + +from numpy import ndarray + +# H5Path data constructors +H5PathContains = NamedTuple("H5PathContains", [("path", Text)]) + +H5PathOptionalItemValue = NamedTuple('H5OptionalItemValue', [('path', Text), + ('default', float)]) # noqa +H5PathWithAttribute = NamedTuple("H5PathWithAttribute", [('attribute', Text), + ('value', bytes)]) + +H5Path = Union[H5PathContains, H5PathOptionalItemValue, H5PathWithAttribute] + +Parameter = NamedTuple("Parameter", [("name", Text), + ("value", float), + ("bounds", Tuple[float, float])]) + + +def _v_attrs(attribute: Text, value: Text, _name: Text, obj) -> h5py.Dataset: + """extract all the images and accumulate them in the acc variable""" + if isinstance(obj, h5py.Dataset): + if attribute in obj.attrs and obj.attrs[attribute] == value: + return obj + + +def _v_item(key: Text, name: Text, obj: h5py.Dataset) -> h5py.Dataset: + if key in name: + return obj + + +def get_shape(h5file: h5py.File, + item: H5Path) -> Tuple: + res = None + if isinstance(item, H5PathContains): + res = h5file.visititems(partial(_v_item, item.path)).shape + elif isinstance(item, H5PathOptionalItemValue): + _item = h5file.visititems(partial(_v_item, item.path)) + res = _item.shape if _item else (1,) + elif isinstance(item, H5PathWithAttribute): + res = h5file.visititems(partial(_v_attrs, + item.attribute, item.value)).shape + return res + + +def get_item_at_index(h5file: h5py.File, + item: H5Path, + index: int) -> Union[float, ndarray]: + res = None + if isinstance(item, H5PathContains): + res = h5file.visititems(partial(_v_item, item.path))[index] + elif isinstance(item, H5PathOptionalItemValue): + _item = h5file.visititems(partial(_v_item, item.path)) + res = _item.value if _item else item.default + elif isinstance(item, H5PathWithAttribute): + res = h5file.visititems(partial(_v_attrs, + item.attribute, item.value))[index] + return res diff --git a/contrib/python/mars.py b/contrib/python/mars.py new file mode 100644 index 0000000..a0ad404 --- /dev/null +++ b/contrib/python/mars.py @@ -0,0 +1,368 @@ +#!/usr/bin/env python3 +""" coding: utf-8 +# Il y a six types de fichiers à traiter. +# +# nb images | tz | poni +# ----------|---------| +# 5 | -1 | x +# 5 | 0 | scan3.poni +# 5 | x(-1) | x +# 5 | x(0) | scan3.poni +# 1 | -1 | x +# 1 | 0 | scan3.poni + +# In[7]: + +# # extraction des scans avec 6 positions en tx, quelque soit tz. +# # LONG +# import glob + +# files = glob.glob(os.path.join(ROOT, "*.nxs")) + +# def is_ok(filename: str) -> bool: +# with h5py.File(filename, mode='r') as f: +# for imgs, tx, tz in zip(get_images(f), get_tx(f), get_tz(f)): +# return True if tx.shape[0] == 6 else False + +# good = [f for f in files if is_ok(f)] +# print(good) +""" + +from typing import Iterator, List, NamedTuple, Text, Tuple, Union + +import os + +from functools import partial +from itertools import chain + +import h5py +import numpy +import pylab +import pyFAI + +from fabio.edfimage import edfimage +from numpy import ndarray + +from pyFAI.goniometer import GeometryTransformation, GoniometerRefinement +from pyFAI.gui import jupyter + +from common import * + + +ROOT = "/home/experiences/instrumentation/picca/jupyter/mars/20160800/" +PUBLISHED = os.path.join(ROOT, "published-data") + +CALIB = os.path.join(ROOT, "scan_3_01.nxs") + +# H5Path data constructors +MetaDataSource = NamedTuple("MetaDataSource", [("images", H5Path), + ("tx", H5Path), + ("tz", H5Path)]) + +MetaData = NamedTuple("MetaData", [("image", ndarray), + ("tx", float), + ("tz", float)]) + +_MultiCalibMarsTxTz = NamedTuple("_MultiCalibMarsTxTz", + [("filename", Text), + ("metasources", MetaDataSource), + ("idxs", List[int]), + ("calibrant", Text), + ("detector", Text), + ("wavelength", float)]) + +class MultiCalibMarsTxTz(_MultiCalibMarsTxTz): + def __len__(self) -> int: + with h5py.File(self.filename, mode='r') as f: + return get_shape(f, self.metasources.images)[0] + + def __item(self, f: h5py.File, index: int) -> MetaData: + return MetaData(get_item_at_index(f, + self.metasources.images, index), + get_item_at_index(f, + self.metasources.tx, index), + get_item_at_index(f, + self.metasources.tz, index)) + + def __item__(self, index: int) -> MetaData: + with h5py.File(self.filename, mode='r') as f: + return self.__item(f, index) + + def frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in self.idxs: + yield self.__item(f, index) + + def all_frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in range(len(self)): + yield self.__item(f, index) + + +def save_as_edf(calib: MultiCalibMarsTxTz, + basedir: Text) -> None: + """Save the multi calib images into edf files in order to do the first + calibration""" + for idx, metadata in zip(calib.idxs, calib.frames()): + base = os.path.splitext(os.path.basename(calib.filename))[0] + output = os.path.join(basedir, base + '_%d.edf' % (idx,)) + edfimage(metadata.image).write(output) + + +def optimize_with_new_images(multicalib: MultiCalibMarsTxTz, + gonioref: GoniometerRefinement, + calibrant: pyFAI.calibrant.Calibrant, + pts_per_deg: float=1) -> None: + """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 + + """ + sg = None + for idx, metadata in enumerate(multicalib.all_frames()): + print() + base = os.path.splitext(os.path.basename(multicalib.filename))[0] + + label = base + "_%d" % (idx,) + if label in gonioref.single_geometries: + continue + print(label) + + sg = gonioref.new_geometry(label, image=metadata.image, + metadata=metadata, calibrant=calibrant) + 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) # noqa + jupyter.display(sg=sg) + + +def get_wavelength(multicalib: MultiCalibMarsTxTz) -> float: + """Return the wavelength""" + return multicalib.wavelength + + +def get_calibrant(multicalib: MultiCalibMarsTxTz) -> pyFAI.calibrant.Calibrant: + """Return the calibrant with the right wavelength""" + calibrant = pyFAI.calibrant.get_calibrant(multicalib.calibrant) + calibrant.wavelength = get_wavelength(multicalib) + return calibrant + + +def get_detector(multicalib: MultiCalibMarsTxTz) -> pyFAI.Detector: + return pyFAI.detector_factory(multicalib.detector) + + +def calibration(json: str) -> None: + """Do a calibration with a bunch of images""" + + wavelength = 4.85945727522e-11 + + multicalib = MultiCalibMarsTxTz(os.path.join(ROOT, "scan_3_01.nxs"), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", 0.0)), # noqa + [2, 5, 8], "LaB6", "xpad_flat", wavelength) + + multicalib2 = MultiCalibMarsTxTz(os.path.join(ROOT, "scan_4_01.nxs"), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", -1.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + + # save all the ref as images in order to do the calibration with + # pyFAI-calib[2]. + save_as_edf(multicalib, PUBLISHED) + + # Definition of the geometry refinement: the parameter order is + # the same as the param_names + calibrant = get_calibrant(multicalib) + detector = get_detector(multicalib) + + distance = 0.258705917299 + poni1_scale = 0.001 + poni1_offset = 0.132825374721 + poni2_scale = 0.0012272727272727272 + poni2_offset = -0.9488181818181818 + rot1 = 0.00388272369359 + rot2 = -0.00942588451226 + rot3 = 7.19961198098e-07 + + parameters = [Parameter("dist", distance, (distance, distance)), + Parameter("poni1_offset", poni1_offset, (0, 0.2)), + Parameter("poni1_scale", poni1_scale, (0, 0.002)), + Parameter("poni2_offset", poni2_offset, (-1, -0.7)), + Parameter("poni2_scale", poni2_scale, (-1, 1)), + Parameter("rot1", rot1, (rot1, rot1)), + Parameter("rot2", rot2, (rot2, rot2)), + Parameter("rot3", rot3, (rot3, rot3))] + + params = {p.name: p.value for p in parameters} + bounds = {p.name: p.bounds for p in parameters} + param_names = [p.name for p in parameters] + + # Let's refine poni1 and poni2 also as function of the distance: + + trans_function = GeometryTransformation(param_names=param_names, + pos_names=["tx", "tz"], + dist_expr="dist", + poni1_expr="tz * poni1_scale + poni1_offset", # noqa + poni2_expr="tx * poni2_scale + poni2_offset", # noqa + rot1_expr="rot1", + rot2_expr="rot2", + rot3_expr="rot3") + + def pos_function(metadata: MetaData) -> Tuple[float, float]: + """Definition of the function reading the detector position from the + header of the image.""" + return metadata.tx, metadata.tz + + gonioref = GoniometerRefinement(params, # initial guess + bounds=bounds, + pos_function=pos_function, + trans_function=trans_function, + detector=detector, + wavelength=wavelength) + + print("Empty refinement object:") + print(gonioref) + + # Let's populate the goniometer refinement object with the known poni + for idx, metadata in zip(multicalib.idxs, multicalib.frames()): + base = os.path.splitext(os.path.basename(multicalib.filename))[0] + + label = base + "_%d" % (idx,) + control_points = os.path.join(PUBLISHED, base + "_%d.npt" % (idx,)) + ai = pyFAI.load(os.path.join(PUBLISHED, base + "_%d.poni" % (idx,))) # noqa + print(ai) + + gonioref.new_geometry(label, metadata.image, metadata, + control_points, calibrant, ai) + + print("Filled refinement object:") + print(gonioref) + print(os.linesep + "\tlabel \t tx") + for k, v in gonioref.single_geometries.items(): + print(k, v.get_position()) + + for g in gonioref.single_geometries.values(): + ai = gonioref.get_ai(g.get_position()) + print(ai) + + for sg in gonioref.single_geometries.values(): + jupyter.display(sg=sg) + + gonioref.refine2() + + for multi in [multicalib, multicalib2]: + optimize_with_new_images(multi, gonioref, calibrant) + + # for idx, sg in enumerate(gonioref.single_geometries.values()): + # sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param) + # jupyter.display(sg=sg) + + gonioref.save(json) + + # pylab.show() + + +def _integrate(json: Text, multicalib: Tuple[MultiCalibMarsTxTz, MultiCalibMarsTxTz]) -> None: + # do not do the computation if the .dat already exist + output = multicalib[0].filename + '.dat' + if os.path.exists(output): + return + + THRESHOLD = 12000 + """Integrate a file with a json calibration file""" + gonio = pyFAI.goniometer.Goniometer.sload(json) + + images = [] + positions = [] + for metadata in chain(multicalib[0].all_frames(), + multicalib[1].all_frames()): + images.append(metadata.image) + positions.append((metadata.tx, metadata.tz)) + mai = gonio.get_mg(positions) + + # compute the mask + detector = get_detector(multicalib[0]) + mask = numpy.array(detector.mask) + lst_mask = [] + for img in images: # remove all pixels above the threshold" + if THRESHOLD is not None: + mask_t = numpy.where(img > THRESHOLD, True, False) + lst_mask.append(numpy.logical_or(mask, mask_t)) + else: + lst_mask.append(mask) + + res = mai.integrate1d(images, 10000, lst_mask=lst_mask) + numpy.savetxt(output, numpy.array(res).T) + + #jupyter.plot1d(res) + #pylab.show() + + +def integrate(json: Text, mcals: List[Tuple[MultiCalibMarsTxTz, MultiCalibMarsTxTz]]) -> None: + for mcal in mcals: + try: + print(mcal[0].filename, mcal[1].filename) + _integrate(json, mcal) + except: + pass + +def main(): + wavelength = 4.85945727522e-11 + lab6 = MultiCalibMarsTxTz(os.path.join(ROOT, "scan_3_01.nxs"), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", 0.0)), # noqa + [2, 5, 8], "LaB6", "xpad_flat", wavelength) + + lab6_2 = MultiCalibMarsTxTz(os.path.join(ROOT, "scan_4_01.nxs"), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", -1.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + + JSON = os.path.join(PUBLISHED, "calibration.json") + #calibration(JSON, [lab6, lab6_2]) + + # integration des echantillons + tz1 = [ MultiCalibMarsTxTz(os.path.join(ROOT, "scan_%d_01.nxs" % (i,)), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa" + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", -1.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + for i in [77, 79, 81, 83, 85, 87, 89, 91]] + + + tz0 = [ MultiCalibMarsTxTz(os.path.join(ROOT, "scan_%d_01.nxs" % (i,)), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa" + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", 0.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + for i in [78, 80, 82, 84, 86, 88, 90, 92]] + + tz3 = [ MultiCalibMarsTxTz(os.path.join(ROOT, "scan_%d_01.nxs" % (i,)), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa" + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", 0.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + for i in range(399, 588, 2) if i not in [523, 527, 581]] + + tz5 = [ MultiCalibMarsTxTz(os.path.join(ROOT, "scan_%d_01.nxs" % (i,)), + MetaDataSource(H5PathWithAttribute("interpretation", b"image"), # noqa" + H5PathContains("scan_data/actuator_1_1"), # noqa + H5PathOptionalItemValue("MARS/D03-1-CX0__DT__DTC_2D-MT_Tz__#1/raw_value", -5.0)), # noqa + [], "LaB6", "xpad_flat", wavelength) + for i in range(400, 589, 2) if i not in [524, 528, 582]] + + # samples = [(lab6, lab6_2)] + list(zip(tz1, tz0)) + list(zip(tz3, tz5)) + samples = list(zip(tz3, tz5)) + integrate(JSON, samples) + pylab.show() + +if __name__ == "__main__": + main() diff --git a/contrib/python/mythen.py b/contrib/python/mythen.py new file mode 100644 index 0000000..cad3a43 --- /dev/null +++ b/contrib/python/mythen.py @@ -0,0 +1,306 @@ +#!/usr/bin/env python3 +""" coding: utf-8 +# Il y a six types de fichiers à traiter. +# +# nb images | tz | poni +# ----------|---------| +# 5 | -1 | x +# 5 | 0 | scan3.poni +# 5 | x(-1) | x +# 5 | x(0) | scan3.poni +# 1 | -1 | x +# 1 | 0 | scan3.poni + +# In[7]: + +# # extraction des scans avec 6 positions en tx, quelque soit tz. +# # LONG +# import glob + +# files = glob.glob(os.path.join(ROOT, "*.nxs")) + +# def is_ok(filename: str) -> bool: +# with h5py.File(filename, mode='r') as f: +# for imgs, tx, tz in zip(get_images(f), get_tx(f), get_tz(f)): +# return True if tx.shape[0] == 6 else False + +# good = [f for f in files if is_ok(f)] +# print(good) +""" + +from typing import Iterator, List, NamedTuple, Text, Tuple, Union + +import os + +from functools import partial +from math import pi, radians + +import h5py +import pylab +import pyFAI + +from fabio.edfimage import edfimage +from numpy import ndarray + +from pyFAI.goniometer import GeometryTransformation, GoniometerRefinement +from pyFAI.gui import jupyter + +from common import * + +ROOT = "/home/akira/Downloads" +PUBLISHED = os.path.join(ROOT, "published-data") + +CALIB = os.path.join(ROOT, "scan_3_01.nxs") + +MetaDataSource = NamedTuple("MetaDataSource", [("spectrum", H5Path), + ("tth", H5Path)]) + +MetaData = NamedTuple("MetaData", [("spectrum", ndarray), + ("tth", float)]) + +_Mythen = NamedTuple("_Mythen", [("filename", Text), + ("metasources", MetaDataSource), + ("idxs", List[int]), + ("module", int), + ("calibrant", Text), + ("detector", Text), + ("wavelength", float)]) + +Parameter = NamedTuple("Parameter", [("name", Text), + ("value", float), + ("bounds", Tuple[float, float])]) + + +def extract_module(spectrum: ndarray, module: int) -> ndarray: + return spectrum[slice(module*1280, (module+1)*1280, 1)] + + +class Mythen(_Mythen): + def __len__(self) -> int: + with h5py.File(self.filename, mode='r') as f: + return get_shape(f, self.metasources.spectrum)[0] + + def __item(self, f: h5py.File, index: int) -> MetaData: + spectrum = extract_module(get_item_at_index(f, self.metasources.spectrum, index), self.module) + spectrum.shape = (1, spectrum.shape[0]) + return MetaData(spectrum, + get_item_at_index(f, + self.metasources.tth, index)) + + def __item__(self, index: int) -> MetaData: + with h5py.File(self.filename, mode='r') as f: + return self.__item(f, index) + + def frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in self.idxs: + yield self.__item(f, index) + + def all_frames(self) -> Iterator[MetaData]: + with h5py.File(self.filename, mode='r') as f: + for index in range(len(self)): + yield self.__item(f, index) + + +def save_as_edf(calib: Mythen, + basedir: Text) -> None: + """Save the multi calib images into edf files in order to do the first + calibration""" + for idx, metadata in zip(calib.idxs, calib.frames()): + base = os.path.splitext(os.path.basename(calib.filename))[0] + output = os.path.join(basedir, base + '_%d.edf' % (idx,)) + edfimage(metadata.spectrum).write(output) + + +def get_wavelength(multicalib: Mythen) -> float: + """Return the wavelength""" + return multicalib.wavelength + + +def get_calibrant(multicalib: Mythen) -> pyFAI.calibrant.Calibrant: + """Return the calibrant with the right wavelength""" + calibrant = pyFAI.calibrant.get_calibrant(multicalib.calibrant) + calibrant.wavelength = get_wavelength(multicalib) + return calibrant + + +def get_detector(multicalib: Mythen) -> pyFAI.Detector: + return pyFAI.detector_factory(multicalib.detector) + + +def optimize_with_new_images(multicalib: Mythen, + gonioref: GoniometerRefinement, + calibrant: pyFAI.calibrant.Calibrant, + pts_per_deg: float=1) -> None: + """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 + + """ + sg = None + for idx, metadata in enumerate(multicalib.all_frames()): + print() + base = os.path.splitext(os.path.basename(multicalib.filename))[0] + + label = base + "_%d" % (idx,) + if label in gonioref.single_geometries: + continue + print(label) + + sg = gonioref.new_geometry(label, image=metadata.spectrum, + metadata=metadata, calibrant=calibrant) + print(sg.extract_cp(method="blob", pts_per_deg=pts_per_deg))# , "watershed")) + print("*"*50) + gonioref.refine2() + if sg: + sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param) # noqa + jupyter.display(sg=sg) + +def calibration(json: str) -> None: + """Do a calibration with a bunch of images""" + + filename = os.path.join(ROOT, "tth2C-16_11_15-55-22_008.nxs") + wavelength = 0.67186e-10 + multicalib = Mythen(filename, + MetaDataSource(H5PathWithAttribute("interpretation", b"spectrum"), + H5PathContains("scan_data/actuator_1_1")), + [18, 21], 0, "LaB6", "mythen", wavelength) + + save_as_edf(multicalib, PUBLISHED) + + # Definition of the geometry refinement: the parameter order is + # the same as the param_names + calibrant = get_calibrant(multicalib) + detector = get_detector(multicalib) + + distance = 0.617400891837 + poni1 = 0 + poni2 = 0.0323014291288 + rot1_scale = -1 * pi / 180. + rot1_offset = radians(74.88) + rot2 = 0 + rot3 = 0 + + parameters = [Parameter("dist", distance, (distance, distance)), + Parameter("poni1", poni1, (poni1, poni1)), + Parameter("poni2", poni2, (poni2, poni2)), + Parameter("rot1_offset", rot1_offset, (rot1_offset, rot1_offset)), + Parameter("rot1_scale", rot1_scale, (rot1_scale, rot1_scale)), + Parameter("rot2", rot2, (rot2, rot2)), + Parameter("rot3", rot3, (rot3, rot3))] + + params = {p.name: p.value for p in parameters} + bounds = {p.name: p.bounds for p in parameters} + param_names = [p.name for p in parameters] + + # Let's refine poni1 and poni2 also as function of the distance: + + trans_function = GeometryTransformation(param_names=param_names, + pos_names=["tth"], + dist_expr="dist", + poni1_expr="poni1", + poni2_expr="poni2", + rot1_expr="tth * rot1_scale + rot1_offset", + rot2_expr="rot2", + rot3_expr="rot3") + + def pos_function(metadata: MetaData) -> Tuple[float]: + """Definition of the function reading the detector position from the + header of the image.""" + return metadata.tth, + + gonioref = GoniometerRefinement(params, # initial guess + bounds=bounds, + pos_function=pos_function, + trans_function=trans_function, + detector=detector, + wavelength=wavelength) + + print("Empty refinement object:") + print(gonioref) + + # Let's populate the goniometer refinement object with the known poni + for idx, metadata in zip(multicalib.idxs, multicalib.frames()): + base = os.path.splitext(os.path.basename(multicalib.filename))[0] + + label = base + "_%d" % (idx,) + control_points = os.path.join(PUBLISHED, base + "_%d.npt" % (idx,)) + ai = pyFAI.load(os.path.join(PUBLISHED, base + "_%d.poni" % (idx,))) # noqa + print(ai) + + gonioref.new_geometry(label, metadata.spectrum, metadata, + control_points, calibrant, ai) + + print("Filled refinement object:") + print(gonioref) + print(os.linesep + "\tlabel \t tx") + for k, v in gonioref.single_geometries.items(): + print(k, v.get_position()) + + for g in gonioref.single_geometries.values(): + ai = gonioref.get_ai(g.get_position()) + print(ai) + + for sg in gonioref.single_geometries.values(): + jupyter.display(sg=sg) + + gonioref.refine2() + + # for multi in [multicalib]: + # optimize_with_new_images(multi, gonioref, calibrant) + + # for idx, sg in enumerate(gonioref.single_geometries.values()): + # sg.geometry_refinement.set_param(gonioref.get_ai(sg.get_position()).param) + # jupyter.display(sg=sg) + + gonioref.save(json) + + # pylab.show() + + +def integrate(json: Text) -> None: + """Integrate a file with a json calibration file""" + gonio = pyFAI.goniometer.Goniometer.sload(json) + filename = os.path.join(ROOT, "tth2C-16_11_15-55-22_008.nxs") + wavelength = 0.67186e-10 + multicalib = Mythen(filename, + MetaDataSource(H5PathWithAttribute("interpretation", b"spectrum"), + H5PathContains("scan_data/actuator_1_1")), + [18, 21], 0, "LaB6", "mythen", wavelength) + + # print(len(multicalib)) + # metadata = multicalib.__item__(21) + + # ai = gonio.get_ai(metadata.tth) + # res = ai.integrate1d(metadata.spectrum, 1280, unit="2th_deg") + # jupyter.plot1d(res) + + images = [] + positions = [] + for metadata in multicalib.all_frames(): + images.append(metadata.spectrum) + positions.append((metadata.tth,)) + mai = gonio.get_mg(positions) + res = mai.integrate1d(images, 30000) + jupyter.plot1d(res) + pylab.show() + + +def main() -> None: + filename = os.path.join(ROOT, "tth2C-16_11_15-55-22_008.nxs") + wavelength = 0.67186e-10 + mythen = Mythen(filename, + MetaDataSource(H5PathWithAttribute("interpretation", b"spectrum"), + H5PathContains("scan_data/actuator_1_1")), + [18, 21], 0, "LaB6", "mythen", wavelength) + + save_as_edf(mythen, PUBLISHED) + + for frame in mythen.all_frames(): + print(frame) + +if __name__ == "__main__": + json = os.path.join(PUBLISHED, "mythen1.json") + calibration(json) + integrate(json) diff --git a/contrib/python/swing/pinhole1.smv b/contrib/python/swing/pinhole1.smv Binary files differnew file mode 100644 index 0000000..cfdee04 --- /dev/null +++ b/contrib/python/swing/pinhole1.smv diff --git a/contrib/python/swing/plot.py b/contrib/python/swing/plot.py new file mode 100644 index 0000000..fce42ac --- /dev/null +++ b/contrib/python/swing/plot.py @@ -0,0 +1,12 @@ +import pylab +import numpy + +FILE = "pinhole1.smv" + +with open(FILE, "rb") as f: + f.seek(512) + img = numpy.fromfile(f, dtype='uint16', sep="") + print(img.shape) + img.shape = ((4096, 4096)) + pylab.imshow(img, interpolate='nearest') + pylab.show() |