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

news/xtype.rst

@@ -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>

src/diffpy/utils/scattering_objects/diffraction_objects.py

@@ -10,6 +10,8 @@
 DQUANTITIES = ["d", "dspace"]
 XQUANTITIES = ANGLEQUANTITIES + DQUANTITIES + QQUANTITIES
 XUNITS = ["degrees", "radians", "rad", "deg", "inv_angs", "inv_nm", "nm-1", "A-1"]
+QMAX = 40
+DMAX = 100
 
 x_grid_emsg = (
  "objects are not on the same x-grid. You may add them using the self.add method "
@@ -283,7 +285,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):
@@ -305,7 +307,7 @@ def q_to_tth(self):
  Parameters
  ----------
  q : array
-An array of :math:`q` values
+The array of :math:`q` values
 
  wavelength : float
  Wavelength of the incoming x-rays
@@ -315,12 +317,16 @@ def q_to_tth(self):
  Returns
  -------
  two_theta : array
-An array of :math:`2\theta` values in radians
+The array of :math:`2\theta` values in radians
  """
  q = self.on_q[0]
  q = np.asarray(q)
  wavelength = float(self.wavelength)
  pre_factor = wavelength / (4 * np.pi)
+ if np.any(np.abs(q * pre_factor) > 1):
+ raise ValueError(
+ "Invalid input for arcsin: some values exceed the range [-1, 1]. " "Check wavelength or q values."
+ )
  return np.rad2deg(2.0 * np.arcsin(q * pre_factor))
 
  def tth_to_q(self):
@@ -344,7 +350,7 @@ def tth_to_q(self):
  Parameters
  ----------
  two_theta : array
-An array of :math:`2\theta` values in units of degrees
+The array of :math:`2\theta` values in units of degrees
 
  wavelength : float
  Wavelength of the incoming x-rays
@@ -354,26 +360,134 @@ def tth_to_q(self):
  Returns
  -------
  q : array
-An array of :math:`q` values in the inverse of the units
+The array of :math:`q` values in the inverse of the units
  of ``wavelength``
  """
  two_theta = np.asarray(np.deg2rad(self.on_tth[0]))
