add workbench autofocus.py

deltatau.py

@@ -1,43 +0,0 @@
-import sys,os,socket
-sys.path.insert(0, os.path.expanduser('..'))
-hostname=socket.gethostname()
-if hostname=='ganymede':
-  sys.path.insert(0, os.path.expanduser('~/Documents/prj/SwissFEL/PBTools/'))
-else:
-  #sys.path.insert(0, os.path.expanduser('/tmp/zamofing_t/PBTools/'))
-  sys.path.insert(0, os.path.expanduser('/sf/cristallina/applications/mx/zamofing_t/PBTools/'))
-
-
-import logging
-_log=logging.getLogger(__name__)
-
-from PyQt5.QtWidgets import (QApplication,)
-from app_config import AppCfg #settings, option, toggle_option
-
-
-from pbtools.misc.pp_comm import PPComm
-from pbtools.misc.gather import Gather
-import shapepath
-
-class Deltatau:
-  def __init__(self):
-    app=QApplication.instance()
-    cfg=app._cfg
-    host=cfg.value(AppCfg.DT_HOST)
-
-    hpp=host.split(':')
-    param={'host':hpp[0]}
-    if len(hpp)>1:
-      param['port']=int(hpp[1])
-    if len(hpp)>2:
-      param['fast_gather_port']=int(hpp[2])
-      _log.info(' -> ssh-tunneling PPComm({host}:{port} {host}:{fast_gather_port})'.format(**param))
-    try:
-      self._comm=comm=PPComm(**param,timeout=2.0)
-      self._gather=gather=Gather(comm)
-    except (socket.timeout,socket.gaierror) as e:
-      _log.critical(f'can not connect to deltatau:"{host}" -> {e}')
-      self._comm=comm=None
-      self._gather=gather=None
-      #return
-    self._shapepath=sp=shapepath.ShapePath(comm, gather, verbose=0xff, sync_mode=1, sync_flag=3)

geometry.py

@@ -247,13 +247,56 @@ class geometry:
 
     _log.debug('least square data:\nK:{}\nAA:{}'.format(K, AA))
 
-  def autofocus(self):
+  @staticmethod
+  def autofocus(cam,mot,rng=(-1,1),n=30):
     # cam   camera object
     # mot   motor object
     # rng   region (min max relative to current position) to seek
     # n     number of images to take in region
     # roi   region of interrest to calculate sharpness
     # mode  mode to calculate sharpness (sum/max-min/hist?  of edge detection in roi)
+    import PIL.Image
+    from scipy import ndimage
+    v=np.ndarray(shape=(len(cam),2))
+    if type(cam) == list:
+      for i, fn in enumerate(cam):
+        img=PIL.Image.open(fn)
+        img=np.asarray(img)
+        s=ndimage.sobel(img)
+        v[i,0]=s.sum()
+        v[i,1]=s.std()
+        #fft=np.log(np.abs(np.fft.fft2(img)))
+        #fft=np.fft.fftshift(fft)
+        #s=np.array(fft.shape,dtype=np.uint16)/2
+        #fft[300:700,400:800]=0
+        #v[i,1]=fft.sum()
+        #if i&0x3==0:
+        #  plt.figure()
+        #  plt.imshow(fft)
+      fig, ax=plt.subplots()
+      mx=v.max(0)
+      mn=v.min(0)
+      v=(v-mn)/(mx-mn)
+      #ax.plot(v[:,0])
+      ax.plot(v)
+      plt.show()
+      pass
+    else:
+      p0=mot.get_rbv()
+      for i,p in enumerate(np.linspace(p0+rng[0],p0+rng[1],n)):
+        mot.move_abs(p,wait=True)
+        pic=cam._pic# get_image()
+        mx=pic.max()
+        if pic.max()>255:
+          scl=2**int(round(np.log2(mx)-8))
+          pic=np.array(pic/scl, dtype=np.uint8)
+        elif pic.dtype!=np.uint8:
+          pic=np.array(pic, dtype=np.uint8)
+
+        img=PIL.Image.fromarray(pic)
+        fn=f'/tmp/image{i:03d}.png'
+        img.save(fn)
+        _log.debug(f'{fn} {pic.dtype} {pic.min()} {pic.max()}')
     pass
 
   def pix2pos(self, p, zoom=None):
