From a95b0aa10e9dacd863a6f1db19ffd6968bc5d5bb Mon Sep 17 00:00:00 2001
From: Marco Cammarata <marcocamma@gmail.com>
Date: Tue, 17 Jan 2017 19:21:02 +0100
Subject: [PATCH] added angle conversion theta/q

---
 xray/utils.py | 24 ++++++++++++++++++++++++
 1 file changed, 24 insertions(+)

diff --git a/xray/utils.py b/xray/utils.py
index 4053f5c..cb21b1c 100644
--- a/xray/utils.py
+++ b/xray/utils.py
@@ -216,3 +216,27 @@ def reshapeToBroadcast(what,ref):
   shape  = [ref.shape[0],] + [1,]*(ref.ndim-1)
   return what.reshape(shape)
 
+def radToQ(theta,**kw):
+  """ theta is the scattering angle (theta_out - theta_in);
+      kw should be have either E or wavelength
+      it returns the scattering vector in the units of wavelength """
+  # Energy or wavelength should be in kw
+  assert "E" in kw or "wavelength" in kw
+  # but not both
+  assert not ("E" in kw and "wavelength" in kw)
+  if "E" in kw: kw["wavelength"] = 12.398/kw["E"]
+  return 4*np.pi/kw["wavelength"]*np.sin(theta)
+
+def degToQ(theta,**kw):
+  theta = theta/180.*np.pi
+  return radToQ(theta,**kw)
+degToQ.__doc__ = radToQ.__doc__
+
+def qToTheta(q,**kw):
+  """ Return scattering angle from q (given E or wavelength) """
+  # Energy or wavelength should be in kw
+  assert "E" in kw or "wavelength" in kw
+  # but not both
+  assert not ("E" in kw and "wavelength" in kw)
+  if "E" in kw: kw["wavelength"] = 12.398/kw["E"]
+  return np.arcsin(q*kw["wavelength"]/4/np.pi)