wip: digital twin

debye_bec/bec_widgets/widgets/digital_twin/calculate_positions.py

@@ -54,10 +54,10 @@ def calc_positions(cfg):
     # Bragg Angle
     # TODO Should the bragg angle be corrected for the symmetric bragg case?
     #      See raytracing script or here: bragg = np.asin(rm.ch / (2.*cfg['dSpacing']*cfg['energyCCM'])) - aCrystal.get_dtheta_symmetric_Bragg(cfg['energyCCM'])
-    if cfg['mo_mode'] == 'Monochromatic':
+    if cfg['mo1_mode'] == 'Monochromatic':
         # Add 2x CM pitch to the bragg angle
-        bragg =  ((2 * cfg['cm_pitch']) + cfg['mo_bragg']) / np.pi * 180
-    elif cfg['mo_mode'] == 'Pinkbeam':
+        bragg =  ((2 * cfg['cm_pitch']) + cfg['mo1_bragg']) / np.pi * 180
+    elif cfg['mo1_mode'] == 'Pinkbeam':
         # Align xtal surfaces parallel to beam
         bragg =  (2 * cfg['cm_pitch']) / np.pi * 180
     else:
@@ -65,13 +65,13 @@ def calc_positions(cfg):
     pos['mo1_bragg_angle'] = {'value': bragg} # Bragg angle in deg
 
     # TRY, Height
-    l = bl.mo1.xtalGap[0]/np.sin(cfg['mo_bragg'])
-    yhor = l*np.cos(2.*(cfg['mo_bragg']+cfg['cm_pitch']))
+    l = bl.mo1.xtalGap[0]/np.sin(cfg['mo1_bragg'])
+    yhor = l*np.cos(2.*(cfg['mo1_bragg']+cfg['cm_pitch']))
     yver = yhor*np.tan(2.*cfg['cm_pitch'])
 
-    if cfg['mo_mode'] == 'Monochromatic':
-        beamOffsetCCM = l*np.sin(2.*(cfg['mo_bragg']+cfg['cm_pitch']))-yver # Resultat ist korrekt!
-    elif cfg['mo_mode'] == 'Pinkbeam':
+    if cfg['mo1_mode'] == 'Monochromatic':
+        beamOffsetCCM = l*np.sin(2.*(cfg['mo1_bragg']+cfg['cm_pitch']))-yver # Resultat ist korrekt!
+    elif cfg['mo1_mode'] == 'Pinkbeam':
         beamOffsetCCM = 0
     else:
         raise Exception('Monochromator mode not supported')
@@ -87,13 +87,13 @@ def calc_positions(cfg):
     # logger.info(f'f = {f}')
     d = bl.mo1.heightOffset # xtal height offset, mm
     # logger.info(f'd = {d}')
-    c = d*csc(cfg['mo_bragg'])-f*cot(cfg['mo_bragg'])
+    c = d*csc(cfg['mo1_bragg'])-f*cot(cfg['mo1_bragg'])
     # logger.info(f'c = {c}')
 
     # Calculate height of center of rotation
-    b = np.sqrt(d**2*csc(cfg['mo_bragg'])**2-2*d*f*cot(cfg['mo_bragg'])*csc(cfg['mo_bragg'])+f**2*cot(cfg['mo_bragg'])**2+f**2)
+    b = np.sqrt(d**2*csc(cfg['mo1_bragg'])**2-2*d*f*cot(cfg['mo1_bragg'])*csc(cfg['mo1_bragg'])+f**2*cot(cfg['mo1_bragg'])**2+f**2)
     # logger.info(f'b = {b}')
-    h = np.cos(np.pi/2-np.arctan(f/c)-cfg['mo_bragg']-2*cfg['cm_pitch'])*b
+    h = np.cos(np.pi/2-np.arctan(f/c)-cfg['mo1_bragg']-2*cfg['cm_pitch'])*b
     # logger.info(f'h = {h}')
     h2 = ((bl.mo1.center[1] - bl.cm.center[1])-np.sqrt(b**2-h**2))*np.tan(2*cfg['cm_pitch'])
     # logger.info(f'mo1 = {bl.mo1.center[1]}')
@@ -103,18 +103,18 @@ def calc_positions(cfg):
     #TODO Mono height not exactly the same as in raytracing
     heightCCM1real = h + h2 # per design, the height should not change if the pitch of the CM is not changed!
     # heightCCM1real = heightCCM1real - 30 # Zero position of stage is at 1430 mm from ground.
-    if cfg['mo_mode'] == 'Monochromatic':
+    if cfg['mo1_mode'] == 'Monochromatic':
         pass