@@ -558,6 +601,11 @@ if __name__=="__main__":
              [2.84446241, 3.54734879, 3.1415]])
       plane=geometry.least_square_plane(pts)
 
+  if args.mode&0x08:
+    import glob
+    imgLst=sorted(glob.glob("scratch/image*.png"))
+    geometry.autofocus(imgLst,None)
+
 
 # pix2pos="[[1.0, 200.0, 400.0, 600.0, 800.0], [[[0.001182928623952055, -2.6941995127711305e-05], [-4.043716694634124e-05, -0.0011894314084263825]], [[0.0007955995220142541, -3.175003727901119e-05], [-2.0896601103372113e-05, -0.0008100805094631365]], [[0.00048302539335378367, -1.1661121407652543e-05], [-2.0673746995751222e-05, -0.0004950857899461772]], [[0.00028775903460576395, -1.3762555219494508e-05], [-9.319936861519268e-06, -0.0002889214488565999]], [[0.0001788819256630411, -6.470841493681516e-06], [-2.0336605088889967e-06, -0.0001831131753499113]]]]"
 

psi_device.py

@@ -0,0 +1,89 @@
+import sys,os,socket
+sys.path.insert(0, os.path.expanduser('..'))
+hostname=socket.gethostname()
+if hostname=='ganymede':
+  sys.path.insert(0, os.path.expanduser('~/Documents/prj/SwissFEL/PBTools/'))
+else:
+  sys.path.insert(0, os.path.expanduser('/sf/cristallina/applications/mx/zamofing_t/PBTools/'))
+  sys.path.insert(0, os.path.expanduser("/photonics/home/gac-cristall/Documents/swissmx_cristallina/slic/"))
+  from slic.core.acquisition import SFAcquisition
+
+
+import logging
+_log=logging.getLogger(__name__)
+
+from PyQt5.QtWidgets import (QApplication,)
+from app_config import AppCfg #settings, option, toggle_option
+
+
+from pbtools.misc.pp_comm import PPComm
+from pbtools.misc.gather import Gather
+import shapepath
+
+try:
+  from slic.core.acquisition import SFAcquisition
+except ImportError as e:
+  _log.warning(e)
+
+class Deltatau:
+  def __init__(self):
+    app=QApplication.instance()
+    cfg=app._cfg
+    host=cfg.value(AppCfg.DT_HOST)
+
+    hpp=host.split(':')
+    param={'host':hpp[0]}
+    if len(hpp)>1:
+      param['port']=int(hpp[1])
+    if len(hpp)>2:
+      param['fast_gather_port']=int(hpp[2])
+      _log.info(' -> ssh-tunneling PPComm({host}:{port} {host}:{fast_gather_port})'.format(**param))
+    try:
+      self._comm=comm=PPComm(**param,timeout=2.0)
+      self._gather=gather=Gather(comm)
+    except (socket.timeout,socket.gaierror) as e:
+      _log.critical(f'can not connect to deltatau:"{host}" -> {e}')
+      self._comm=comm=None
+      self._gather=gather=None
+      #return
+    self._shapepath=sp=shapepath.ShapePath(comm, gather, verbose=0xff, sync_mode=1, sync_flag=3)
+
+class Jungfrau:
+  def __init__(self):
+    #python /sf/jungfrau/applications/daq_client/daq_client.py -h
+    #python /sf/jungfrau/applications/daq_client/daq_client.py -p p19739 -t no_beam_test
+    # -c/sf/cristallina/config/channel_lists/channel_list_bs-e/sf/cristallina/config/channel_lists/channel_list_ca
+    # -f/sf/cristallina/config/jungfrau/jf_1d5M.json
+    # --start_pulseid 15382163895 --stop_pulseid 15382163905
+    #rsync -vai gac-cristall@saresc-cons-03:/sf/jungfrau/applications/daq_client/daq_client.py .
+
+
+    # setup slic parameters
+    detectors = [ { "name" : "JF17T16V01",
+                            "adc_to_energy" : True,
+                            "compression" : True,
+                            "factor" : 11.33,
+                            "geometry" : True,
+                            "double_pixel_action" : "mask",
+                            "remove_raw_files" : False
+                        } ] #["JF17T16V01"]
+    bs_channels = [ # have timing signitures
+        "SARES30-LSCP1-CRISTA1:CH0:1",
+    #    "SARES30-CAMS156-SMX-OAV:FPICTURE"
+    ]
+    epics_channels = [ # no time signitures
+    ]
+
+    try:
+      self._daq=SFAcquisition(
+        "cristallina", "p20516",
+        default_detectors=detectors, default_channels=bs_channels, default_pvs=epics_channels,
+        rate_multiplicator=1, append_user_tag_to_data_dir=True
+      )
+    except NameError as e:
+      _log.warning(f'Jungfrau not connected: {e}')
+      self._daq=None
+
+  def acquire(self,run_name, n_pulses, wait=False):
+    if self._daq is not None:
+      self._daq.acquire( run_name, n_pulses=n_pulses, wait=False)

