From 745e50cb167be93313a2417f6d809d63a1de1a42 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 11:48:43 +0100 Subject: [PATCH 01/19] correction exit and scattering angles --- src/pyFAI/units.py | 62 ++++++++++++++++++++++++++++++++++------------ 1 file changed, 46 insertions(+), 16 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 712044886..d476b1349 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -349,7 +349,7 @@ def eq_q(x, y, z, wavelength): return 4.0e-9 * numpy.pi * numpy.sin(eq_2th(x, y, z) / 2.0) / wavelength -def eq_exitangle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): +def eq_scattering_angle_vertical(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): """Calculates the vertical exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position @@ -361,7 +361,29 @@ def eq_exitangle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, s return numpy.arctan2(y, numpy.sqrt(z ** 2 + x ** 2)) -def eq_exitangle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): +def eq_exit_angle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): + """Calculates the vertical exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction + + :param x: horizontal position, towards the center of the ring, from sample position + :param y: vertical position, to the roof, from sample position + :param z: distance from sample along the beam + :param wavelength: in meter + :return: vertical exit angle in radians + """ + rot_incident_angle = numpy.array([[1,0,0], + [0,numpy.cos(incident_angle), numpy.sin(-incident_angle)], + [0, numpy.sin(incident_angle), numpy.cos(incident_angle)]], + ) + rot_tilt_angle = numpy.array([[numpy.cos(tilt_angle), numpy.sin(-tilt_angle), 0], + [numpy.sin(tilt_angle), numpy.cos(tilt_angle), 0], + [0, 0, 1]], + ) + rotated_xyz = numpy.tensordot(rot_incident_angle, numpy.stack((x,y,z)), axes=1) + xp, yp, zp = numpy.tensordot(rot_tilt_angle, rotated_xyz, axes=1) + return numpy.arctan2(yp, numpy.sqrt(zp ** 2 + xp ** 2)) + + +def eq_scattering_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): """Calculates the horizontal exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position @@ -382,9 +404,9 @@ def q_lab_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sam :param wavelength: in meter :return: horizontal scattering vector in inverse nm """ - exit_angle = eq_exitangle(x=x, y=y, z=z, wavelength=wavelength) - exit_angle_horz = eq_exitangle_horz(x=x, y=y, z=z, wavelength=wavelength) - return 2.0e-9 / wavelength * numpy.pi * numpy.cos(exit_angle) * numpy.sin(exit_angle_horz) + scattering_angle_vertical = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) + scattering_angle_horz = eq_scattering_angle_horz(x=x, y=y, z=z, wavelength=wavelength) + return 2.0e-9 / wavelength * numpy.pi * numpy.cos(scattering_angle_vertical) * numpy.sin(scattering_angle_horz) def q_lab_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): @@ -396,7 +418,7 @@ def q_lab_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sam :param wavelength: in meter :return: vertical scattering vector in inverse nm """ - exit_angle = eq_exitangle(x=x, y=y, z=z, wavelength=wavelength) + exit_angle = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) return 2.0e-9 / wavelength * numpy.pi * numpy.sin(exit_angle) @@ -409,8 +431,8 @@ def q_lab_beam(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sam :param wavelength: in meter :return: beam scattering vector in inverse nm """ - exit_angle = eq_exitangle(x=x, y=y, z=z, wavelength=wavelength) - exit_angle_horz = eq_exitangle_horz(x=x, y=y, z=z, wavelength=wavelength) + exit_angle = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) + exit_angle_horz = eq_scattering_angle_horz(x=x, y=y, z=z, wavelength=wavelength) return 2.0e-9 / wavelength * numpy.pi * (numpy.cos(exit_angle) * numpy.cos(exit_angle_horz) - 1) @@ -877,19 +899,27 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o unit_symbol="nm^{-1}", positive=False) -register_radial_fiber_unit("exitangle_rad", +register_radial_fiber_unit("scattering_angle_vert", scale=1.0, - label=r"Exit scattering angle (rad)", - equation=eq_exitangle, - short_name="exitangle", + label=r"Vertical scattering angle (rad)", + equation=eq_scattering_angle_vertical, + short_name="scatangle_vert", + unit_symbol="rad", + positive=False) + +register_radial_fiber_unit("scattering_angle_horz", + scale=1.0, + label=r"Horizontal scattering angle(rad)", + equation=eq_scattering_angle_horz, + short_name="scatangle_horz", unit_symbol="rad", positive=False) -register_radial_fiber_unit("horz_exitangle_rad", +register_radial_fiber_unit("exit_angle", scale=1.0, - label=r"Exit scattering angle (rad) in the horizontal axis", - equation=eq_exitangle_horz, - short_name="exitangle_horz", + label=r"Exit angle(rad)", + equation=eq_exit_angle, + short_name="exitangle", unit_symbol="rad", positive=False) From e7203eb0d53ccb59386ed3bcce48183efbf01e77 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 17 Jan 2025 10:52:33 +0000 Subject: [PATCH 02/19] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/pyFAI/units.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index d476b1349..e8cf8e25b 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -373,7 +373,7 @@ def eq_exit_angle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, rot_incident_angle = numpy.array([[1,0,0], [0,numpy.cos(incident_angle), numpy.sin(-incident_angle)], [0, numpy.sin(incident_angle), numpy.cos(incident_angle)]], - ) + ) rot_tilt_angle = numpy.array([[numpy.cos(tilt_angle), numpy.sin(-tilt_angle), 0], [numpy.sin(tilt_angle), numpy.cos(tilt_angle), 0], [0, 0, 1]], From b8369b6f6990c1ed2ea50e10639ecb28294ca7d9 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 11:53:07 +0100 Subject: [PATCH 03/19] compatible --- src/pyFAI/units.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index d476b1349..3700db97f 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -361,7 +361,7 @@ def eq_scattering_angle_vertical(x, y, z, wavelength=None, incident_angle=0.0, t return numpy.arctan2(y, numpy.sqrt(z ** 2 + x ** 2)) -def eq_exit_angle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): +def eq_exitangle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): """Calculates the vertical exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position @@ -918,7 +918,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle", scale=1.0, label=r"Exit angle(rad)", - equation=eq_exit_angle, + equation=eq_exitangle, short_name="exitangle", unit_symbol="rad", positive=False) From b9494ee7a39e7c8eaa50f6a668eb835d0544d5a5 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 13:32:38 +0100 Subject: [PATCH 04/19] horizontal exit angle --- src/pyFAI/units.py | 57 ++++++++++++++++++++++++++++++++++++++-------- 1 file changed, 47 insertions(+), 10 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index cfdf1f4fc..5317bd767 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -350,7 +350,7 @@ def eq_q(x, y, z, wavelength): def eq_scattering_angle_vertical(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): - """Calculates the vertical exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction + """Calculates the vertical scattering angle (relative to direct beam axis), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position :param y: vertical position, to the roof, from sample position @@ -361,13 +361,27 @@ def eq_scattering_angle_vertical(x, y, z, wavelength=None, incident_angle=0.0, t return numpy.arctan2(y, numpy.sqrt(z ** 2 + x ** 2)) -def eq_exitangle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): - """Calculates the vertical exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction +def eq_scattering_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): + """Calculates the horizontal scattering angle (relative to direct beam axis), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position :param y: vertical position, to the roof, from sample position :param z: distance from sample along the beam :param wavelength: in meter + :return: horizontal exit angle in radians + """ + return numpy.arctan2(x, z) + + +def eq_exit_angle_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): + """Calculates the vertical exit angle relative to the horizon (for thin films), used for GI/Fiber diffraction + + :param x: horizontal position, towards the center of the ring, from sample position + :param y: vertical position, to the roof, from sample position + :param z: distance from sample along the beam + :param incident_angle: tilting of the sample towards the beam (analog to rot2): in radians + :param tilt_angle: tilting of the sample orthogonal to the beam direction (analog to rot3): in radians + :param wavelength: in meter :return: vertical exit angle in radians """ rot_incident_angle = numpy.array([[1,0,0], @@ -383,16 +397,31 @@ def eq_exitangle(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, s return numpy.arctan2(yp, numpy.sqrt(zp ** 2 + xp ** 2)) -def eq_scattering_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): - """Calculates the horizontal exit scattering angle (relative to direct beam axis), used for GI/Fiber diffraction +def eq_exit_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): + """Calculates the horizontal exit angle relative to the horizon (for thin films), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position :param y: vertical position, to the roof, from sample position :param z: distance from sample along the beam + :param incident_angle: tilting of the sample towards the beam (analog to rot2): in radians + :param tilt_angle: tilting of the sample orthogonal to the beam direction (analog to rot3): in radians :param wavelength: in meter :return: horizontal exit angle in radians """ - return numpy.arctan2(x, z) + rot_incident_angle = numpy.array([[1,0,0], + [0,numpy.cos(incident_angle), numpy.sin(-incident_angle)], + [0, numpy.sin(incident_angle), numpy.cos(incident_angle)]], + ) + rot_tilt_angle = numpy.array([[numpy.cos(tilt_angle), numpy.sin(-tilt_angle), 0], + [numpy.sin(tilt_angle), numpy.cos(tilt_angle), 0], + [0, 0, 1]], + ) + rotated_xyz = numpy.tensordot(rot_incident_angle, numpy.stack((x,y,z)), axes=1) + xp, yp, zp = numpy.tensordot(rot_tilt_angle, rotated_xyz, axes=1) + return numpy.arctan2(xp, zp) + + +eq_exitangle = eq_exit_angle_vert def q_lab_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): @@ -909,17 +938,25 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("scattering_angle_horz", scale=1.0, - label=r"Horizontal scattering angle(rad)", + label=r"Horizontal scattering angle (rad)", equation=eq_scattering_angle_horz, short_name="scatangle_horz", unit_symbol="rad", positive=False) -register_radial_fiber_unit("exit_angle", +register_radial_fiber_unit("exit_angle_vert", + scale=1.0, + label=r"Vertical exit angle (rad)", + equation=eq_exit_angle_vert, + short_name="exitangle_vert", + unit_symbol="rad", + positive=False) + +register_radial_fiber_unit("exit_angle_horz", scale=1.0, label=r"Exit angle(rad)", - equation=eq_exitangle, - short_name="exitangle", + equation=eq_exit_angle_horz, + short_name="exitangle_horz", unit_symbol="rad", positive=False) From 15b1ba3bae2ce85c59c00c5c722414d4c9421292 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 13:34:30 +0100 Subject: [PATCH 05/19] unifity terms --- src/pyFAI/units.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 5317bd767..5f075808d 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -447,8 +447,8 @@ def q_lab_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sam :param wavelength: in meter :return: vertical scattering vector in inverse nm """ - exit_angle = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) - return 2.0e-9 / wavelength * numpy.pi * numpy.sin(exit_angle) + scattering_angle_vertical = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) + return 2.0e-9 / wavelength * numpy.pi * numpy.sin(scattering_angle_vertical) def q_lab_beam(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): @@ -460,9 +460,9 @@ def q_lab_beam(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sam :param wavelength: in meter :return: beam scattering vector in inverse nm """ - exit_angle = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) - exit_angle_horz = eq_scattering_angle_horz(x=x, y=y, z=z, wavelength=wavelength) - return 2.0e-9 / wavelength * numpy.pi * (numpy.cos(exit_angle) * numpy.cos(exit_angle_horz) - 1) + scattering_angle_vertical = eq_scattering_angle_vertical(x=x, y=y, z=z, wavelength=wavelength) + scattering_angle_horz = eq_scattering_angle_horz(x=x, y=y, z=z, wavelength=wavelength) + return 2.0e-9 / wavelength * numpy.pi * (numpy.cos(scattering_angle_vertical) * numpy.cos(scattering_angle_horz) - 1) def q_lab(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): From 053ebb2cfeed158ea8f271cac4a169dcc3f03142 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 13:41:26 +0100 Subject: [PATCH 06/19] exit angles in degrees --- src/pyFAI/units.py | 22 +++++++++++++++++++--- 1 file changed, 19 insertions(+), 3 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 5f075808d..c6109051e 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -374,7 +374,7 @@ def eq_scattering_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_ def eq_exit_angle_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): - """Calculates the vertical exit angle relative to the horizon (for thin films), used for GI/Fiber diffraction + """Calculates the vertical exit angle in radians relative to the horizon (for thin films), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position :param y: vertical position, to the roof, from sample position @@ -398,7 +398,7 @@ def eq_exit_angle_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle= def eq_exit_angle_horz(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): - """Calculates the horizontal exit angle relative to the horizon (for thin films), used for GI/Fiber diffraction + """Calculates the horizontal exit angle in radians relative to the horizon (for thin films), used for GI/Fiber diffraction :param x: horizontal position, towards the center of the ring, from sample position :param y: vertical position, to the roof, from sample position @@ -954,12 +954,28 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle_horz", scale=1.0, - label=r"Exit angle(rad)", + label=r"Horizontal exit angle (rad)", equation=eq_exit_angle_horz, short_name="exitangle_horz", unit_symbol="rad", positive=False) +register_radial_fiber_unit("exit_angle_vert_deg", + scale=180.0 / numpy.pi, + label=r"Vertical exit angle (deg)", + equation=eq_exit_angle_vert, + short_name="exitangle_vert_deg", + unit_symbol="deg", + positive=False) + +register_radial_fiber_unit("exit_angle_horz_deg", + scale=180.0 / numpy.pi, + label=r"Horizontal exit angle (deg)", + equation=eq_exit_angle_horz, + short_name="exitangle_horz_deg", + unit_symbol="deg", + positive=False) + register_radial_fiber_unit("qxgi_nm^-1", scale=1.0, label=r"Scattering vector $q_x$ ($nm^{-1}$)", From 923b5ad6e7224caf379b6a51e2432a375b499ebe Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 13:44:14 +0100 Subject: [PATCH 07/19] default exit angle is degrees --- src/pyFAI/units.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index c6109051e..1f1712fa5 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -944,35 +944,35 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o unit_symbol="rad", positive=False) -register_radial_fiber_unit("exit_angle_vert", +register_radial_fiber_unit("exit_angle_vert_rad", scale=1.0, label=r"Vertical exit angle (rad)", equation=eq_exit_angle_vert, - short_name="exitangle_vert", + short_name="exitangle_vert_rad", unit_symbol="rad", positive=False) -register_radial_fiber_unit("exit_angle_horz", +register_radial_fiber_unit("exit_angle_horz_rad", scale=1.0, label=r"Horizontal exit angle (rad)", equation=eq_exit_angle_horz, - short_name="exitangle_horz", + short_name="exitangle_horz_rad", unit_symbol="rad", positive=False) -register_radial_fiber_unit("exit_angle_vert_deg", +register_radial_fiber_unit("exit_angle_vert", scale=180.0 / numpy.pi, label=r"Vertical exit angle (deg)", equation=eq_exit_angle_vert, - short_name="exitangle_vert_deg", + short_name="exitangle_vert", unit_symbol="deg", positive=False) -register_radial_fiber_unit("exit_angle_horz_deg", +register_radial_fiber_unit("exit_angle_horz", scale=180.0 / numpy.pi, label=r"Horizontal exit angle (deg)", equation=eq_exit_angle_horz, - short_name="exitangle_horz_deg", + short_name="exitangle_horz", unit_symbol="deg", positive=False) From ee000a6ddde4a777245c0edc1d67c1253eed4798 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 14:40:08 +0100 Subject: [PATCH 08/19] attemp for numexpr --- src/pyFAI/units.py | 24 ++++++++++++++++-------- 1 file changed, 16 insertions(+), 8 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 1f1712fa5..7bdfca777 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -638,9 +638,6 @@ def eq_qhorz_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_ :return: component of the scattering vector along the horizontal direction in inverse nm """ hpos, vpos = rotate_sample_orientation(x=x, y=y, sample_orientation=sample_orientation) -# The above code is using multi-line comments in Python, which are denoted by three consecutive pound -# signs ( - return eq_qhorz(hpos=hpos, vpos=vpos, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle) @@ -744,11 +741,20 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula_qx = f"4.0e-9*π/λ*sin(arctan2(x, z)/2.0)" # TODO: wrong, fix me formula_qy = f"4.0e-9*π/λ*sin(arctan2(y, z)/2.0)" # TODO: wrong, fix me -formula_exit_angle = "arctan2(y, sqrt(z*z + x*x))" -formula_exit_angle_horz = "arctan2(x,z)" -formula_qbeam_lab = f"2.0e-9/λ*π*(cos({formula_exit_angle})*cos({formula_exit_angle_horz}) - 1)" -formula_qhorz_lab = f"2.0e-9/λ*π*cos({formula_exit_angle})*sin({formula_exit_angle_horz})" -formula_qvert_lab = f"2.0e-9/λ*π*sin({formula_exit_angle})" +formula_scattering_angle_vert = "arctan2(y, sqrt(z*z + x*x))" +formula_scattering_angle_horz = "arctan2(x,z)" +# formula_x_rot_iangle = "x" +# formula_y_rot_iangle = "y * cos(η) - z * sin(η)" +# formula_z_rot_iangle = "y * sin(η) + z * cos(η)" +# formula_x_rot_tangle = f"({formula_x_rot_iangle}) * cos(χ) - ({formula_y_rot_iangle}) * sin(χ)" +# formula_y_rot_tangle = f"({formula_x_rot_iangle}) * sin(χ) + ({formula_y_rot_iangle}) * cos(χ)" +# formula_z_rot_tangle = formula_z_rot_iangle +# formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} ** {formula_z_rot_tangle} + {formula_x_rot_tangle} ** {formula_x_rot_tangle}))" +# formula_exit_angle_horz = f"arctan2(({formula_x_rot_tangle}), ({formula_z_rot_tangle}))" + +formula_qbeam_lab = f"2.0e-9/λ*π*(cos({formula_scattering_angle_vert})*cos({formula_scattering_angle_horz}) - 1)" +formula_qhorz_lab = f"2.0e-9/λ*π*cos({formula_scattering_angle_vert})*sin({formula_scattering_angle_horz})" +formula_qvert_lab = f"2.0e-9/λ*π*sin({formula_scattering_angle_vert})" formula_qbeam_rot = f"cos(η)*({formula_qbeam_lab})+sin(η)*({formula_qvert_lab})" formula_qhorz_rot = f"cos(χ)*({formula_qhorz_lab})-sin(χ)*sin(η)*({formula_qbeam_lab})+sin(χ)*cos(η)*({formula_qvert_lab})" formula_qvert_rot = f"-sin(χ)*({formula_qhorz_lab})-cos(χ)*sin(η)*({formula_qbeam_lab})+cos(χ)*cos(η)*({formula_qvert_lab})" @@ -931,6 +937,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("scattering_angle_vert", scale=1.0, label=r"Vertical scattering angle (rad)", + # formula=formula_scattering_angle_vert, equation=eq_scattering_angle_vertical, short_name="scatangle_vert", unit_symbol="rad", @@ -939,6 +946,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("scattering_angle_horz", scale=1.0, label=r"Horizontal scattering angle (rad)", + # formula=formula_scattering_angle_horz, equation=eq_scattering_angle_horz, short_name="scatangle_horz", unit_symbol="rad", From a9c2d991a564d3be510bc07552e3588dbf1c9c93 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 16:13:56 +0100 Subject: [PATCH 09/19] all numexpr --- src/pyFAI/units.py | 24 ++++++++++++++---------- 1 file changed, 14 insertions(+), 10 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 7bdfca777..992471367 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -743,14 +743,14 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula_scattering_angle_vert = "arctan2(y, sqrt(z*z + x*x))" formula_scattering_angle_horz = "arctan2(x,z)" -# formula_x_rot_iangle = "x" -# formula_y_rot_iangle = "y * cos(η) - z * sin(η)" -# formula_z_rot_iangle = "y * sin(η) + z * cos(η)" -# formula_x_rot_tangle = f"({formula_x_rot_iangle}) * cos(χ) - ({formula_y_rot_iangle}) * sin(χ)" -# formula_y_rot_tangle = f"({formula_x_rot_iangle}) * sin(χ) + ({formula_y_rot_iangle}) * cos(χ)" -# formula_z_rot_tangle = formula_z_rot_iangle -# formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} ** {formula_z_rot_tangle} + {formula_x_rot_tangle} ** {formula_x_rot_tangle}))" -# formula_exit_angle_horz = f"arctan2(({formula_x_rot_tangle}), ({formula_z_rot_tangle}))" +formula_x_rot_iangle = "x" +formula_y_rot_iangle = "y * cos(η) - z * sin(η)" +formula_z_rot_iangle = "y * sin(η) + z * cos(η)" +formula_x_rot_tangle = f"({formula_x_rot_iangle}) * cos(χ) - ({formula_y_rot_iangle}) * sin(χ)" +formula_y_rot_tangle = f"({formula_x_rot_iangle}) * sin(χ) + ({formula_y_rot_iangle}) * cos(χ)" +formula_z_rot_tangle = formula_z_rot_iangle +formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} ** {formula_z_rot_tangle} + {formula_x_rot_tangle} ** {formula_x_rot_tangle}))" +formula_exit_angle_horz = f"arctan2(({formula_x_rot_tangle}), ({formula_z_rot_tangle}))" formula_qbeam_lab = f"2.0e-9/λ*π*(cos({formula_scattering_angle_vert})*cos({formula_scattering_angle_horz}) - 1)" formula_qhorz_lab = f"2.0e-9/λ*π*cos({formula_scattering_angle_vert})*sin({formula_scattering_angle_horz})" @@ -937,7 +937,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("scattering_angle_vert", scale=1.0, label=r"Vertical scattering angle (rad)", - # formula=formula_scattering_angle_vert, + formula=formula_scattering_angle_vert, equation=eq_scattering_angle_vertical, short_name="scatangle_vert", unit_symbol="rad", @@ -946,7 +946,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("scattering_angle_horz", scale=1.0, label=r"Horizontal scattering angle (rad)", - # formula=formula_scattering_angle_horz, + formula=formula_scattering_angle_horz, equation=eq_scattering_angle_horz, short_name="scatangle_horz", unit_symbol="rad", @@ -955,6 +955,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle_vert_rad", scale=1.0, label=r"Vertical exit angle (rad)", + formula=formula_exit_angle_vert, equation=eq_exit_angle_vert, short_name="exitangle_vert_rad", unit_symbol="rad", @@ -963,6 +964,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle_horz_rad", scale=1.0, label=r"Horizontal exit angle (rad)", + formula=formula_exit_angle_horz, equation=eq_exit_angle_horz, short_name="exitangle_horz_rad", unit_symbol="rad", @@ -971,6 +973,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle_vert", scale=180.0 / numpy.pi, label=r"Vertical exit angle (deg)", + formula=formula_exit_angle_vert, equation=eq_exit_angle_vert, short_name="exitangle_vert", unit_symbol="deg", @@ -979,6 +982,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("exit_angle_horz", scale=180.0 / numpy.pi, label=r"Horizontal exit angle (deg)", + formula=formula_exit_angle_horz, equation=eq_exit_angle_horz, short_name="exitangle_horz", unit_symbol="deg", From b2ea9e84dc36c7eebbef1ee72ad05370e88df8bb Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 16:15:53 +0100 Subject: [PATCH 10/19] compatibility --- src/pyFAI/units.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 992471367..7d8f58555 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -751,6 +751,8 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula_z_rot_tangle = formula_z_rot_iangle formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} ** {formula_z_rot_tangle} + {formula_x_rot_tangle} ** {formula_x_rot_tangle}))" formula_exit_angle_horz = f"arctan2(({formula_x_rot_tangle}), ({formula_z_rot_tangle}))" +formula_exit_angle = formula_scattering_angle_vert +formula_exit_angle_horz = formula_scattering_angle_horz formula_qbeam_lab = f"2.0e-9/λ*π*(cos({formula_scattering_angle_vert})*cos({formula_scattering_angle_horz}) - 1)" formula_qhorz_lab = f"2.0e-9/λ*π*cos({formula_scattering_angle_vert})*sin({formula_scattering_angle_horz})" From 169eddf7d032b047f25846096a8841e539e1cef0 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 17:04:20 +0100 Subject: [PATCH 11/19] correct numexpr formulas --- src/pyFAI/units.py | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 7d8f58555..f5fc69c65 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -385,12 +385,12 @@ def eq_exit_angle_vert(x, y, z, wavelength=None, incident_angle=0.0, tilt_angle= :return: vertical exit angle in radians """ rot_incident_angle = numpy.array([[1,0,0], - [0,numpy.cos(incident_angle), numpy.sin(-incident_angle)], - [0, numpy.sin(incident_angle), numpy.cos(incident_angle)]], + [0,numpy.cos(incident_angle), numpy.sin(-incident_angle)], + [0, numpy.sin(incident_angle), numpy.cos(incident_angle)]], ) rot_tilt_angle = numpy.array([[numpy.cos(tilt_angle), numpy.sin(-tilt_angle), 0], - [numpy.sin(tilt_angle), numpy.cos(tilt_angle), 0], - [0, 0, 1]], + [numpy.sin(tilt_angle), numpy.cos(tilt_angle), 0], + [0, 0, 1]], ) rotated_xyz = numpy.tensordot(rot_incident_angle, numpy.stack((x,y,z)), axes=1) xp, yp, zp = numpy.tensordot(rot_tilt_angle, rotated_xyz, axes=1) @@ -741,16 +741,16 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula_qx = f"4.0e-9*π/λ*sin(arctan2(x, z)/2.0)" # TODO: wrong, fix me formula_qy = f"4.0e-9*π/λ*sin(arctan2(y, z)/2.0)" # TODO: wrong, fix me -formula_scattering_angle_vert = "arctan2(y, sqrt(z*z + x*x))" +formula_scattering_angle_vert = "arctan2(y, sqrt(z*z+x*x))" formula_scattering_angle_horz = "arctan2(x,z)" formula_x_rot_iangle = "x" -formula_y_rot_iangle = "y * cos(η) - z * sin(η)" -formula_z_rot_iangle = "y * sin(η) + z * cos(η)" -formula_x_rot_tangle = f"({formula_x_rot_iangle}) * cos(χ) - ({formula_y_rot_iangle}) * sin(χ)" -formula_y_rot_tangle = f"({formula_x_rot_iangle}) * sin(χ) + ({formula_y_rot_iangle}) * cos(χ)" +formula_y_rot_iangle = "(y*cos(η)-z*sin(η))" +formula_z_rot_iangle = "(y*sin(η)+z*cos(η))" +formula_x_rot_tangle = f"({formula_x_rot_iangle} * cos(χ) - {formula_y_rot_iangle} * sin(χ))" +formula_y_rot_tangle = f"({formula_x_rot_iangle} * sin(χ) + {formula_y_rot_iangle} * cos(χ))" formula_z_rot_tangle = formula_z_rot_iangle -formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} ** {formula_z_rot_tangle} + {formula_x_rot_tangle} ** {formula_x_rot_tangle}))" -formula_exit_angle_horz = f"arctan2(({formula_x_rot_tangle}), ({formula_z_rot_tangle}))" +formula_exit_angle_vert = f"arctan2({formula_y_rot_tangle}, sqrt({formula_z_rot_tangle} * {formula_z_rot_tangle} + {formula_x_rot_tangle} * {formula_x_rot_tangle}))" +formula_exit_angle_horz = f"arctan2({formula_x_rot_tangle}, ({formula_z_rot_tangle}))" formula_exit_angle = formula_scattering_angle_vert formula_exit_angle_horz = formula_scattering_angle_horz From 4a8beeb8c8e5037cca7079fc3703c08da85b57cb Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 19:24:17 +0100 Subject: [PATCH 12/19] polar angle unit --- src/pyFAI/units.py | 52 +++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 51 insertions(+), 1 deletion(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index f5fc69c65..00ce1d462 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -312,6 +312,18 @@ def register_azimuthal_unit(name, scale=1, label=None, equation=None, formula=No corner, delta, short_name, unit_symbol, positive, period) ANY_UNITS.update(AZIMUTHAL_UNITS) +def register_azimuthal_fiber_unit(name, scale=1, label=None, equation=None, formula=None, + incident_angle=0.0, tilt_angle=0.0, sample_orientation=1, + center=None, corner=None, delta=None, short_name=None, + unit_symbol=None, positive=False, period=None): + AZIMUTHAL_UNITS[name] = UnitFiber(name=name, scale=scale, label=label, + equation=equation, formula=formula, + incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation, + center=center, corner=corner, delta=delta, + short_name=short_name, unit_symbol=unit_symbol, + positive=positive, period=period, + ) + ANY_UNITS.update(AZIMUTHAL_UNITS) def eq_r(x, y, z=None, wavelength=None): """Calculates the radius in meter @@ -731,6 +743,22 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o eq_qip(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 + eq_qoop(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 ) + +def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): + """Calculates the polar angle from the vertical axis (fiber or thin-film main axis) + + :param x: horizontal position, towards the center of the ring, from sample position + :param y: vertical position, to the roof, from sample position + :param z: distance from sample along the beam + :param wavelength: in meter + :param incident_angle: tilting of the sample towards the beam (analog to rot2): in radians + :param tilt_angle: tilting of the sample orthogonal to the beam direction (analog to rot3): in radians + :param int sample_orientation: 1-8, orientation of the fiber axis according to EXIF orientation values (see def rotate_sample_orientation) + :return: component of the scattering vector in the plane YZ, in inverse nm + """ + qoop = eq_qoop(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) + qip = eq_qip(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) + return numpy.arctan2(qoop, qip) formula_r = "sqrt(x * x + y * y)" formula_2th = f"arctan2({formula_r}, z)" @@ -762,6 +790,8 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula_qvert_rot = f"-sin(χ)*({formula_qhorz_lab})-cos(χ)*sin(η)*({formula_qbeam_lab})+cos(χ)*cos(η)*({formula_qvert_lab})" formula_qip = f"sqrt(({formula_qbeam_rot})**2+({formula_qhorz_rot})**2)*((({formula_qhorz_rot} > 0) * 2) - 1)" formula_qoop = formula_qvert_rot +formula_qtot = f"sqrt(({formula_qip}) * ({formula_qip}) + ({formula_qoop}) * ({formula_qoop}))" +formula_chi_gi = f"arctan2(({formula_qip}), ({formula_qoop}))" register_radial_unit("r_mm", center="rArray", @@ -1038,6 +1068,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o register_radial_fiber_unit("qtot_nm^-1", scale=1.0, label=r"Scattering vector $q_{xyz}$ ($nm^{-1}$)", + formula=formula_qtot, equation=eq_q_total, short_name="q", unit_symbol="nm^{-1}", @@ -1116,6 +1147,20 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o formula=formula_chi, positive=False, period=360) +register_azimuthal_fiber_unit(name="chigi_rad", + scale=1.0, + label=r"Polar angle $\chi$ ($rad$)", + formula=formula_chi_gi, + equation=eq_chi_gi, + positive=False, + period=2.*pi) +register_azimuthal_fiber_unit(name="chigi_deg", + scale=180. / pi, + label=r"Polar angle $\chi$ ($^{o}$)", + formula=formula_chi_gi, + equation=eq_chi_gi, + positive=False, + period=360) AZIMUTHAL_UNITS["qx_nm^-1"] = RADIAL_UNITS["qx_nm^-1"] AZIMUTHAL_UNITS["qy_nm^-1"] = RADIAL_UNITS["qy_nm^-1"] @@ -1179,7 +1224,12 @@ def get_unit_fiber(name, incident_angle:float =0.0, tilt_angle:float =0.0, sampl :param float tilt angle: roll angle. Its rotation axis is orthogonal to the beam, the horizontal axis of the lab frame :param int sample_orientation: 1-8, orientation of the fiber axis according to EXIF orientation values (see def rotate_sample_orientation) """ - unit = copy.deepcopy(RADIAL_UNITS.get(name, None)) + if name in RADIAL_UNITS: + unit = copy.deepcopy(RADIAL_UNITS.get(name, None)) + elif name in AZIMUTHAL_UNITS: + unit = copy.deepcopy(AZIMUTHAL_UNITS.get(name, None)) + else: + unit = None if isinstance(unit, UnitFiber): unit.set_incident_angle(incident_angle) From aa1930ca567ed1986f275f3f9ae46969758b9574 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 17 Jan 2025 18:24:37 +0000 Subject: [PATCH 13/19] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- src/pyFAI/units.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 00ce1d462..72e842a73 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -316,8 +316,8 @@ def register_azimuthal_fiber_unit(name, scale=1, label=None, equation=None, form incident_angle=0.0, tilt_angle=0.0, sample_orientation=1, center=None, corner=None, delta=None, short_name=None, unit_symbol=None, positive=False, period=None): - AZIMUTHAL_UNITS[name] = UnitFiber(name=name, scale=scale, label=label, - equation=equation, formula=formula, + AZIMUTHAL_UNITS[name] = UnitFiber(name=name, scale=scale, label=label, + equation=equation, formula=formula, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation, center=center, corner=corner, delta=delta, short_name=short_name, unit_symbol=unit_symbol, @@ -743,7 +743,7 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o eq_qip(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 + eq_qoop(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 ) - + def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): """Calculates the polar angle from the vertical axis (fiber or thin-film main axis) From 92d3576b91c35119e7367335c0c1c1b504018a49 Mon Sep 17 00:00:00 2001 From: Edgar Date: Fri, 17 Jan 2025 19:43:00 +0100 Subject: [PATCH 14/19] fix polar --- src/pyFAI/units.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index 00ce1d462..278e42199 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -758,7 +758,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or """ qoop = eq_qoop(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) qip = eq_qip(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) - return numpy.arctan2(qoop, qip) + return numpy.arctan2(qip, qoop) formula_r = "sqrt(x * x + y * y)" formula_2th = f"arctan2({formula_r}, z)" From ce44e6c7b00156d6e408918407ff0cbafc17c644 Mon Sep 17 00:00:00 2001 From: edgar1993a Date: Mon, 20 Jan 2025 13:09:55 +0100 Subject: [PATCH 15/19] explicit _rad or _deg --- src/pyFAI/units.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index fdf28d293..da2572f98 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -966,7 +966,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or unit_symbol="nm^{-1}", positive=False) -register_radial_fiber_unit("scattering_angle_vert", +register_radial_fiber_unit("scattering_angle_vert_rad", scale=1.0, label=r"Vertical scattering angle (rad)", formula=formula_scattering_angle_vert, @@ -975,7 +975,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or unit_symbol="rad", positive=False) -register_radial_fiber_unit("scattering_angle_horz", +register_radial_fiber_unit("scattering_angle_horz_rad", scale=1.0, label=r"Horizontal scattering angle (rad)", formula=formula_scattering_angle_horz, @@ -1002,7 +1002,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or unit_symbol="rad", positive=False) -register_radial_fiber_unit("exit_angle_vert", +register_radial_fiber_unit("exit_angle_vert_deg", scale=180.0 / numpy.pi, label=r"Vertical exit angle (deg)", formula=formula_exit_angle_vert, @@ -1011,7 +1011,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or unit_symbol="deg", positive=False) -register_radial_fiber_unit("exit_angle_horz", +register_radial_fiber_unit("exit_angle_horz_deg", scale=180.0 / numpy.pi, label=r"Horizontal exit angle (deg)", formula=formula_exit_angle_horz, From 62ab80df38a99d7af2e8cf93e379c26a9cf14c81 Mon Sep 17 00:00:00 2001 From: edgar1993a Date: Thu, 23 Jan 2025 15:25:28 +0100 Subject: [PATCH 16/19] qtot fast --- src/pyFAI/units.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/src/pyFAI/units.py b/src/pyFAI/units.py index b061a6d10..0ffd1f670 100644 --- a/src/pyFAI/units.py +++ b/src/pyFAI/units.py @@ -739,10 +739,8 @@ def eq_q_total(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_o :param int sample_orientation: 1-8, orientation of the fiber axis according to EXIF orientation values (see def rotate_sample_orientation) :return: component of the scattering vector in the plane YZ, in inverse nm """ - return numpy.sqrt( - eq_qip(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 + - eq_qoop(x=x, y=y, z=z, wavelength=wavelength, incident_angle=incident_angle, tilt_angle=tilt_angle, sample_orientation=sample_orientation) ** 2 - ) + hpos, vpos = rotate_sample_orientation(x=x, y=y, sample_orientation=sample_orientation) + return 4.0e-9 * numpy.pi * numpy.sin(eq_2th(hpos, vpos, z) / 2.0) / wavelength def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_orientation=1): """Calculates the polar angle from the vertical axis (fiber or thin-film main axis) @@ -790,7 +788,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or formula_qvert_rot = f"-sin(χ)*({formula_qhorz_lab})-cos(χ)*sin(η)*({formula_qbeam_lab})+cos(χ)*cos(η)*({formula_qvert_lab})" formula_qip = f"sqrt(({formula_qbeam_rot})**2+({formula_qhorz_rot})**2)*((({formula_qhorz_rot} > 0) * 2) - 1)" formula_qoop = formula_qvert_rot -formula_qtot = f"sqrt(({formula_qip}) * ({formula_qip}) + ({formula_qoop}) * ({formula_qoop}))" +formula_qtot = formula_q formula_chi_gi = f"arctan2(({formula_qip}), ({formula_qoop}))" register_radial_unit("r_mm", @@ -1101,6 +1099,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or register_radial_fiber_unit("qip_A^-1", scale=0.1, label=r"Scattering vector $q_{IP}$ ($A^{-1}$)", + formula=formula_qip, equation=eq_qip, short_name="qip", unit_symbol="A^{-1}", @@ -1109,6 +1108,7 @@ def eq_chi_gi(x, y, z, wavelength, incident_angle=0.0, tilt_angle=0.0, sample_or register_radial_fiber_unit("qoop_A^-1", scale=0.1, label=r"Scattering vector $q_{OOP}$ ($A^{-1}$)", + formula=formula_qoop, equation=eq_qoop, short_name="qoop", unit_symbol="A^{-1}", From 5e554ec32ce62c339673d8c4916c7b88f5ecc02c Mon Sep 17 00:00:00 2001 From: edgar1993a Date: Thu, 23 Jan 2025 16:45:47 +0100 Subject: [PATCH 17/19] fix res.aximuthal --- src/pyFAI/integrator/fiber.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pyFAI/integrator/fiber.py b/src/pyFAI/integrator/fiber.py index d5ffa9658..f12760b24 100644 --- a/src/pyFAI/integrator/fiber.py +++ b/src/pyFAI/integrator/fiber.py @@ -276,7 +276,7 @@ def integrate_fiber(self, data, else: sum_variance = None sigma = None - result = Integrate1dResult(res.azimuthal * unit_scale, intensity, sigma) + result = Integrate1dResult(res.azimuthal, intensity, sigma) result._set_method_called("integrate_radial") result._set_unit(output_unit) result._set_sum_normalization(sum_normalization) From c046d6576fa1d09b2f77c40e2180bddb90770106 Mon Sep 17 00:00:00 2001 From: Edgar Date: Sun, 26 Jan 2025 22:21:02 +0100 Subject: [PATCH 18/19] script sandbox --- sandbox/test_gi_different_units.py | 190 +++++++++++++++++++++++++++++ 1 file changed, 190 insertions(+) create mode 100644 sandbox/test_gi_different_units.py diff --git a/sandbox/test_gi_different_units.py b/sandbox/test_gi_different_units.py new file mode 100644 index 000000000..57102224d --- /dev/null +++ b/sandbox/test_gi_different_units.py @@ -0,0 +1,190 @@ +from pyFAI.gui.jupyter import plot1d, plot2d, subplots +import matplotlib.pyplot as plt +from pyFAI.test.utilstest import UtilsTest +from pyFAI import load +import fabio +from pyFAI import detector_factory +from pyFAI.calibrant import get_calibrant +from pyFAI.integrator.fiber import FiberIntegrator + +if __name__ == "__main__": + + + # EXAMPLE WITH CLASSIC QIP-QOOP NM^-1 + fi = load(UtilsTest.getimage("LaB6_3.poni"), type_="pyFAI.integrator.fiber.FiberIntegrator") + data = fabio.open(UtilsTest.getimage("Y6.edf")).data + air = fabio.open(UtilsTest.getimage("air.edf")).data + data_clean = data - 0.03 * air + sample_orientation = 6 + + res2d_0 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, + ) + + res2d_patch_1 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, + ip_range=[-1,1], oop_range=[0,20], + ) + res2d_patch_2 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, + ip_range=[-17.5,17.5], oop_range=[5,6.5], + ) + + + res1d_0 = fi.integrate1d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, + ip_range=[-1,1], oop_range=[0,20], + npt_ip=100, npt_oop=500, + vertical_integration=True, + ) + res1d_1 = fi.integrate1d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, + ip_range=[-17.5,17.5], oop_range=[5,6.5], + npt_ip=500, npt_oop=100, + vertical_integration=False, + ) + + fig, ax = subplots(nrows=4, ncols=4) + plot2d(result=res2d_0, ax=ax[0,0]) + img = ax[0,0].get_images()[0] + img.set_cmap("viridis") + img.set_clim(20,200) + + plot2d(result=res2d_patch_1, ax=ax[0,1]) + plot2d(result=res2d_patch_2, ax=ax[0,1]) + plot2d(result=res2d_0, ax=ax[0,1]) + + ax[0,1].get_images()[0].set_cmap("viridis") + ax[0,1].get_images()[1].set_cmap("viridis") + ax[0,1].get_images()[2].set_cmap("gray") + ax[0,1].get_images()[2].set_alpha(0.5) + + plot1d(result=res1d_0, ax=ax[0,2]) + plot1d(result=res1d_1, ax=ax[0,3]) + + + + # EXAMPLE WITH POLAR UNITS RAD + det = detector_factory("Eiger_4M") + cal = get_calibrant("LaB6") + fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) + data = cal.fake_calibration_image(ai=fi) * 1000 + 100 + + sample_orientation = 1 + + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + ) + + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_nm^-1", unit_oop="chigi_rad", + ) + res1d_0 = fi.integrate1d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_nm^-1", unit_oop="chigi_rad", + ip_range=[0,40], oop_range=[-0.2,0.2], + npt_ip=1000, npt_oop=500, + vertical_integration=False, + ) + + res2d_patch = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_nm^-1", unit_oop="chigi_rad", + ip_range=[0,40], oop_range=[-0.2,0.2], + ) + + plot2d(result=res2d, ax=ax[1,0]) + plot2d(result=res2d_polar, ax=ax[1,1]) + plot2d(result=res2d_patch, ax=ax[1,2]) + plot2d(result=res2d_polar, ax=ax[1,2]) + + img = ax[1,0].get_images()[0] + img.set_cmap("viridis") + img = ax[1,1].get_images()[0] + img.set_cmap("viridis") + + ax[1,2].get_images()[0].set_cmap("viridis") + ax[1,2].get_images()[1].set_cmap("gray") + ax[1,2].get_images()[1].set_alpha(0.5) + plot1d(result=res1d_0, ax=ax[1,3]) + + + + # EXAMPLE WITH POLAR UNITS DEG-QA-1 + det = detector_factory("Eiger_4M") + cal = get_calibrant("LaB6") + fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) + data = cal.fake_calibration_image(ai=fi) * 1000 + 100 + + sample_orientation = 1 + + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + ) + + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_A^-1", unit_oop="chigi_deg", + ) + res1d_0 = fi.integrate1d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_A^-1", unit_oop="chigi_deg", + ip_range=[0,4], oop_range=[-10,10], + npt_ip=1000, npt_oop=500, + vertical_integration=False, + ) + + res2d_patch = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_ip="qtot_A^-1", unit_oop="chigi_deg", + ip_range=[0,4], oop_range=[-10,10], + ) + + plot2d(result=res2d, ax=ax[2,0]) + plot2d(result=res2d_polar, ax=ax[2,1]) + plot2d(result=res2d_patch, ax=ax[2,2]) + plot2d(result=res2d_polar, ax=ax[2,2]) + + img = ax[2,0].get_images()[0] + img.set_cmap("viridis") + img = ax[2,1].get_images()[0] + img.set_cmap("viridis") + + ax[2,2].get_images()[0].set_cmap("viridis") + ax[2,2].get_images()[1].set_cmap("gray") + ax[2,2].get_images()[1].set_alpha(0.5) + plot1d(result=res1d_0, ax=ax[2,3]) + + + + + # EXAMPLE WITH POLAR UNITS DEG-QA-1 NOW VERTICAL + det = detector_factory("Eiger_4M") + cal = get_calibrant("LaB6") + fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) + data = cal.fake_calibration_image(ai=fi) * 1000 + 100 + + sample_orientation = 1 + + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + ) + + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_oop="qtot_A^-1", unit_ip="chigi_deg", + ) + res1d_0 = fi.integrate1d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_oop="qtot_A^-1", unit_ip="chigi_deg", + oop_range=[0,4], ip_range=[-10,10], + npt_oop=1000, npt_ip=500, + vertical_integration=True, + ) + + res2d_patch = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, + unit_oop="qtot_A^-1", unit_ip="chigi_deg", + oop_range=[0,4], ip_range=[-10,10], + ) + + plot2d(result=res2d, ax=ax[3,0]) + plot2d(result=res2d_polar, ax=ax[3,1]) + plot2d(result=res2d_patch, ax=ax[3,2]) + plot2d(result=res2d_polar, ax=ax[3,2]) + + img = ax[3,0].get_images()[0] + img.set_cmap("viridis") + img = ax[3,1].get_images()[0] + img.set_cmap("viridis") + + ax[3,2].get_images()[0].set_cmap("viridis") + ax[3,2].get_images()[1].set_cmap("gray") + ax[3,2].get_images()[1].set_alpha(0.5) + plot1d(result=res1d_0, ax=ax[3,3]) + + plt.show() \ No newline at end of file From 51e0307545479144ad13c33f19e5425b1ae4d9d9 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Sun, 26 Jan 2025 21:21:59 +0000 Subject: [PATCH 19/19] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- sandbox/test_gi_different_units.py | 82 +++++++++++++++--------------- 1 file changed, 41 insertions(+), 41 deletions(-) diff --git a/sandbox/test_gi_different_units.py b/sandbox/test_gi_different_units.py index 57102224d..827117086 100644 --- a/sandbox/test_gi_different_units.py +++ b/sandbox/test_gi_different_units.py @@ -8,26 +8,26 @@ from pyFAI.integrator.fiber import FiberIntegrator if __name__ == "__main__": - - + + # EXAMPLE WITH CLASSIC QIP-QOOP NM^-1 fi = load(UtilsTest.getimage("LaB6_3.poni"), type_="pyFAI.integrator.fiber.FiberIntegrator") data = fabio.open(UtilsTest.getimage("Y6.edf")).data air = fabio.open(UtilsTest.getimage("air.edf")).data data_clean = data - 0.03 * air sample_orientation = 6 - + res2d_0 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, ) - + res2d_patch_1 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, ip_range=[-1,1], oop_range=[0,20], ) res2d_patch_2 = fi.integrate2d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, ip_range=[-17.5,17.5], oop_range=[5,6.5], ) - - + + res1d_0 = fi.integrate1d_grazing_incidence(data=data_clean, sample_orientation=sample_orientation, ip_range=[-1,1], oop_range=[0,20], npt_ip=100, npt_oop=500, @@ -37,39 +37,39 @@ ip_range=[-17.5,17.5], oop_range=[5,6.5], npt_ip=500, npt_oop=100, vertical_integration=False, - ) - + ) + fig, ax = subplots(nrows=4, ncols=4) plot2d(result=res2d_0, ax=ax[0,0]) img = ax[0,0].get_images()[0] img.set_cmap("viridis") img.set_clim(20,200) - + plot2d(result=res2d_patch_1, ax=ax[0,1]) plot2d(result=res2d_patch_2, ax=ax[0,1]) plot2d(result=res2d_0, ax=ax[0,1]) - + ax[0,1].get_images()[0].set_cmap("viridis") ax[0,1].get_images()[1].set_cmap("viridis") ax[0,1].get_images()[2].set_cmap("gray") ax[0,1].get_images()[2].set_alpha(0.5) - + plot1d(result=res1d_0, ax=ax[0,2]) plot1d(result=res1d_1, ax=ax[0,3]) - - - + + + # EXAMPLE WITH POLAR UNITS RAD det = detector_factory("Eiger_4M") cal = get_calibrant("LaB6") fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) data = cal.fake_calibration_image(ai=fi) * 1000 + 100 - + sample_orientation = 1 - + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, ) - + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, unit_ip="qtot_nm^-1", unit_oop="chigi_rad", ) @@ -84,35 +84,35 @@ unit_ip="qtot_nm^-1", unit_oop="chigi_rad", ip_range=[0,40], oop_range=[-0.2,0.2], ) - + plot2d(result=res2d, ax=ax[1,0]) plot2d(result=res2d_polar, ax=ax[1,1]) plot2d(result=res2d_patch, ax=ax[1,2]) plot2d(result=res2d_polar, ax=ax[1,2]) - + img = ax[1,0].get_images()[0] img.set_cmap("viridis") img = ax[1,1].get_images()[0] img.set_cmap("viridis") - + ax[1,2].get_images()[0].set_cmap("viridis") ax[1,2].get_images()[1].set_cmap("gray") ax[1,2].get_images()[1].set_alpha(0.5) plot1d(result=res1d_0, ax=ax[1,3]) - - - + + + # EXAMPLE WITH POLAR UNITS DEG-QA-1 det = detector_factory("Eiger_4M") cal = get_calibrant("LaB6") fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) data = cal.fake_calibration_image(ai=fi) * 1000 + 100 - + sample_orientation = 1 - + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, ) - + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, unit_ip="qtot_A^-1", unit_oop="chigi_deg", ) @@ -127,36 +127,36 @@ unit_ip="qtot_A^-1", unit_oop="chigi_deg", ip_range=[0,4], oop_range=[-10,10], ) - + plot2d(result=res2d, ax=ax[2,0]) plot2d(result=res2d_polar, ax=ax[2,1]) plot2d(result=res2d_patch, ax=ax[2,2]) plot2d(result=res2d_polar, ax=ax[2,2]) - + img = ax[2,0].get_images()[0] img.set_cmap("viridis") img = ax[2,1].get_images()[0] img.set_cmap("viridis") - + ax[2,2].get_images()[0].set_cmap("viridis") ax[2,2].get_images()[1].set_cmap("gray") ax[2,2].get_images()[1].set_alpha(0.5) plot1d(result=res1d_0, ax=ax[2,3]) - - - - + + + + # EXAMPLE WITH POLAR UNITS DEG-QA-1 NOW VERTICAL det = detector_factory("Eiger_4M") cal = get_calibrant("LaB6") fi = FiberIntegrator(detector=det, wavelength=1e-10, dist=0.1, poni1=0.05, poni2=0.05) data = cal.fake_calibration_image(ai=fi) * 1000 + 100 - + sample_orientation = 1 - + res2d= fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, ) - + res2d_polar = fi.integrate2d_grazing_incidence(data=data, sample_orientation=sample_orientation, unit_oop="qtot_A^-1", unit_ip="chigi_deg", ) @@ -171,20 +171,20 @@ unit_oop="qtot_A^-1", unit_ip="chigi_deg", oop_range=[0,4], ip_range=[-10,10], ) - + plot2d(result=res2d, ax=ax[3,0]) plot2d(result=res2d_polar, ax=ax[3,1]) plot2d(result=res2d_patch, ax=ax[3,2]) plot2d(result=res2d_polar, ax=ax[3,2]) - + img = ax[3,0].get_images()[0] img.set_cmap("viridis") img = ax[3,1].get_images()[0] img.set_cmap("viridis") - + ax[3,2].get_images()[0].set_cmap("viridis") ax[3,2].get_images()[1].set_cmap("gray") ax[3,2].get_images()[1].set_alpha(0.5) plot1d(result=res1d_0, ax=ax[3,3]) - - plt.show() \ No newline at end of file + + plt.show()