Skip to content

add functions to allow conversions between d and tth, and d and q #154

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.

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

Closed
wants to merge 13 commits into from
Closed
23 changes: 23 additions & 0 deletions news/xtype.rst
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
**Added:**

* functions to allow conversions between d and tth, and d and q.

**Changed:**

* <news item>

**Deprecated:**

* <news item>

**Removed:**

* <news item>

**Fixed:**

* <news item>

**Security:**

* <news item>
109 changes: 108 additions & 1 deletion src/diffpy/utils/scattering_objects/diffraction_objects.py
Original file line number Diff line number Diff line change
Expand Up @@ -283,7 +283,7 @@ def insert_scattering_quantity(
elif xtype.lower() in ANGLEQUANTITIES:
self.on_tth = [np.array(xarray), np.array(yarray)]
elif xtype.lower() in DQUANTITIES:
self.on_tth = [np.array(xarray), np.array(yarray)]
self.on_d = [np.array(xarray), np.array(yarray)]
self.set_all_arrays()

def q_to_tth(self):
Expand Down Expand Up @@ -362,18 +362,125 @@ def tth_to_q(self):
pre_factor = (4 * np.pi) / wavelength
return pre_factor * np.sin(two_theta / 2)

def q_to_d(self):
r"""
Helper function to convert q to d using :math:`d = \frac{2 \pi}{q}`

adds a small value (epsilon = 1e-10) to `q` where `q` is close to zero to avoid division by zero
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this seems like not a great way to handle this. Let's think about this.


Parameters
----------
q : array
An array of :math:`q` values

Returns
-------
d : array
An array of :math:`d` values in the inverse of the units of ``wavelength``
"""
epsilon = 1e-10
q = np.asarray(self.on_q[0])
q = np.where(np.abs(q) < epsilon, q + epsilon, q)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

let's think a bit how best to do this. Also, we may need to think if we want the d array to be so non-linear. How are we handling the interpolation when we go from q to tth? We can use that as inspiration.

return (2 * np.pi) / q

def d_to_q(self):
r"""
Helper function to convert d to q using :math:`q = \frac{2 \pi}{d}`

adds a small value (epsilon = 1e-10) to `d` where `d` is close to zero to avoid division by zero

Parameters
----------
d : array
An array of :math:`d` values

Returns
-------
q : array
An array of :math:`q` values in the inverse of the units of ``wavelength``
"""
epsilon = 1e-10
d = np.asarray(self.on_d[0])
d = np.where(np.abs(d) < epsilon, d + epsilon, d)
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Added d to a small number to avoid dividing by 0, which would cause q to become some arbitrary large number (as commented later in the test function). Not sure if this is a good way to address this problem - do we want to just raise error in these cases?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

pls see my comments

return (2 * np.pi) / d

def tth_to_d(self):
r"""
Helper function to convert two-theta to d

uses the formula .. math:: d = \frac{\lambda}{2 \sin\left(\frac{2\theta}{2}\right)},
and adds a small value (epsilon = 1e-10) to sin where sin is close to zero to avoid division by zero

Parameters
----------
two_theta : array
An array of :math:`2\theta` values in units of degrees

wavelength : float
Wavelength of the incoming x-rays

Returns
-------
d : array
An array of :math:`d` values in the inverse of the units
of ``wavelength``
"""
epsilon = 1e-10
two_theta = np.asarray(np.deg2rad(self.on_tth[0]))
wavelength = float(self.wavelength)
sin_two_theta = np.sin(two_theta / 2)
sin_two_theta = np.where(np.abs(sin_two_theta) < epsilon, sin_two_theta + epsilon, sin_two_theta)
return wavelength / (2 * sin_two_theta)

def d_to_tth(self):
r"""
Helper function to convert d to two-theta

uses the formula .. math:: 2\theta = 2 \arcsin\left(\frac{\lambda}{2d}\right),
and adds a small value (epsilon = 1e-10) to `d` where `d` is close to zero to avoid division by zero

Parameters
----------
d : array
An array of :math:`d` values

wavelength : float
Wavelength of the incoming x-rays

Returns
-------
two_theta : array
An array of :math:`2\theta` values in radians
"""
epsilon = 1e-10
d = np.asarray(self.on_d[0])
d = np.where(np.abs(d) < epsilon, d + epsilon, d)
wavelength = float(self.wavelength)
return np.rad2deg(2.0 * np.arcsin(wavelength / (2 * d)))

def set_all_arrays(self):
master_array, xtype = self._get_original_array()
if xtype == "q":
self.on_tth[0] = self.q_to_tth()
self.on_tth[1] = master_array[1]
self.on_d[0] = self.q_to_d()
self.on_d[1] = master_array[1]
if xtype == "tth":
self.on_q[0] = self.tth_to_q()
self.on_q[1] = master_array[1]
self.on_d[0] = self.tth_to_d()
self.on_d[1] = master_array[1]
if xtype == "d":
self.on_tth[0] = self.d_to_tth()
self.on_tth[1] = master_array[1]
self.on_q[0] = self.d_to_q()
self.on_q[1] = master_array[1]
self.tthmin = self.on_tth[0][0]
self.tthmax = self.on_tth[0][-1]
self.qmin = self.on_q[0][0]
self.qmax = self.on_q[0][-1]
self.dmin = self.on_d[0][0]
self.dmax = self.on_d[0][-1]

def _get_original_array(self):
if self.input_xtype in QQUANTITIES:
Expand Down
79 changes: 79 additions & 0 deletions tests/test_diffraction_objects.py
Original file line number Diff line number Diff line change
Expand Up @@ -231,6 +231,85 @@ def test_diffraction_objects_equality(inputs1, inputs2, expected):
assert (diffraction_object1 == diffraction_object2) == expected


def test_q_to_tth():
actual = Diffraction_object(wavelength=0.71)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

let's have a quick convo about what we want to test. Just one well behavied number is good, but any edge cases? IIn which case which ones? then maybe we can figure out a way to create the inputs once and reuse them to test all the transforms together to save time? That may also be a mistake.

setattr(actual, "on_q", [[0, 4.58087], [1, 1]])
actual_tth = actual.q_to_tth()
expected_tth = [0, 30]
assert np.allclose(actual_tth, expected_tth)


def test_tth_to_q():
actual = Diffraction_object(wavelength=0.71)
setattr(actual, "on_tth", [[0, 30], [1, 1]])
actual_q = actual.tth_to_q()
expected_q = [0, 4.58087]
assert np.allclose(actual_q, expected_q)


def test_q_to_d():
actual = Diffraction_object(wavelength=0.71)
setattr(actual, "on_q", [[0, 4.58087], [1, 1]])
actual_d = actual.q_to_d()
expected_d = [62831853071.8, 1.37161]
assert np.allclose(actual_d, expected_d)


def test_d_to_q():
actual = Diffraction_object(wavelength=0.71)
setattr(actual, "on_d", [[1e10, 1.37161], [1, 1]])
actual_q = actual.d_to_q()
expected_q = [0, 4.58087]
assert np.allclose(actual_q, expected_q)


def test_tth_to_d():
actual = Diffraction_object(wavelength=0.71)
setattr(actual, "on_tth", [[0, 30], [1, 1]])
actual_d = actual.tth_to_d()
expected_d = [3550000000, 1.37161]
assert np.allclose(actual_d, expected_d)


def test_d_to_tth():
actual = Diffraction_object(wavelength=0.71)
setattr(actual, "on_d", [[1e10, 1.37161], [1, 1]])
actual_tth = actual.d_to_tth()
expected_tth = [0, 30]
assert np.allclose(actual_tth, expected_tth)


params_array = [
(["q", "on_q", [4.58087, 8.84956], [1, 2]]),
(["tth", "on_tth", [30, 60], [1, 2]]),
(["d", "on_d", [1.37161, 0.71], [1, 2]]),
]


@pytest.mark.parametrize("inputs", params_array)
def test_set_all_arrays(inputs):
input_xtype, on_xtype, xarray, yarray = inputs
expected_values = {
"on_tth": [np.array([30, 60]), np.array([1, 2])],
"on_q": [np.array([4.58087, 8.84956]), np.array([1, 2])],
"on_d": [np.array([1.37161, 0.71]), np.array([1, 2])],
"tthmin": 30,
"tthmax": 60,
"qmin": 4.58087,
"qmax": 8.84956,
"dmin": 1.37161,
"dmax": 0.71,
}

actual = Diffraction_object(wavelength=0.71)
setattr(actual, "input_xtype", input_xtype)
setattr(actual, on_xtype, [xarray, yarray])
actual.set_all_arrays()
for attr, expected in expected_values.items():
actual_value = getattr(actual, attr)
assert np.allclose(actual_value, expected)


def test_dump(tmp_path, mocker):
x, y = np.linspace(0, 5, 6), np.linspace(0, 5, 6)
directory = Path(tmp_path)
Expand Down