swissmx.py

@@ -100,8 +100,7 @@ import epics
 from epics.ca import pend_event
 import camera,backlight,zoom,illumination,geometry
 ts.log('Import part 7/8:')
-import deltatau
+import psi_device
-import detector
 ts.log('Import part 8/8:')
 
 def tdstamp():
@@ -362,7 +361,7 @@ class WndSwissMx(QMainWindow, Ui_MainWindow):
     grp.setTransform(tr)  # assign transform
 
     #--- image ---
-    self._goImg=img=pg.ImageItem()
+    self._goImg=img=pg.ImageItem() #border=pg.mkPen('r',width=2))
     grp.addItem(img)
 
     #--- opctical center ----
@@ -1080,7 +1079,7 @@ class WndSwissMx(QMainWindow, Ui_MainWindow):
 
 
         #create trace
-        self.add_camera_trace()
+        self.add_camera_trace(pFix)
 
         _log.debug(f'dist to beam ({pFix[0]:>0.6g} {pFix[1]:>0.6g}mm)')
         fx_motor=self.tweakers["fast_x"]
@@ -1301,7 +1300,7 @@ class WndSwissMx(QMainWindow, Ui_MainWindow):
     grp.setTransform(tr)
     grpTrc.setTransform(tr)
 
-  def add_camera_trace(self):
+  def add_camera_trace(self,delta_mm):
     # add up to tracelen image when moveing motors to new positions
     try:
       tracelen=self._goTrace._tracelen
@@ -1316,7 +1315,8 @@ class WndSwissMx(QMainWindow, Ui_MainWindow):
     trImg=self._goImgGrp.transform()
     trTrk=self._goTracked.transform()
     tr=trImg*trTrk.inverted()[0]  # TODO: check if oder correct !!!
-    img=pg.ImageItem(self._goImg.image)
+    cm=pg.ColorMap((0,255),((0,0,64),(255,255,128)))#,opacity=.5)
+    img=pg.ImageItem(self._goImg.image,lut=cm.getLookupTable(0,255))
     img.setTransform(tr)
     self._goTrace.addItem(img) # automatically added to self.vb
 
@@ -1598,13 +1598,18 @@ Author Thierry Zamofing (thierry.zamofing@psi.ch)
       plt.stem(x, y)
       plt.show(block=False)
 
-    step=1
+    step=2
     if step==0:
       vb=self.vb
       vb.autoRange(items=(self._goImg,))
     elif step==1:
       testMatplotlib()
-    if step==3:
+    elif step==2:
+      app=QApplication.instance()
+      cfg=app._cfg
+      geo=app._geometry
+      geo.autofocus(app._camera, self.tweakers['base_z'])
+    elif step==3:
       grp=pg.ItemGroup()
       vb.addItem(grp)
       obj=UsrGO.Marker((100, 100), (100, 100), mode=1)
@@ -2019,11 +2024,11 @@ Author Thierry Zamofing (thierry.zamofing@psi.ch)
     try:
       dt=app._deltatau
     except AttributeError:
-      app._deltatau=dt=deltatau.Deltatau()
+      app._deltatau=dt=psi_device.Deltatau()
     try:
       jf=app._jungfrau
     except AttributeError:
-      app._jungfrau=jf=detector.Jungfrau()
+      app._jungfrau=jf=psi_device.Jungfrau()
 
     sp=dt._shapepath
     sp.verbose=dt_misc['verbose']
@@ -2042,7 +2047,8 @@ Author Thierry Zamofing (thierry.zamofing@psi.ch)
     sp.homing()  # homing if needed
     sp.run()  # start motion program
     sp.wait_armed()  # wait until motors are at first position
