Finished math for forward/inverse asymmotron calculations.

To do this, I changed the α parameter from being a pseudo position to just be a soft signal.
spc-group · Oct 23, 2024 · 4826aa7 · 4826aa7
1 parent c31c5d6
commit 4826aa7
Show file tree

Hide file tree

Showing 2 changed files with 16 additions and 111 deletions.
diff --git a/src/haven/devices/analyzer.py b/src/haven/devices/analyzer.py
@@ -136,10 +136,10 @@ def __init__(
     d_spacing: Signal = Cpt(Signal, name="d_spacing")  # In Å
     rowland_diameter: Signal = Cpt(Signal, name="rowland_diameter")  # In mm
     wedge_angle: Signal = Cpt(Signal, name="wedge_angle")  # In radians
+    alpha: Signal = Cpt(Signal, name="alpha")
 
     # Pseudo axes
     energy: PseudoSingle = Cpt(PseudoSingle, name="energy", limits=(0, 1000))
-    alpha: PseudoSingle = Cpt(PseudoSingle, name="alpha", limits=(0, 180))
 
     # Real axes
     x: EpicsMotor = FCpt(EpicsMotor, "{x_motor_pv}", name="x")
@@ -149,11 +149,12 @@ def __init__(
     def forward(self, pseudo_pos):
         """Run a forward (pseudo -> real) calculation"""
         # Convert distance to microns and degrees to radians
-        energy, alpha = pseudo_pos.energy, pseudo_pos.alpha
+        energy = pseudo_pos.energy
         # Step 0: convert energy to bragg angle
         bragg = energy_to_bragg(energy, d=self.d_spacing.get())
         # Step 1: Convert energy params to geometry params
         D = self.rowland_diameter.get(use_monitor=True)
+        alpha = self.alpha.get()
         theta_M = bragg + alpha
         rho = D * np.sin(theta_M)
         # Step 2: Convert geometry params to motor positions
@@ -173,106 +174,13 @@ def inverse(self, real_pos):
         x = real_pos.x
         z = real_pos.z
         # Step 1: Convert motor positions to geometry parameters
-        theta_M = np.arctan((x + z * np.sin(beta)) / (z * np.cos(beta)))
+        theta_M = np.arctan2((x + z * np.sin(beta)) , (z * np.cos(beta)))
         rho = z * np.cos(beta) / np.cos(theta_M)
-        print(f"{x=}, {z=}, {beta=}, {theta_M=}, {rho}")
-        return self.PseudoPosition(energy=0, alpha=0)
-
-
-# Rewrite the following four equations so that they calculate D, theta and alpha based on inputting x, y, z1 and z.
-
-#         x = D * (sin(theta + alpha)) ** 2
-#         y = D * ((sin(theta + alpha)) ** 2 - (sin(theta - alpha)) ** 2)
-#         z1 = D * sin(theta - alpha) * cos(theta + alpha)
-#         z = z1 + D * sin(theta - alpha) * cos(theta - alpha)
-
-# To calculate D, theta, and alpha based on x, y, z1, and z, we can rearrange the given equations as follows:
-
-# D = x / (sin(theta + alpha)) ** 2
-
-# sin(theta - alpha) = z1 / (D * cos(theta + alpha))
-
-# cos(theta - alpha) = (z - z1) / (D * sin(theta - alpha))
-
-# Substitute the above two equations into the expression for y:
-
-# y = D * ((sin(theta + alpha)) ** 2 - (sin(theta - alpha)) ** 2)
-
-# y = D * ((sin(theta + alpha)) ** 2 - (z1 / (D * cos(theta + alpha))) ** 2)
-
-# y = D * ((sin(theta + alpha)) ** 2 - (z1 ** 2 / D ** 2) * (cos(theta + alpha)) ** (-2))
-
-# Solve for cos(theta + alpha) in the above equation:
-
-# cos(theta + alpha) = (z1 ** 2 / (D ** 2 * y + z1 ** 2)) ** 0.5
-
-# Substitute the value of cos(theta + alpha) into the equation for sin(theta - alpha):
-
-# sin(theta - alpha) = z1 / (D * (z1 ** 2 / (D ** 2 * y + z1 ** 2)) ** 0.5)
-
-# Simplify the above equation:
-
-# sin(theta - alpha) = (z1 ** 2 / (D ** 2 * y + z1 ** 2)) ** 0.5
-
-# Finally, substitute the values of cos(theta + alpha) and sin(theta - alpha) into the equation for D:
-
-# D = x / ((sin(theta + alpha)) ** 2)
-
-# D = x / ((1 - (z1 ** 2 / (D ** 2 * y + z1 ** 2))) ** 2)
-
-# This equation can be solved numerically to obtain the value of D. Once D is known, we can use the equations for cos(theta + alpha) and sin(theta - alpha) to calculate theta and alpha.
-
-
-# class LERIXSpectrometer(Device):
-#     rowland = Cpt(
-#         RowlandPositioner,
-#         x_motor_pv="vme_crate_ioc:m1",
-#         y_motor_pv="vme_crate_ioc:m2",
-#         z_motor_pv="vme_crate_ioc:m3",
-#         z1_motor_pv="vme_crate_ioc:m4",
-#         name="rowland",
-#     )
-
-
-# async def make_lerix_device(name: str, x_pv: str, y_pv: str, z_pv: str, z1_pv: str):
-#     dev = RowlandPositioner(
-#         name=name,
-#         x_motor_pv=x_pv,
-#         y_motor_pv=y_pv,
-#         z_motor_pv=z_pv,
-#         z1_motor_pv=z1_pv,
-#         labels={"lerix_spectrometers"},
-#     )
-#     pvs = ", ".join((x_pv, y_pv, z_pv, z1_pv))
-#     try:
-#         await await_for_connection(dev)
-#     except TimeoutError as exc:
-#         log.warning(f"Could not connect to LERIX spectrometer: {name} ({pvs})")
-#     else:
-#         log.info(f"Created area detector: {name} ({pvs})")
-#         return dev
-
-
-# def load_lerix_spectrometers(config=None):
-#     """Create devices for the LERIX spectrometer."""
-#     if config is None:
-#         config = load_config()
-#     # Create spectrometers
-#     devices = []
-#     for name, cfg in config.get("lerix", {}).items():
-#         rowland = cfg["rowland"]
-#         devices.append(
-#             make_device(
-#                 RowlandPositioner,
-#                 name=name,
-#                 x_motor_pv=rowland["x_motor_pv"],
-#                 y_motor_pv=rowland["y_motor_pv"],
-#                 z_motor_pv=rowland["z_motor_pv"],
-#                 z1_motor_pv=rowland["z1_motor_pv"],
-#                 labels={"lerix_spectromoters"},
-#             )
-#         )
-#     return devices
+        alpha = self.alpha.get()
+        bragg = theta_M - alpha
+        print(f"{x=}, {z=}, {beta=}, {theta_M=}, {rho=}, {bragg=}, {alpha=}")
+        energy = bragg_to_energy(bragg, d=self.d_spacing.get())
+        return self.PseudoPosition(energy=energy)
 
 
 # -----------------------------------------------------------------------------

diff --git a/src/haven/tests/test_analyzer.py b/src/haven/tests/test_analyzer.py
@@ -4,7 +4,7 @@
 import pytest
 from ophyd.sim import make_fake_device
 
-from haven.instrument import analyzer
+from haven.devices import analyzer
 
 um_per_mm = 1000
 
@@ -53,14 +53,10 @@ def test_wavelength_to_bragg(theta, d_spacing, wavelength):
     theta = np.radians(theta)
     d_spacing *= 1e-10
     wavelength *= 1e-10
-<<<<<<< Updated upstream
-    assert pytest.approx(analyzer.wavelength_to_bragg(wavelength, d=d_spacing), rel=0.001) == theta
-=======
     assert (
         pytest.approx(analyzer.wavelength_to_bragg(wavelength, d=d_spacing), rel=0.001)
         == theta
     )
->>>>>>> Stashed changes
 
 
 analyzer_values = [
@@ -87,26 +83,27 @@ def xtal(sim_registry):
 @pytest.mark.parametrize("bragg,alpha,beta,z,x", analyzer_values)
 def test_rowland_circle_forward(xtal, bragg, alpha, beta, x, z):
     xtal.wedge_angle.set(np.radians(beta)).wait()
+    xtal.alpha.set(np.radians(alpha)).wait()
     d = xtal.d_spacing.get()
     bragg = np.radians(bragg)
     energy = analyzer.bragg_to_energy(bragg, d=d)
     # Check the result is correct (convert cm -> m)
     expected = (x / 100, z / 100)
-    actual = xtal.forward(energy, np.radians(alpha))
+    actual = xtal.forward(energy)
     assert actual == pytest.approx(expected, rel=0.01)
 
 
 @pytest.mark.parametrize("bragg,alpha,beta,z,x", analyzer_values)
 def test_rowland_circle_inverse(xtal, bragg, alpha, beta, x, z):
     xtal.wedge_angle.set(np.radians(beta)).wait()
-    # Calculate the expected answer (convert cm -> m)
+    xtal.alpha.set(np.radians(alpha)).wait()
+    # Calculate the expected answer
     bragg = np.radians(bragg)
     d = xtal.d_spacing.get()
-    energy = analyzer.bragg_to_energy(bragg, d=d)
-    expected = (energy, alpha)
+    expected_energy = analyzer.bragg_to_energy(bragg, d=d)
     # Compare to the calculated inverse
     actual = xtal.inverse(x, z)
-    assert actual == pytest.approx(expected)
+    assert actual[0] == pytest.approx(expected_energy, abs=0.2)
 
 
 # -----------------------------------------------------------------------------