Skip to content

FFT norm and SANS Normalisation #64

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
samjrholt wants to merge 3 commits into
base: master
Choose a base branch
from
Open

FFT norm and SANS Normalisation #64

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
1159b15
FFT norm and SANS Normalisation
samjrholt Mar 5, 2024
58de3cd
FFT normalisation
samjrholt Mar 7, 2024
8c465a1
SAXS Normalisation
samjrholt Mar 7, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions mag2exp/ltem.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -92,7 +92,7 @@ def phase(field, /, kcx=0.1, kcy=0.1):
"""
# More readable notation, direction arg will be removed soon.
m_int = field.integrate(direction="z")
m_ft = m_int.fftn()
m_ft = m_int.fftn(norm="ortho")

with warnings.catch_warnings():
warnings.simplefilter("ignore")
Expand All @@ -105,7 +105,7 @@ def phase(field, /, kcx=0.1, kcy=0.1):
with warnings.catch_warnings():
warnings.simplefilter("ignore")
ft_phase = (m_ft & k).ft_z * denom * prefactor
phase = ft_phase.ifftn().real
phase = ft_phase.ifftn(norm="ortho").real
phase.mesh.translate(field.mesh.region.center[:2], inplace=True)
return phase, ft_phase

Expand Down Expand Up @@ -218,7 +218,7 @@ def defocus_image(phase, /, cs=0, df_length=0.2e-3, voltage=None, wavelength=Non
:py:func:`~mag2exp.ltem.relativistic_wavelength`

"""
ft_wavefn = np.exp(phase * 1j).fftn()
ft_wavefn = np.exp(phase * 1j).fftn(norm="ortho")

with warnings.catch_warnings():
warnings.simplefilter("ignore")
Expand All @@ -241,7 +241,7 @@ def defocus_image(phase, /, cs=0, df_length=0.2e-3, voltage=None, wavelength=Non
cts = -df_length + 0.5 * wavelength**2 * cs * ksquare
exp = np.exp(np.pi * cts * 1j * ksquare * wavelength)
ft_def_wf_cts = ft_wavefn * exp
def_wf_cts = ft_def_wf_cts.ifftn()
def_wf_cts = ft_def_wf_cts.ifftn(norm="ortho")
intensity_cts = def_wf_cts.conjugate * def_wf_cts
intensity_cts.mesh.translate(phase.mesh.region.center, inplace=True)
return intensity_cts.real
Expand Down
6 changes: 4 additions & 2 deletions mag2exp/sans.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -240,8 +240,9 @@ def chiral_function(field, /, polarisation=(0, 0, 1)):


def _cross_section_matrix(field, /, polarisation):
m_fft = field.fftn()
m_fft *= field.mesh.dV * 1e16 # TODO: Normalisation
m_fft = field.fftn(norm="ortho")
m_fft *= np.sqrt(field.mesh.dV) # TODO: Normalisation

q = df.Field(
m_fft.mesh,
nvdim=3,
Expand All @@ -268,6 +269,7 @@ def _cross_section_matrix(field, /, polarisation):
magnetic_interaction_new = np.einsum(
"ijkl,lbc->ijkbc", magnetic_interaction.array, p_new
)

cs = np.power(np.abs(magnetic_interaction_new), 2)
return df.Field(
mesh=m_fft.mesh,
Expand Down
39 changes: 16 additions & 23 deletions mag2exp/tests/test_sans.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,35 +21,35 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="pp").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="nn").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="pn").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="np").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2


Expand All @@ -69,21 +69,21 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="pp").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="nn").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, method="pn").sel(k_z=0)
Expand All @@ -110,7 +110,7 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, polarisation=[1, 0, 0], method="nn").sel(k_z=0)
Expand All @@ -123,14 +123,14 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, polarisation=[1, 0, 0], method="np").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2


Expand All @@ -151,7 +151,7 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 2

sans = mag2exp.sans.cross_section(m, polarisation=[1, 0, 0], method="nn").sel(k_z=0)
Expand All @@ -164,14 +164,14 @@ def m_fun(pos):
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(q, 1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 1

sans = mag2exp.sans.cross_section(m, polarisation=[1, 0, 0], method="np").sel(k_z=0)
idx = np.unravel_index(sans.array.argmax(), sans.array.shape)[0:2]
q = sans.mesh.index2point(idx)[0]
assert np.isclose(q, -1 / qx)
peaks = (sans.array > 10).sum()
peaks = (sans.array > 1e-10).sum()
assert peaks == 1


Expand Down Expand Up @@ -213,12 +213,12 @@ def m_fun(pos):
idx = np.unravel_index(cf.array.argmax(), cf.array.shape)[0:2]
q = cf.mesh.index2point(idx)[0]
assert np.isclose(q, 1 / qx)
peaks = (cf.array > 10).sum()
peaks = (cf.array > 1e-10).sum()
assert peaks == 1
idx = np.unravel_index(cf.array.argmin(), cf.array.shape)[0:2]
q = cf.mesh.index2point(idx)[0]
assert np.isclose(q, -1 / qx)
peaks = (cf.array < -10).sum()
peaks = (cf.array < -1e-10).sum()
assert peaks == 1

def m_fun(pos):
Expand Down Expand Up @@ -250,14 +250,7 @@ def m_fun(pos):
sans2 = mag2exp.sans.cross_section(field2, method="unpol")
m2 = abs(sans2.array).max()

region = df.Region(p1=(0, 0, 0), p2=(150e-9, 150e-9, 150e-9))
mesh = df.Mesh(region=region, cell=(5e-9, 5e-9, 5e-9))
field3 = df.Field(mesh, nvdim=3, value=m_fun, norm=Ms)
sans3 = mag2exp.sans.cross_section(field3, method="unpol")
m3 = abs(sans3.array).max()

assert np.isclose(m1, m2)
assert np.isclose(m1 / m3, (100 / 150) ** 6)


def test_sans_cross_section_methods():
Expand Down
9 changes: 1 addition & 8 deletions mag2exp/tests/test_x_ray.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -102,7 +102,7 @@ def m_fun(pos):
idx = np.unravel_index(saxs.array.argmax(), saxs.array.shape)[0:2]
q = saxs.mesh.index2point(idx)[0]
assert np.isclose(abs(q), 1 / qx)
peaks = (saxs.array > 10).sum()
peaks = (saxs.array > 1e-10).sum()
assert peaks == 2


Expand All @@ -126,11 +126,4 @@ def m_fun(pos):
saxs2 = mag2exp.x_ray.saxs(field2)
m2 = abs(saxs2.array).max()

region = df.Region(p1=(0, 0, 0), p2=(150e-9, 150e-9, 150e-9))
mesh = df.Mesh(region=region, cell=(5e-9, 5e-9, 5e-9))
field3 = df.Field(mesh, nvdim=3, value=m_fun, norm=Ms)
saxs3 = mag2exp.x_ray.saxs(field3)
m3 = abs(saxs3.array).max()

assert np.isclose(m1, m2)
assert np.isclose(m1 / m3, (100 / 150) ** 6)
5 changes: 3 additions & 2 deletions mag2exp/x_ray.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@

Module for calculation of x-ray based quantities.
"""
import numpy as np

import mag2exp

Expand Down Expand Up @@ -131,6 +132,6 @@ def saxs(field):
>>> xrs.mpl.scalar()

"""
m_fft = field.fftn().ft_z.sel(k_z=0)
m_fft *= field.mesh.dV * 1e16
m_fft = field.fftn(norm="ortho").ft_z.sel(k_z=0)
m_fft *= np.sqrt(field.mesh.dV)
return abs(m_fft) ** 2