-    elif cfg['mo_mode'] == 'Pinkbeam':
+    elif cfg['mo1_mode'] == 'Pinkbeam':
         heightCCM1real = heightCCM1real - 13 # Move down to let beam pass between both crystal without touching copper cooler
     else:
         raise Exception('Monochromator mode not supported')
     pos['mo1_try'] = {'value': heightCCM1real}
 
     # TRX, Crystal selection
-    if cfg['mo_mode'] == 'Monochromatic':
+    if cfg['mo1_mode'] == 'Monochromatic':
         try:
-            xtal = cfg['mo_xtal'].translate(str.maketrans('', '', '()')) # Remove brackets from xtal name to conform with parameters
+            xtal = cfg['mo1_xtal'].translate(str.maketrans('', '', '()')) # Remove brackets from xtal name to conform with parameters
             index = bl.mo1.xtal.index(xtal)
         except:
             raise ValueError(f"Requested xtal {xtal} not found in parameters!")

debye_bec/bec_widgets/widgets/digital_twin/calculate_sideview.py

@@ -0,0 +1,42 @@
+import numpy as np
+import debye_bec.bec_widgets.widgets.digital_twin.x01da_parameters as bl
+
+def calc_sideview(cfg):
+
+    # Calculate height of beam after CM
+    height = 2 * bl.cm.center[1] * np.tan(cfg['v_acc'])
+
+    # beam height (Y=height, Z=along beam)
+    beam = {}
+    beam['Z'] = []
+    beam['Y'] = []
+    beam['Z'].append(0) # Source
+    beam['Y'].append(bl.sourceHeight)
+    beam['Z'].append(bl.cm.center[1]) # CM 
+    beam['Y'].append(bl.sourceHeight)
+    if cfg['mo1_mode'] in 'Monochromatic':
+        diag = bl.mo1.xtalGap[0]/np.sin(cfg['mo1_bragg']) # Calculations for Mono
+        dy = diag*np.sin(2*(cfg['cm_pitch']+cfg['mo1_bragg']))
+        dz = diag*np.cos(2*(cfg['cm_pitch']+cfg['mo1_bragg']))
+        beam['Z'].append(bl.mo1.center[1]-dz/2) # Mono 1.1
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.mo1.center[1]-dz/2-bl.cm.center[1]))
+        beam['Z'].append(bl.mo1.center[1]+dz/2) # Mono 1.2
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.mo1.center[1]-dz/2-bl.cm.center[1])+dy)
+        beam['Z'].append(bl.fm.center[1]) # FM
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.fm.center[1]-bl.cm.center[1]-dz)+dy)
+        beam['Z'].append(cfg['smpl']) # Experiment
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.fm.center[1]-bl.cm.center[1]-dz)+dy+np.tan(2*(cfg['cm_pitch']-cfg['fm_pitch']))*(cfg['smpl']-bl.fm.center[1]))
+    elif cfg['mo1_mode'] == 'Pinkbeam':
+        beam['Z'].append(bl.fm.center[1]) # FM
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.fm.center[1]-bl.cm.center[1]))
+        beam['Z'].append(cfg['smpl']) # Experiment
+        beam['Y'].append(bl.sourceHeight+np.tan(2*cfg['cm_pitch'])*(bl.fm.center[1]-bl.cm.center[1])+np.tan(2*(cfg['cm_pitch']-cfg['fm_pitch']))*(cfg['smpl']-bl.fm.center[1]))
+
+    dy_fm_ex = beam['Y'][-1] - beam['Y'][-2]
+    dz_fm_ex = beam['Z'][-1] - beam['Z'][-2]
+    dz_fm_win = bl.ehWindow.center[1] - beam['Z'][-2]
+    h_at_win = beam['Y'][-2] + dy_fm_ex / dz_fm_ex * dz_fm_win
+
+    beam['heightWindow'] = h_at_win
+
+    return beam

debye_bec/bec_widgets/widgets/digital_twin/digital_twin.py

@@ -8,7 +8,7 @@ from qtpy.QtWidgets import (
     QApplication, QLayout
 )
 # pylint: disable=E0611
-from qtpy.QtCore import Qt
+from qtpy.QtCore import Qt, QTimer
 from qtpy.QtGui import QColor, QFont
 import pyqtgraph as pg
 
@@ -20,6 +20,7 @@ from bec_widgets.utils.error_popups import SafeSlot
 from debye_bec.bec_widgets.widgets.qt_widgets import InputNumberField, ComboBox, Group, NumberIndicator
 
 from debye_bec.bec_widgets.widgets.digital_twin.calculate_positions import calc_positions
+from debye_bec.bec_widgets.widgets.digital_twin.calculate_sideview import calc_sideview
 
 import debye_bec.bec_widgets.widgets.digital_twin.x01da_parameters as bl
 
