auto-co.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "af2f3a19",
   "metadata": {},
   "source": [
    "This program will prepare a basic automation routine for functionalizing CO molecules on a clean Cu(111) surface. This program is provided to demonstrate of the automation capabilities for CO-AFM utilizing a CreaTec STM/AFM system. \n",
    "\n",
    "Programs from other publications may provide fully-autonomous construction and tip-preparation for STM tips, but here we wish to focus on the basics for CO functionalization:\n",
    " - loading the model,\n",
    " - utilizing CV tools to help segment your images,\n",
    " - functionalizing the tip,\n",
    " - using the model to determine the quality of CO functionalization.\n",
    "\n",
    "This is a tool designed for SPM practitioners who have some understanding of Python programming, machine learning, and hardware automation. Our hope is that hardware manufacturers will integrate some of these techniques into their own programs to aid experimental researchers speed up their own processes. The reason for this program is that automation STM tip preparation has been shown in numerous examples, but functionalization automation is relatively unexplored. \n",
    "\n",
    "Keep in mind that automated CO functionalization does not mean taking into account every corner case imaginable. Surface variations, impurities, and cleaning targets will need to be integrated into future models and programs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "bb45eaf9",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Initialization of libraries ###\n",
    "import numpy as np\n",
    "import win32com.client\n",
    "import tensorflow as tf\n",
    "from tensorflow import keras\n",
    "import cv2\n",
    "import time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "322c3b03",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Load CreaTec scan read parameters ###\n",
    "# Channel Selection #\n",
    "chantopo     = 1\n",
    "chancurrent  = 2\n",
    "\n",
    "# Unit Selection #\n",
    "unit_default = 0\n",
    "unit_volt    = 1\n",
    "unit_dac     = 2\n",
    "unit_ampere  = 3\n",
    "unit_nm      = 4\n",
    "unit_Hz      = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "358d7e5c",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Load CreaTec scan parameters ###\n",
    "# This assumes a 512x512 pixel size with 400x400 Å size and a scanning speed equivalent to about 5-6 minutes per image. Why?\n",
    "# This is approximately the limit of reasonable CO imaging and gives you a large area to test with. Too large of an image\n",
    "# may yield odd CO images. The program will attempt to scale smaller images as well, but your results may be poor.\n",
    "def SetScanParameters():\n",
    "    # Set X x Y image size (Pixels) #\n",
    "    NumX = 512\n",
    "    NumY = 512\n",
    "    stm.setparam('Num.X', NumX)\n",
    "    stm.setparam('Num.Y', NumY)\n",
    "\n",
    "    # Set the X x Y image size (Å)\n",
    "    DeltaX = 119\n",
    "    DeltaY = 119\n",
    "    stm.setparam('Delta X [DAC]', DeltaX)\n",
    "    stm.setparam('Delta Y [DAC]', DeltaY)\n",
    "\n",
    "    # Set the speed #\n",
    "    stm.setparam('DX/DDeltaX', 33)\n",
    "    \n",
    "    # Set the Topography and Current channels #\n",
    "    stm.setparam('Channelselectval', 3)\n",
    "\n",
    "    # Set to Forward + Backward\n",
    "    stm.setparam('ScanXMode', 0)\n",
    "\n",
    "    # Set to Constant Current mode #\n",
    "    stm.setparam('CHMode', 0)\n",
    "\n",
    "    # Set Bias Voltage (mV) #\n",
    "    stm.setparam('BiasVolt.[mV]', 100)\n",
    "\n",
    "    # Set Setpoint Current (pA) #\n",
    "    stm.setparam('FBLogIset', 50)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "id": "43e10a8e",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Load point spectra parameters ###\n",
    "def CO_pickup():\n",
    "    # Set the voltage #\n",
    "    stm.setparam('Vpoint0.V', 0)\n",
    "    stm.setparam('Vpoint0.V', 2600)\n",
    "    stm.setparam('Vpoint0.V', 2600)\n",
    "    stm.setparam('Vpoint0.V', 2600)\n",
    "    stm.setparam('Vpoint0.V', 0)\n",
    "\n",
    "    # Set the time #\n",
    "    stm.setparam('Vpoint0.t', 0)\n",
    "    stm.setparam('Vpoint0.t', 20)\n",
    "    stm.setparam('Vpoint0.t', 500)\n",
    "    stm.setparam('Vpoint0.t', 980)\n",
    "    stm.setparam('Vpoint0.t', 1000)\n",
    "\n",
    "    # Set the parameters #\n",
    "    stm.setparam('Vertchannelselectval',4097)\n",
    "    stm.setparam('Vertmangain', 6)\n",
    "    stm.setparam('Vertmandelay', 100) # Need to calculate this value (clock cycle)\n",
    "    stm.setparam('VertSpecBack', 0)\n",
    "    stm.setparam('VertSpecAvrgnr', 1)\n",
    "    stm.setparam('VertAvrgdelay', 100)\n",
    "    stm.setparam('VertRepeatCounter', 1)\n",
    "    stm.setparam('VertLineCount', 2)\n",
    "    stm.setparam('VertFBMode', 4)\n",
    "    stm.setparam('VertFBLogiset', 0.400)\n",
    "    stm.setparam('VertLatddx', 15.223)\n",
    "    stm.setparam('VertLatdelay', 166)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "id": "5a82ad02",
   "metadata": {},
   "outputs": [],
   "source": [
    "def CO_dropoff():\n",
    "    # Set the voltage #\n",
    "    stm.setparam('Vpoint0.V', 0)\n",
    "    stm.setparam('Vpoint0.V', 3300)\n",
    "    stm.setparam('Vpoint0.V', 3300)\n",
    "    stm.setparam('Vpoint0.V', 3300)\n",
    "    stm.setparam('Vpoint0.V', 0)\n",
    "\n",
    "    # Set the time #\n",
    "    stm.setparam('Vpoint0.t', 0)\n",
    "    stm.setparam('Vpoint0.t', 20)\n",
    "    stm.setparam('Vpoint0.t', 500)\n",
    "    stm.setparam('Vpoint0.t', 980)\n",
    "    stm.setparam('Vpoint0.t', 1000)\n",
    "\n",
    "    # Set the parameters #\n",
    "    stm.setparam('Vertchannelselectval', 4097)\n",
    "    stm.setparam('Vertmangain', 6)\n",
    "    stm.setparam('Vertmandelay', 100) # Need to calculate this value (clock cycle)\n",
    "    stm.setparam('VertSpecBack', 0)\n",
    "    stm.setparam('VertSpecAvrgnr', 1)\n",
    "    stm.setparam('VertAvrgdelay', 100)\n",
    "    stm.setparam('VertRepeatCounter', 1)\n",
    "    stm.setparam('VertLineCount', 2)\n",
    "    stm.setparam('VertFBMode', 4)\n",
    "    stm.setparam('VertFBLogiset', 0.400)\n",
    "    stm.setparam('VertLatddx', 15.223)\n",
    "    stm.setparam('VertLatdelay', 166)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "0bca8b7f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Start a Scan\n",
    "def RunScan():\n",
    "    stm.scanstart()\n",
    "    while stm.scanstatus > 0:\n",
    "        time.sleep(0.1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "id": "ff288dd8",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Read the scan image ###\n",
    "def ReadSTMImage(channel, unit):\n",
    "    image_STM = stm.scandata(channel,unit)\n",
    "    image = np.asarray(image_STM)\n",
    "    return image"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "id": "f5ffd9f0",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Normalize the image for import into OpenCV ###\n",
    "def NormalizeSTMImage(image):\n",
    "    # Normalize the image which is needed for later usage in OpenCV tools.\n",
    "    img_n = cv2.normalize(src=image, dst=None, alpha=0, beta=255, norm_type=cv2.NORM_MINMAX, dtype=cv2.CV_8U)\n",
    "    return img_n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 179,
   "id": "eb260b7c",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Size filter to remove points with areas significantly larger features ###\n",
    "def FilterCO(kp):\n",
    "    # Gather the CO sizes\n",
    "    s = [int(kp[i].size) for i in range(len(kp))]\n",
    "    \n",
    "    # Extract COs based on the feature size. This is needed for SURF, where the feature detection and CO points are somewhat\n",
    "    # correlated. A better way to do this is HoughCircle in OpenCV, but this has a habit of detecting COs at the edge.\n",
    "    kp_XY = [kp[i].pt for i in range(len(kp)) if kp[i].size <= np.median(s)]\n",
    "    \n",
    "    # Extract the CO sizes. This value isn't really used, but was in testing.\n",
    "    kp_size = [kp[i].size for i in range(len(kp)) if kp[i].size <= np.median(s)]\n",
    "    return s, kp_XY, kp_size"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 230,
   "id": "adbd0a51",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Generate mask for COs for import into empty area function\n",
    "def GenerateCOMask(image, kp_XY):\n",
    "    # Extract dimensions from image\n",
    "    height, width = image.shape[:2]\n",
    "\n",
    "    # Generate a blank mask of all zeros\n",
    "    binaryMap = np.zeros((width, height))\n",
    "\n",
    "    # Extract the points from the filtered CO list and convert to X,Y lists\n",
    "    CO_list = [(int(element[0]), int(element[1])) for element in kp_XY]\n",
    "    \n",
    "    # Generate a circle mask around each point of similar size to the original radius\n",
    "    for i in range(len(CO_list)):\n",
    "        cv2.circle(binaryMap,(CO_list[i][0],CO_list[i][1]),int(np.median(s)/2),(255,255,255),thickness=-1)\n",
    "\n",
    "    # Output a binary map\n",
    "    th, dst = cv2.threshold(binaryMap, 0, 255, cv2.THRESH_BINARY)\n",
    "    \n",
    "    return th, dst"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 229,
   "id": "53872bdc",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Crop and Store COs ###\n",
    "def CropCO(norm_image,kp_XY):\n",
    "    # Set an empty array to fill with images\n",
    "    CO_crop = np.array([])\n",
    "\n",
    "    # Crop each CO with input size to feed into model\n",
    "    for i in range(len(circles)):\n",
    "        \n",
    "        # Set a zoom to gather the entire CO (center +surrounding halo)\n",
    "        zoom = 1.5\n",
    "        # Convert the kp coordinates for easier reference\n",
    "        kp_coord = kp_XY[i].pt\n",
    "        \n",
    "        # Set the radius to include zoom value. Expands the surrounding rectangle (square)\n",
    "        radius = (kp_XY[i].size/zoom).astype(\"int\")\n",
    "        \n",
    "        # Crop at the target coordinates\n",
    "        CO_crop.append(norm_image[kp_coord[0]-radius:kp_coord[1]+radius, kp_coord[0]-radius:kp_coord[1]+radius])\n",
    "    return CO_crop"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 226,
   "id": "b7be245c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def FindEmptySquareArea(mat, ZERO=0):\n",
    "    #Find the largest square of ZERO's in the matrix `mat`.\n",
    "    #Source: https://stackoverflow.com/a/1726667\n",
    "    \n",
    "    # Extract array shape\n",
    "    nrows, ncols = mat.shape\n",
    "    \n",
    "    # Check for null array / empty matrix or rows\n",
    "    if not (nrows and ncols):\n",
    "        return 0\n",
    "    \n",
    "    # Initiate null array based on row and column sizes\n",
    "    counts = np.zeros((nrows, ncols))\n",
    "    \n",
    "    # For each row\n",
    "    for i in reversed(range(nrows)):\n",
    "        # The matrix must be rectangular (actually square)\n",
    "        assert len(mat[i]) == ncols\n",
    "        \n",
    "        # For each element in the row\n",
    "        for j in reversed(range(ncols)):\n",
    "            \n",
    "            # If the element is zero, check around it\n",
    "            if mat[i][j] == ZERO:\n",
    "                counts[i][j] = (1 + min(\n",
    "                    counts[i][j+1],   # east\n",
    "                    counts[i+1][j],   # south\n",
    "                    counts[i+1][j+1]  # south-east\n",
    "                    )) if i < (nrows - 1) and j < (ncols - 1) else 1  # edges\n",
    "    mx = -1\n",
    "    lx = -1\n",
    "    ly = -1\n",
    "    for row in range(len(counts)):\n",
    "        for col in range(len(counts[row])):\n",
    "            if counts[row, col] > mx:\n",
    "                mx = counts[row, col]\n",
    "                ly = row\n",
    "                lx = col\n",
    "\n",
    "    return mx, ly, lx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 228,
   "id": "98fa096d",
   "metadata": {},
   "outputs": [],
   "source": [
    "def FindEmptySquareAreaPosition(img):\n",
    "    # Returns the Y,X center of the largest empty square in an image and size.\n",
    "    # Prints a nice reference image for sanity check.\n",
    "    mx, ly, lx = FindEmptySquareArea(img, 0)\n",
    "    \n",
    "    # Calculate the center\n",
    "    ly = ly+mx/2\n",
    "    lx = lx+mx/2\n",
    "    \n",
    "    # Visual output of the square\n",
    "    # Normalize the image\n",
    "    img = (img-img.min())/(img.max()-img.min())\n",
    "    img *= 255\n",
    "    ly = int(ly)\n",
    "    lx = int(lx)\n",
    "\n",
    "    # Draw a box around the edge of the square\n",
    "    img[int(ly),int(lx)]=255\n",
    "    img[int(ly-1),int(lx-1)]=255\n",
    "    img[int(ly-1),int(lx+1)]=255\n",
    "    img[int(ly+1),int(lx-1)]=255\n",
    "    img[int(ly+1),int(lx+1)]=255\n",
    "    for i in range(int(mx)):\n",
    "        img[int(ly+i-mx/2),int(lx-mx/2)] = 255\n",
    "        img[int(ly+i-mx/2),int(lx+mx-mx/2)] = 255\n",
    "        img[int(ly-mx/2),int(lx+i-mx/2)] = 255\n",
    "        img[int(ly+mx-mx/2),int(lx+i-mx/2)] = 255\n",
    "        \n",
    "    # Print the image (quick reference)\n",
    "    #plt.imshow(img, cmap = \"gray\");\n",
    "    return ly, lx, mx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "c4bf3505",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Initialize CreaTec STMAFM program ###\n",
    "stm = win32com.client.Dispatch(\"pstmafm.stmafmrem\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "8f993acd",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Import Keras model ###\n",
    "reconstructed_model = keras.models.load_model(\"/pretrained_weights/model.h5\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "840bf609",
   "metadata": {},
   "outputs": [],
   "source": [
    "### Simplified program to run autonomous CO functionalization on a single scan area ###\n",
    "\n",
    "score = 0\n",
    "target = 0.9\n",
    "CO_count = 0\n",
    "pulse = 0\n",
    "\n",
    "while score < target:\n",
    "    # Load STM Scan Parameters #\n",
    "    SetScanParameters()\n",
    "    \n",
    "    # Start Scan #\n",
    "    RunScan()\n",
    "    \n",
    "    # Load STM Image #\n",
    "    stm_image = ReadSTMImage(chantopo, unit_volt)\n",
    "    \n",
    "    # Convert to Image Format Compatible with SURF\n",
    "    norm_image = NormalizeSTMImage(stm_image)\n",
    "    \n",
    "    # Load SURF parameters \n",
    "    surf = cv2.xfeatures2d.SURF_create(3000)\n",
    "    \n",
    "    # Find keypoints and descriptors\n",
    "    kp, des = surf.detectAndCompute(norm_image,None)\n",
    "    \n",
    "    # If no COs, alert operator\n",
    "    if len(kp) == 0:\n",
    "        print(\"No objects detected. Please check scan area or image dimensions.\")\n",
    "        break\n",
    "    \n",
    "    # Filter out non-CO positions\n",
    "    s, kp_XY, kp_size = FilterCO(kp)\n",
    "    \n",
    "    # Plot the SURF identified features\n",
    "    surf_image = cv2.drawKeypoints(norm_image,kp_XY,None,(255,0,0),4)\n",
    "    #plt.imshow(surf_image);\n",
    "    \n",
    "    # Extract the points from the filtered CO list and convert to X,Y lists\n",
    "    CO_list = [(int(element[0]), int(element[1])) for element in kp_XY]\n",
    "    \n",
    "    ### Sort the centerpoints from top to bottom, left to right ###\n",
    "    CO_list_sorted = CO_list[np.lexsort((CO_list[:,0], CO_list[:,1]))]\n",
    "    \n",
    "    # Generate a CO mask to identify empty areas\n",
    "    _, CO_mask = GenerateCOMask(norm_image, CO_list_sorted)\n",
    "    \n",
    "    # Crop and Store CO images\n",
    "    # Cropped images need to be scaled to 32x32 images for analysis by model\n",
    "    CO_crop = CropCO(norm_image, circles)\n",
    "    CO_crop = np.array(CO_crop, dtype=\"object\")\n",
    "    resized_CO = np.array([cv2.resize(np.copy(CO_crop), (32, 32)) for CO in CO_crop])\n",
    "    \n",
    "    # Analyze with model\n",
    "    add_norm_CO(resized_CO)\n",
    "    outputs = model.predict_on_batch(resized_CO)\n",
    "    probs = np.squeeze(np.array(outputs))\n",
    "    score = np.median(probs)\n",
    "    \n",
    "    if CO_count >= 10:\n",
    "        # Load dropoff parameters\n",
    "        CO_dropoff()\n",
    "        \n",
    "        # Target tip to an empty area\n",
    "        ly, lx, mx = FindEmptySquareAreaPosition(CO_mask)\n",
    "        \n",
    "        # Run vertical manipulation\n",
    "        stm.btn_vertspec(lx, ly)\n",
    "        \n",
    "        # Brief sleep due to STMAFM program quirks\n",
    "        time.sleep(0.05)\n",
    "        \n",
    "        # Track the number of tip pulses\n",
    "        pulse += pulse\n",
    "        \n",
    "        # Reset counter back to 0\n",
    "        CO_count = 0\n",
    "    else:\n",
    "        # Score comparison\n",
    "        # Typically if the score is less than 0.005, this is, from initial tests, no CO.\n",
    "        # Not always true, but further refinements of the program needed to detect blank or extremely tilted COs.\n",
    "        if score < 0.005:\n",
    "            \n",
    "            # Load pickup parameters\n",
    "            CO_pickup()\n",
    "            \n",
    "            # Run vertical manipulation at first point in list\n",
    "            stm.btn_vertspec(CO_list_sorted[0][0], CO_list_sorted[0][1])\n",
    "            \n",
    "            # Brief sleep due to STMAFM program quirks\n",
    "            time.sleep(0.05)\n",
    "            \n",
    "            # Track the number of attempted CO pickups\n",
    "            CO_count += CO_count\n",
    "        \n",
    "        # For scores between 0.005 and the target, CO or O-terminated is already attached. Drop-off is needed.\n",
    "        elif score >= 0.005 and score < target:\n",
    "            \n",
    "            # Load dropoff parameters\n",
    "            CO_dropoff()\n",
    "\n",
    "            # Target tip to an empty area\n",
    "            ly, lx, mx = FindEmptySquareAreaPosition(CO_mask)\n",
    "\n",
    "            # Run vertical manipulation\n",
    "            stm.btn_vertspec(lx, ly)\n",
    "\n",
    "            # Brief sleep due to STMAFM program quirks\n",
    "            time.sleep(0.05)\n",
    "\n",
    "            # Track the number of tip pulses\n",
    "            pulse += pulse\n",
    "\n",
    "        # Final condition is greater than or equal to target, which means that the model identifies the image as CO-terminated\n",
    "        # within the target parameters. \n",
    "        else:\n",
    "            print(\"CO termination completed\")\n",
    "            print(\"Target: \", target)\n",
    "            print(\"Score:\", score)\n",
    "            print(\"# of Pulses: \", pulse)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "createc2",
   "language": "python",
   "name": "createc2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}