-    sp.trigger(0.5)  # send a start trigger (if needed) ater given time
+    sp.trigger(0.5)  # send a start trigger (if needed) after given time
+    jf.acquire('filename.tmp',sp.points.shape[0])
     if not dt._comm is None:
       while True:
         p=sp.progress()

workbench/Readme.md

@@ -0,0 +1,8 @@
+testApps
+--------
+
+this directory contains test applications to analyse data and try basig functionalities.
+
+'scratch' contains quick and dirty basic tests
+
+'testApps' contains clean small applications to tests e.g. image analysis stuff.

workbench/autofocus.py

@@ -0,0 +1,94 @@
+# -*- coding: utf-8 -*-
+"""
+Demonstrate a simple data-slicing task: given 3D data (displayed at top), select
+a 2D plane and interpolate data along that plane to generate a slice image
+(displayed at bottom).
+
+
+"""
+
+## Add path to library (just for examples; you do not need this)
+
+import numpy as np
+from pyqtgraph.Qt import QtCore, QtGui
+import pyqtgraph as pg
+
+app=QtGui.QApplication([])
+
+## Create window with two ImageView widgets
+win=QtGui.QMainWindow()
+win.resize(800, 800)
+win.setWindowTitle('pyqtgraph example: DataSlicing')
+cw=QtGui.QWidget()
+win.setCentralWidget(cw)
+l=QtGui.QGridLayout()
+cw.setLayout(l)
+imv1=pg.ImageView()
+imv2=pg.ImageView()
+l.addWidget(imv1, 0, 0)
+l.addWidget(imv2, 1, 0)
+sld=QtGui.QSlider(QtCore.Qt.Horizontal)
+sld.setMinimum(10)
+sld.setMaximum(30)
+sld.setValue(20)
+sld.setTickPosition(QtGui.QSlider.TicksBelow)
+sld.setTickInterval(5)
+
+l.addWidget(sld, 2, 0)
+win.show()
+
+roi=pg.LineSegmentROI([[10, 64], [120, 64]], pen='r')
+imv1.addItem(roi)
+
+x1=np.linspace(-30, 10, 128)[:, np.newaxis, np.newaxis]
+x2=np.linspace(-20, 20, 128)[:, np.newaxis, np.newaxis]
+y=np.linspace(-30, 10, 128)[np.newaxis, :, np.newaxis]
+z=np.linspace(-20, 20, 128)[np.newaxis, np.newaxis, :]
+d1=np.sqrt(x1**2+y**2+z**2)
+d2=2*np.sqrt(x1[::-1]**2+y**2+z**2)
+d3=4*np.sqrt(x2**2+y[:, ::-1]**2+z**2)
+data=(np.sin(d1)/d1**2)+(np.sin(d2)/d2**2)+(np.sin(d3)/d3**2)
+
+import PIL.Image
+from scipy import ndimage
+import glob
+imgLst=sorted(glob.glob("image*.png"))
+v=np.ndarray(shape=(len(imgLst), 2))
+
+#for i, fn in enumerate(imgLst):
+#    img=PIL.Image.open(fn)
+#    img=np.asarray(img)
+#    s=ndimage.sobel(img)
+#    v[i, 0]=s.sum()
+#    v[i, 1]=s.std()
+  #fig, ax=plt.subplots()
+  #mx=v.max(0)
+  #mn=v.min(0)
+  #v=(v-mn)/(mx-mn)
+  # ax.plot(v[:,0])
+  #ax.plot(v)
+  #plt.show()
+#  pass
+
+
+def update():
+  global data, imv1, imv2, imgLst
+  d2=roi.getArrayRegion(data, imv1.imageItem, axes=(1, 2))
+  imv2.setImage(d2)
+
+
+roi.sigRegionChanged.connect(update)
+
+## Display the data
+imv1.setImage(data)
+imv1.setHistogramRange(-0.01, 0.01)
+imv1.setLevels(-0.003, 0.003)
+
+update()
+
+## Start Qt event loop unless running in interactive mode.
+if __name__=='__main__':
+  import sys
+
+  if (sys.flags.interactive!=1) or not hasattr(QtCore, 'PYQT_VERSION'):
+    QtGui.QApplication.instance().exec_()