@@ -79,29 +80,76 @@ class DigitalTwin(BECWidget, QWidget):
         self.calc_bragg_angle()
         self.calc_ideal_fm_pitch()
         self.calc_crit_angle()
+        self.calc_assistant_sideview()
+        self.calc_reality_sideview()
 
-    @SafeSlot()
-    def calc_positions(self, *args):
-        # logger.info(f'Got field {field} and number {qt_obj} and number {number} and args {args}')
+        # Timer: update plot every 1 second
+        self._timer = QTimer(self)
+        self._timer.setInterval(1000)
+        self._timer.timeout.connect(self.calc_reality_sideview)
+        # TODO: Check if I need to stop the timer if the widget is closed?
+        self._timer.start()
 
+    def get_assistant_config(self):
         config = { # Config in SI units!
                 'energy' : self.input.energy.value(),
                 'h_acc' : self.input.sldi_hacc.value() * 1e-3,
                 'v_acc' : self.input.sldi_vacc.value() * 1e-3,
                 'cm_pitch' : -self.input.cm_pitch.value() * 1e-3,
                 'cm_stripe' : self.input.cm_stripe.currentText(),
-                'mo_mode' : self.input.mo1_mode.currentText(),
-                'mo_xtal' : self.input.mo1_xtal.currentText(),
-                'mo_bragg' : self.bragg_angle,
+                'mo1_mode' : self.input.mo1_mode.currentText(),
+                'mo1_xtal' : self.input.mo1_xtal.currentText(),
+                'mo1_bragg' : self.bragg_angle,
                 'fm_pitch' : -self.input.fm_pitch.value() * 1e-3,
                 'fm_stripe' : self.input.fm_stripe.currentText(),
                 'fm_gain_height' : 1,
                 'smpl' : self.input.smpl.value(),
             }
+        # logger.info(f'Config created: {config}')
+        return config
+    
+    # TODO Needs to run in a loop in a separate thread due to the long time it takes to get the values from self.dev...
+    def get_reality_config(self):
+        if abs(self.dev.mo1_trx.position) > 5:
+            mo1_mode = 'Monochromatic'
+        else:
+            mo1_mode = 'Pinkbeam'
+        # TODO: stripe detection, mo1_bragg and acceptance
+        config = { # Config in SI units!
+                'energy' : self.dev.mo1_bragg.position,
+                'h_acc' : self.input.sldi_hacc.value() * 1e-3,
+                'v_acc' : self.input.sldi_vacc.value() * 1e-3,
+                'cm_pitch' : -self.dev.cm_rotx.position * 1e-3,
+                'cm_stripe' : self.input.cm_stripe.currentText(),
+                'mo1_mode' : mo1_mode,
+                'mo1_xtal' : self.dev.mo1_bragg.crystal.current_xtal_string.get(),
+                'mo1_bragg' : self.dev.mo1_bragg.angle.get(),
+                'fm_pitch' : -self.dev.fm_rotx.position * 1e-3,
+                'fm_stripe' : self.input.fm_stripe.currentText(),
+                'fm_gain_height' : 1,
+                'smpl' : self.dev.ot_es1_trz.position,
+            }
+        # logger.info(f'Config created: {config}')
+        return config
 
-        logger.info(f'Config created: {config}')
-        out = calc_positions(config)
-        logger.info(f'Got positions: {out}')
+    @SafeSlot()
+    def calc_assistant_sideview(self, *args):
+        beam = calc_sideview(self.get_assistant_config())
+        self.plot_widget.data['assistant'][0] = beam['Z']
+        self.plot_widget.data['assistant'][1] = beam['Y']
+        self.plot_widget.update_curves()
+
+    @SafeSlot()
+    def calc_reality_sideview(self):
+        logger.info('Update reality plot')
+        beam = calc_sideview(self.get_reality_config())
+        self.plot_widget.data['reality'][0] = beam['Z']
+        self.plot_widget.data['reality'][1] = beam['Y']
+        self.plot_widget.update_curves()
+
+    @SafeSlot()
+    def calc_positions(self, *args):
+        out = calc_positions(self.get_assistant_config())
 
         self.positions.sldi_gapx.setValue(out['sldi_gapx']['value'])
         self.positions.sldi_gapy.setValue(out['sldi_gapy']['value'])
@@ -124,6 +172,9 @@ class DigitalTwin(BECWidget, QWidget):
         self.positions.ot_rotx.setValue(out['ot_rotx']['value'])
         self.positions.ot_es1_trz.setValue(out['ot_es1_trz']['value'])
 
+        # TODO move to somewhere else!
+        self.calc_assistant_sideview()
+
     @SafeSlot()
     def calc_bragg_angle(self, *args):
         """