+ if np.any(two_theta > np.pi):
+ raise ValueError("Two theta exceeds 180 degrees.")
  wavelength = float(self.wavelength)
  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}`, set dmax to DMAX=100
+
+ Parameters
+ ----------
+ q : array
+ The array of :math:`q` values
+
+ Returns
+ -------
+ d : array
+ The array of :math:`d` values in the inverse of the units of ``wavelength``
+ """
+ q = np.asarray(self.on_q[0])
+ d = np.where(q != 0, (2 * np.pi) / q, np.inf)
+ d = np.minimum(d, DMAX)
+ return d[::-1]
+
+ def d_to_q(self):
+ r"""
+ Helper function to convert d to q using :math:`q = \frac{2 \pi}{d}`, set qmax to QMAX=40
+
+ Parameters
+ ----------
+ d : array
+ The array of :math:`d` values
+
+ Returns
+ -------
+ q : array
+ The array of :math:`q` values in the inverse of the units of ``wavelength``
+ """
+ d = np.asarray(self.on_d[0])
+ q = np.where(d != 0, (2 * np.pi) / d, np.inf)
+ q = np.minimum(q, QMAX)
+ return q[::-1]
+
+ 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)}, set dmax to DMAX=100
+
+ Parameters
+ ----------
+ two_theta : array
+ The array of :math:`2\theta` values in units of degrees
+
+ wavelength : float
+ Wavelength of the incoming x-rays
+
+ Returns
+ -------
+ d : array
+ The array of :math:`d` values in the inverse of the units
+ of ``wavelength``
+ """
+ two_theta = np.asarray(np.deg2rad(self.on_tth[0]))
+ if np.any(two_theta > np.pi):
+ raise ValueError("Two theta exceeds 180 degrees.")
+ wavelength = float(self.wavelength)
+ sin_two_theta = np.sin(two_theta / 2)
+ d = np.where(sin_two_theta != 0, wavelength / (2 * sin_two_theta), np.inf)
+ d = np.minimum(d, DMAX)
+ return d[::-1]
+
+ 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), set tth to 180 when d=0
+
+ Parameters
+ ----------
+ d : array
+ The array of :math:`d` values
+
+ wavelength : float
+ Wavelength of the incoming x-rays
+
+ Returns
+ -------
+ two_theta : array
+ The array of :math:`2\theta` values in radians
+ """
+ d = np.asarray(self.on_d[0])
+ wavelength = float(self.wavelength)
+ pre_factor = wavelength / (2 * d)
+ if np.any(np.abs(pre_factor) > 1):
+ raise ValueError(
+ "Invalid input for arcsin: some values exceed the range [-1, 1]. " "Check wavelength or d values."
+ )
+ return np.rad2deg(2.0 * np.arcsin(pre_factor))[::-1]
+
  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:

tests/test_diffraction_objects.py

@@ -231,6 +231,128 @@ def test_diffraction_objects_equality(inputs1, inputs2, expected):
  assert (diffraction_object1 == diffraction_object2) == expected
 
 
+def test_q_to_tth():
+ # valid q values including edge cases when q=0 or 1, expected tth values are 2 * arcsin(q)
+ actual = Diffraction_object(wavelength=4 * np.pi)
+ setattr(actual, "on_q", [[0, 0.2, 0.4, 0.6, 0.8, 1], [1, 2, 3, 4, 5, 6]])
+ actual_tth = actual.q_to_tth()
+ expected_tth = [0, 23.07392, 47.15636, 73.73980, 106.26020, 180]
+ assert np.allclose(actual_tth, expected_tth)
+
+
+def test_q_to_tth_bad():
+ # invalid q values when arcsin value is not in the range of [-1, 1]
+ actual = Diffraction_object(wavelength=4 * np.pi)
+ setattr(actual, "on_q", [[0.6, 0.8, 1, 1.2], [1, 2, 3, 4]])
+ with pytest.raises(ValueError):
+ actual.q_to_tth()
+
+
+def test_tth_to_q():
+ # valid tth values including edge cases when tth=0 or 180 degree,
+ # expected q vales are sin15, sin30, sin45, sin60, sin90
+ actual = Diffraction_object(wavelength=4 * np.pi)
+ setattr(actual, "on_tth", [[0, 30, 60, 90, 120, 180], [1, 2, 3, 4, 5, 6]])
+ actual_q = actual.tth_to_q()
+ expected_q = [0, 0.258819, 0.5, 0.707107, 0.866025, 1]
+ assert np.allclose(actual_q, expected_q)
+
+
+def test_tth_to_q_bad():
+ # invalid tth values > 180 degree
+ actual = Diffraction_object(wavelength=4 * np.pi)
+ setattr(actual, "on_tth", [[0, 30, 60, 90, 120, 181], [1, 2, 3, 4, 5, 6]])
+ with pytest.raises(ValueError):
+ actual.tth_to_q()
+
+
+def test_q_to_d():
+ actual = Diffraction_object()
+ setattr(actual, "on_q", [[0, np.pi, 2 * np.pi, 3 * np.pi, 4 * np.pi, 5 * np.pi], [1, 2, 3, 4, 5, 6]])
+ actual_d = actual.q_to_d()
+ # expected d values are DMAX=100, 2/1, 2/2, 2/3, 2/4, 2/5, in reverse order
+ expected_d = [0.4, 0.5, 0.66667, 1, 2, 100]
+ assert np.allclose(actual_d, expected_d)
+
+
+def test_d_to_q():
+ actual = Diffraction_object()
+ setattr(actual, "on_d", [[0, np.pi, 2 * np.pi, 3 * np.pi, 4 * np.pi, 5 * np.pi], [1, 2, 3, 4, 5, 6]])
+ actual_q = actual.d_to_q()
+ # expected q values are QMAX=40, 2/1, 2/2, 2/3, 2/4, 2/5, in reverse order
+ expected_q = [0.4, 0.5, 0.66667, 1, 2, 40]
+ assert np.allclose(actual_q, expected_q)
+
+
+def test_tth_to_d():
+ # valid tth values including edge cases when tth=0 or 180 degree,
+ # expected d values are DMAX=100, 1/sin15, 1/sin30, 1/sin45, 1/sin60, 1/sin90, in reverse order
+ actual = Diffraction_object(wavelength=2)
+ setattr(actual, "on_tth", [[0, 30, 60, 90, 120, 180], [1, 2, 3, 4, 5, 6]])
+ actual_d = actual.tth_to_d()
+ expected_d = [1, 1.1547, 1.41421, 2, 3.8637, 100]
+ assert np.allclose(actual_d, expected_d)
+
+
+def test_tth_to_d_bad():
+ # invalid tth values > 180 degree
+ actual = Diffraction_object(wavelength=2)
+ setattr(actual, "on_tth", [[0, 30, 60, 90, 120, 181], [1, 2, 3, 4, 5, 6]])
+ with pytest.raises(ValueError):
+ actual.tth_to_d()
+
+
+def test_d_to_tth():
+ actual = Diffraction_object(wavelength=2)
+ setattr(actual, "on_d", [[0, 2, 4, 6, 8, 100, 200], [1, 2, 3, 4, 5, 6, 7]])
+ actual_tth = actual.d_to_tth()
+ # expected tth values are 2*arcsin(1/d), in reverse order
+ # when d is 0, we set tth to 180
+ # we allow d to exceed DMAX=100 for user inputs
+ expected_tth = [0.57296, 1.14593, 14.36151, 19.18814, 28.95502, 60, 180]
+ assert np.allclose(actual_tth, expected_tth)
+
+
+def test_d_to_tth_bad():
+ # invalid d values when arcsin value is not in the range of [-1, 1]
+ actual = Diffraction_object(wavelength=2)
+ setattr(actual, "on_d", [[0, 10, 20, 30], [1, 2, 3, 4]])
+ with pytest.raises(ValueError):
+ actual.d_to_tth()
+
+
+# Inputs cover a good valid range and provide the same expected values, ensuring conversion precision.
+params_array = [
+ (["tth", "on_tth", [30, 60, 90, 120, 150], [1, 2, 3, 4, 5]]),
+ (["q", "on_q", [1.626208, 3.141593, 4.442883, 5.441398, 6.069091], [1, 2, 3, 4, 5]]),
+ (["d", "on_d", [1.035276, 1.154701, 1.414214, 2, 3.863703], [1, 2, 3, 4, 5]]),
+]
+
+
+@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, 90, 120, 150]), np.array([1, 2, 3, 4, 5])],
+ "on_q": [np.array([1.626208, 3.141593, 4.442883, 5.441398, 6.069091]), np.array([1, 2, 3, 4, 5])],
+ "on_d": [np.array([1.035276, 1.154701, 1.414214, 2, 3.863703]), np.array([1, 2, 3, 4, 5])],
+ "tthmin": 30,
+ "tthmax": 150,
+ "qmin": 1.626208,
+ "qmax": 6.069091,
+ "dmin": 1.035276,
+ "dmax": 3.863703,
+ }
+
+ actual = Diffraction_object(wavelength=2)
+ 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)