8_Full_series_of_files-orig.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# RUN Main import block and TODO list\n",
    "\n",
    "# TODO: see how uri calculated the ridges\n",
    "\n",
    "# TODO: Perform Histogram equalization - start with it\n",
    "# TODO: \n",
    "# take integral from the Highest peak+-0.005 divide by integral of the entire graph \n",
    "# This will be the peakness measure for the PSD ==> The desired ridge index\n",
    "# TODO:\n",
    "# take integral from the Highest peak+-0.005 divide by integral of the entire graph - it's the peakness measure for the PSD\n",
    "# must select a peak above a min threshold in order to ignore noisy frequency\n",
    "# must ignore peaks above a certain threshold in order to detect meaningful frequency\n",
    "# run the PSD in moving windows every 200 px (deduced from the below PSD pointing to a freq of 1/0.02=50-> times 4= 200px)\n",
    "# and medianf the result of the windows\n",
    "# TODO:\n",
    "# Another alternative: (with Yariv)\n",
    "# Run PSD column by column - get the phase, freq, peakness and reconstruct an artificial ridge slice\n",
    "# from this - reconstruct a \"clean\" artificial ridge image\n",
    "\n",
    "%matplotlib inline\n",
    "import matplotlib.image as img\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import urllib.request\n",
    "import os\n",
    "import glob\n",
    "from scipy import ndimage\n",
    "from scipy import signal\n",
    "import cv2\n",
    "\n",
    "from PIL import Image, ImageDraw, ImageFont\n",
    "\n",
    "import mahotas as mh\n",
    "from mahotas import polygon\n",
    "# import pymorph as pm\n",
    "\n",
    "import networkx as nx\n",
    "\n",
    "from scipy import ndimage as nd\n",
    "import skimage.transform as transform\n",
    "import skimage.morphology as mp\n",
    "import skimage.io as sio\n",
    "import scipy.misc as sm\n",
    "from skimage.filters import threshold_otsu, threshold_adaptive\n",
    "from skimage.feature import hessian_matrix, hessian_matrix_eigvals\n",
    "from skimage import exposure\n",
    "from skimage import data, img_as_float\n",
    "from sklearn.metrics import precision_recall_fscore_support\n",
    "\n",
    "from bisect import bisect_left\n",
    "import math\n",
    "import warnings\n",
    "import csv\n",
    "import tensorflow as tf\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "code_folding": [
     0,
     113,
     131,
     148,
     160,
     195,
     201,
     206,
     211,
     218,
     223,
     228,
     236,
     243,
     252
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# RUN Utility functions\n",
    "\n",
    "# One time init\n",
    "# with open('results.csv', 'w') as csvfile:\n",
    "#     csvout = csv.writer(csvfile)\n",
    "#     csvout.writerow([\"File\", \"Model\", \"Gap\", \"Slice_size\", \"Count\", \"Precision\", \"Recall\", \"F-score\", \"True Count\", \"Error Rate\"])\n",
    "\n",
    "#BASIC CROP FRAME\n",
    "X_START = 1000\n",
    "X_END = 6000\n",
    "Y_START = 800\n",
    "Y_END = 4300\n",
    "BG_2_OBJ_RATIO = 0.92\n",
    "CUBE_SIZE = 250\n",
    "\n",
    "# Simple crop by x/y ranges\n",
    "def crop(image, ymin, ymax, xmin, xmax):\n",
    "    return image[ymin:ymax, xmin:xmax]\n",
    "\n",
    "# returns a logical matrix of values beyond a threshld\n",
    "def thresholded(image, val): \n",
    "    return np.logical_and(*[image[...] > val  for t in enumerate([0, 0])])\n",
    "\n",
    "def find_min_max_without_orphand_pixels(nonzero_dimension, crop_filter=20):\n",
    "    sorted = np.sort(nonzero_dimension)\n",
    "    prev=-1\n",
    "    min_val = sorted[0]\n",
    "    for i, x in enumerate(sorted[:100]):\n",
    "        if prev >= 0 and x - prev > crop_filter:\n",
    "            min_val = x\n",
    "        prev = x\n",
    "    prev=-1\n",
    "    max_val = sorted[-1]\n",
    "    for i, x in enumerate(sorted[-100:]):\n",
    "        if prev >= 0 and x - prev > crop_filter:\n",
    "            max_val = prev\n",
    "            break\n",
    "        prev = x\n",
    "    \n",
    "    return min_val, max_val\n",
    "\n",
    "def calc_min_max_coordinates(image, crop_val=50):\n",
    "    temp = thresholded(image, crop_val)\n",
    "    temp = temp * 1\n",
    "    temp = np.nonzero(temp)\n",
    "    ymin, ymax = find_min_max_without_orphand_pixels(temp[0])\n",
    "    xmin,xmax = find_min_max_without_orphand_pixels(temp[1])\n",
    "    return ymin, ymax, xmin, xmax\n",
    "\n",
    "def calc_min_max_coordinates_dynamic(image, cutoff=1):\n",
    "    temp = exposure.equalize_adapthist(image, clip_limit=0.03)\n",
    "    flat = np.sort(np.matrix.getA1(temp))\n",
    "    sum_all = np.sum(flat)\n",
    "    index = np.argmin(flat.cumsum() < (sum_all * cutoff))\n",
    "\n",
    "    temp = thresholded(temp, flat[index])\n",
    "    temp = temp * 1\n",
    "    temp = np.nonzero(temp)\n",
    "    ymin, ymax = find_min_max_without_orphand_pixels(temp[0])\n",
    "    xmin,xmax = find_min_max_without_orphand_pixels(temp[1])\n",
    "    return ymin, ymax, xmin, xmax\n",
    "\n",
    "# initial static crop and a seondary dynamic crop based on signal2noise ratio\n",
    "def crop_full_scan(image):\n",
    "    temp = crop(image, Y_START, Y_END, X_START, X_END)\n",
    "    ymin, ymax, xmin, xmax = calc_min_max_coordinates_dynamic(temp, cutoff=BG_2_OBJ_RATIO)\n",
    "    temp = crop(image, Y_START+ymin, Y_START+ymax, X_START+xmin, X_START+xmax)\n",
    "    return temp\n",
    "\n",
    "def crop_thresholded(image):\n",
    "    temp = crop(image, 0, image.shape[0]-1, 0, image.shape[1]-1)\n",
    "    ymin, ymax, xmin, xmax = calc_min_max_coordinates(temp)\n",
    "    temp = crop(image, ymin, ymax, xmin, xmax)\n",
    "    return temp\n",
    "\n",
    "def read_and_crop(image_name):\n",
    "    if \"il239838\" in os.getcwd():\n",
    "        image = img.imread(\"/Users/il239838/Downloads/private/Thesis/Papyrus/\" + image_name)\n",
    "    else:\n",
    "        f = urllib.request.urlopen(\"https://dl.dropboxusercontent.com/s/31b96942qdcn73k/\" + image_name)\n",
    "        image = img.imread(f, format='jpeg')\n",
    "\n",
    "    # Smart-crop the image to get rid of all the noise and redundant area\n",
    "    # return crop_full_scan(image)\n",
    "    cropped = crop_full_scan(image)\n",
    "    return exposure.equalize_adapthist(cropped, clip_limit=0.03)\n",
    "\n",
    "\n",
    "# TODO: fix performance!!! http://scikit-image.org/docs/dev/user_guide/tutorial_parallelization.html\n",
    "def combine_3_images_to_RGB(red, green, blue):\n",
    "    new_image = np.empty((blue.shape[0],blue.shape[1],3))\n",
    "    for x in range(0, blue.shape[0]):\n",
    "        for y in range(0, blue.shape[1]):\n",
    "            new_image[x,y,0] = red[x,y]\n",
    "            new_image[x,y,1] = green[x,y]\n",
    "            new_image[x,y,2] = blue[x,y]\n",
    "    return new_image\n",
    "\n",
    "def slice_image_left_edge(original, width=200, rotate=0):\n",
    "    rot = ndimage.rotate(original, rotate)\n",
    "    # Slice the left slice of the so-called \"blue\" image\n",
    "    left_edge_orig = crop(rot, 1, 1400, 1, width)\n",
    "    left_edge_orig = crop_thresholded(left_edge_orig)\n",
    "\n",
    "    # Copy to a new array so we don't thrash the origin\n",
    "    left_edge = np.empty_like (left_edge_orig)\n",
    "    np.copyto(left_edge, left_edge_orig)\n",
    "\n",
    "    # Zero down low level \"noise\" values\n",
    "    low_values_indices = left_edge < 30  # Where values are low\n",
    "    left_edge[low_values_indices] = 0  # All low values set to 0\n",
    "    return left_edge\n",
    "\n",
    "def get_best_angle_rotation(original, crop=True, width=200):\n",
    "    min_var = 99999999999\n",
    "    best_angle = -10\n",
    "    for x in range(-5,5):\n",
    "        if crop:            \n",
    "            rot_edge = slice_image_left_edge(original, width, x)\n",
    "        else:\n",
    "            rot_edge = ndimage.rotate(original, x)\n",
    "        left_var = np.var(rot_edge, axis=1)\n",
    "        # left_var = np.apply_along_axis(lambda v: np.var(v[np.nonzero(v)]), 1, rot_edge)\n",
    "        var_sum = np.sum(left_var)\n",
    "        if (var_sum < min_var):\n",
    "            min_var = var_sum\n",
    "            best_angle = x\n",
    "    print (\"best_angle=\"+str(best_angle))\n",
    "    return best_angle\n",
    "\n",
    "#     import pdb; pdb.set_trace()\n",
    "def calc_neighbors(slice_map, col, row):\n",
    "    # import pdb; pdb.set_trace()\n",
    "    if ((col-1, row) in slice_map and slice_map[(col-1, row)] != None):\n",
    "        slice_map[(col, row)][\"left\"] = slice_map[(col-1, row)]\n",
    "        slice_map[(col-1, row)][\"right\"] = slice_map[(col, row)]\n",
    "    if ((col+1, row) in slice_map and slice_map[(col+1, row)] != None):\n",
    "        slice_map[(col, row)][\"right\"] = slice_map[(col+1, row)]\n",
    "        slice_map[(col+1, row)][\"left\"] = slice_map[(col, row)]\n",
    "    if ((col, row-1) in slice_map and slice_map[(col, row-1)] != None):\n",
    "        slice_map[(col, row)][\"top\"] = slice_map[(col, row-1)]\n",
    "        slice_map[(col, row-1)][\"bottom\"] = slice_map[(col, row)]\n",
    "    if ((col, row+1) in slice_map and slice_map[(col, row+1)] != None):\n",
    "        slice_map[(col, row)][\"bottom\"] = slice_map[(col, row+1)]\n",
    "        slice_map[(col, row+1)][\"top\"] = slice_map[(col, row)]\n",
    "    \n",
    "\n",
    "\n",
    "def create_cube(raw, x, y):\n",
    "    cube = {}\n",
    "    cube[\"cube\"] = raw\n",
    "    cube[\"top_row\"] = x\n",
    "    cube[\"left_col\"] = y\n",
    "    cube[\"right_col\"] = y + CUBE_SIZE\n",
    "    return cube\n",
    "    \n",
    "\n",
    "ZERO_CUBE = create_cube(np.zeros((CUBE_SIZE, CUBE_SIZE), dtype=np.int), -1, -2)\n",
    "\n",
    "# slice an image to cubes with 250X250 pixel size\n",
    "def slice_to_static_slices(cropped_original):\n",
    "    structure = {}\n",
    "    # cropped_original = cropped_original / 256 # divide by 256 to \"normalize\" between 0 and 1\n",
    "    x, y = cropped_original.shape\n",
    "    print (x,y)\n",
    "    n = 0\n",
    "    # every 250 pixels on the x axis == rows\n",
    "    while ((n + 1) * CUBE_SIZE < x):\n",
    "        m = 0\n",
    "        # every 250 pixels on the y axis == cols\n",
    "        while ((m + 1) * CUBE_SIZE < y):\n",
    "            # cut a cube of 250X250\n",
    "            cube = (crop(cropped_original, n * CUBE_SIZE, (n + 1) * CUBE_SIZE, m * CUBE_SIZE, (m + 1) * CUBE_SIZE))\n",
    "            # keep only cubes for which half of the pixels have some \"color\"\n",
    "            # import pdb; pdb.set_trace()\n",
    "            # print(np.median(cube))\n",
    "            if np.median(cube) > 0.2: # aligned with the normalization 0.2 correlates to 50\n",
    "                # keep the cube\n",
    "                structure[(m, n)] = create_cube(cube, n * CUBE_SIZE, m * CUBE_SIZE)\n",
    "                structure[(m, n)][\"col\"] = m\n",
    "                structure[(m, n)][\"row\"] = n\n",
    "            else:\n",
    "                structure[(m, n)] = None\n",
    "            m += 1\n",
    "        n += 1\n",
    "    \n",
    "    # this loop has to be performed only after we've established all the None cubes\n",
    "    for cube in structure.values():\n",
    "        # set the reference to neighbor cubes\n",
    "        if cube != None:\n",
    "            calc_neighbors(structure, cube[\"col\"], cube[\"row\"])\n",
    "        \n",
    "    # return the data structure with all the cubes and the counters of the rows and columns\n",
    "    return structure, m, n\n",
    "\n",
    "def pad_above(original, above, amount):\n",
    "    res = np.insert(original[\"cube\"], np.zeros(amount), above[\"cube\"][-amount:], axis=0)\n",
    "    res = np.delete(res, np.arange(CUBE_SIZE,CUBE_SIZE+amount), axis=0)\n",
    "    return create_cube(res, original[\"top_row\"] - amount, original[\"left_col\"])\n",
    "  \n",
    "\n",
    "def pad_below(original, below, amount):\n",
    "    res = np.insert(original[\"cube\"], np.full(amount, CUBE_SIZE), below[\"cube\"][:amount], axis=0)\n",
    "    res = np.delete(res, np.arange(0, amount), axis=0)\n",
    "    return create_cube(res, original[\"top_row\"] + amount, original[\"left_col\"])\n",
    "\n",
    "def pad_left(original, left, amount):\n",
    "    res = np.insert(original[\"cube\"], np.zeros(amount, dtype=int), left[\"cube\"][:,-amount:], axis=1)\n",
    "    res = np.delete(res, np.arange(CUBE_SIZE, CUBE_SIZE+amount), axis=1)\n",
    "    return create_cube(res, original[\"top_row\"], original[\"left_col\"] - amount)\n",
    "\n",
    "def pad_right(original, right, amount):\n",
    "    res = np.insert(original[\"cube\"], [CUBE_SIZE], right[\"cube\"][:,:amount], axis=1)\n",
    "    res = np.delete(res, np.arange(0, amount), axis=1)\n",
    "    return create_cube(res, original[\"top_row\"], original[\"left_col\"] + amount)\n",
    "    \n",
    "\n",
    "# \"Shave\" the right edge of the cube with <gap> pixels and pad with zeros on the left\n",
    "def shave_right(original, amount):\n",
    "    return pad_left(original, ZERO_CUBE, amount)\n",
    "    \n",
    "\n",
    "# \"Shave\" the left edge of the cube with <gap> pixels and pad with zeros on the right    \n",
    "def shave_left(original, amount):\n",
    "    return pad_right(original, ZERO_CUBE, amount)\n",
    "    \n",
    "\n",
    "# concatenate cubes \n",
    "def concatenate_cubes(left, right, slice_size):\n",
    "    con = np.concatenate((left[\"cube\"][:,-slice_size:], right[\"cube\"][:,:slice_size]), axis=1)\n",
    "    x_delta = right[\"top_row\"] - left[\"top_row\"]\n",
    "    y_delta = right[\"left_col\"] - left[\"right_col\"] \n",
    "    return con, x_delta, y_delta\n",
    "    \n",
    "\n",
    "# concatenate cubes while artificially creating a gap between them. Pad the other end of the cube with zeros\n",
    "def concatenate_cubes_zero_pad_gaps(left_orig, right_orig, gap):\n",
    "    left = left_orig if gap == 0 else shave_right(left_orig, gap)\n",
    "    right = right_orig if gap == 0 else shave_left(right_orig, gap)\n",
    "    return concatenate_cubes(left, right)    \n",
    "\n",
    "# concatenate cubes while artificially creating a gap between them. Pad the other end of the cobe with the nearby\n",
    "# continuation of the cubes\n",
    "def concatenate_cubes_with_gap(left_orig, right_orig, gap, left_pad, right_pad, slice_size):\n",
    "    # import pdb; pdb.set_trace()\n",
    "    left = left_orig if gap == 0 else pad_left(left_orig, left_pad, gap)\n",
    "    right = right_orig if gap == 0 else pad_right(right_orig, right_pad, gap)\n",
    "    return concatenate_cubes(left, right, slice_size)        \n",
    "\n",
    "# convert the data structure of cubes into a train set of 2 arrays of images and labels\n",
    "# each image is a concatanation of 2 images from the original cubes set, covering all combinations of images\n",
    "# effectively creating Nx(N-1) images\n",
    "def build_train_set(rows, cols, cubes, gap, slice_size):\n",
    "    # import pdb; pdb.set_trace()\n",
    "    warnings.filterwarnings(\"ignore\")\n",
    "\n",
    "    train_imgs = []\n",
    "    train_lbls = []\n",
    "    train_x_delta = []\n",
    "    train_y_delta = []\n",
    "    # iterate over the rows and cols, essentially going over the grid of sliced cubes\n",
    "    for row in range(0, rows):\n",
    "        for col in range(0, cols):\n",
    "            # if this cube exists (could have been removed previously due to lack of data)\n",
    "            if (cubes[(col, row)] != None):\n",
    "                # for each \"current\" image in the iteration\n",
    "                curr = cubes[(col, row)]\n",
    "                # iterate over all the cubes to find all the \"other\" (adjacent) cubes\n",
    "                for adj_row in range(0, rows):\n",
    "                    for adj_col in range(0, cols):\n",
    "                        if (adj_row != row or adj_col != col):\n",
    "                            if (cubes[(adj_col, adj_row)] != None):\n",
    "                                adj = cubes[(adj_col, adj_row)]\n",
    "                                # append the adjacent image to the current image\n",
    "                                # pass the filling cubes on the right and left to pad against the gap\n",
    "                                if (gap == 0 or (\"left\" in curr.keys() and \"right\" in adj.keys())):\n",
    "                                    if (gap == 0):\n",
    "                                        conc, x_delta, y_delta = concatenate_cubes(curr, adj, slice_size)\n",
    "                                    else:\n",
    "                                        conc, x_delta, y_delta = concatenate_cubes_with_gap(curr, adj, gap, curr[\"left\"], adj[\"right\"], slice_size)\n",
    "                                    train_imgs.append(conc)\n",
    "                                    train_x_delta.append(x_delta)\n",
    "                                    train_y_delta.append(y_delta)\n",
    "                                    # if the adj image is on the same row and on the right of the curr image - it will be marked as match    \n",
    "                                    if (adj_row == row and adj_col == (col + 1)):\n",
    "                                        # mark the image as matched\n",
    "                                        train_lbls.append([0,1])\n",
    "                                        # need to enrich the set with a few more tru positive samples - so we offset \n",
    "                                        # the matched images up ad down a few times and create more matches\n",
    "                                        if (\"top\" in curr.keys() and \"top\"in adj.keys()):\n",
    "                                            for i in range(5, 101, 5):\n",
    "                                                curr1 = pad_above(curr, curr[\"top\"],i)\n",
    "                                                adj1 = pad_above(adj, adj[\"top\"],i)\n",
    "                                                if (gap == 0 or (\"left\" in curr.keys() and \"right\" in adj.keys() and \"top\" in curr[\"left\"].keys() and \"top\"in curr[\"right\"].keys())):\n",
    "                                                    if (gap == 0):\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes(curr1, adj1, slice_size)\n",
    "                                                    else:\n",
    "                                                        curr1Left = pad_above(curr[\"left\"], curr[\"left\"][\"top\"], i) # FIXIT?\n",
    "                                                        adj1Right = pad_above(adj[\"right\"], curr[\"right\"][\"top\"], i) # FIXIT?\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes_with_gap(curr1, adj1, gap, curr1Left, adj1Right, slice_size)\n",
    "                                                    train_imgs.append(conc)\n",
    "                                                    train_x_delta.append(x_delta)\n",
    "                                                    train_y_delta.append(y_delta)\n",
    "                                                    # mark the image as matched\n",
    "                                                    train_lbls.append([0,1])\n",
    "                                        if (\"bottom\" in curr.keys() and \"bottom\"in adj.keys()):\n",
    "                                            for i in range(5, 101, 5):\n",
    "                                                curr1 = pad_below(curr, curr[\"bottom\"],i)\n",
    "                                                adj1 = pad_below(adj, adj[\"bottom\"],i)\n",
    "                                                if (gap == 0 or (\"left\" in curr.keys() and \"right\" in adj.keys() and \"bottom\" in curr[\"left\"].keys() and \"bottom\"in curr[\"right\"].keys())):\n",
    "                                                    if (gap == 0):\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes(curr1, adj1, slice_size)\n",
    "                                                    else:\n",
    "                                                        curr1Left = pad_below(curr[\"left\"], curr[\"left\"][\"bottom\"], i) # FIXIT?\n",
    "                                                        adj1Right = pad_below(adj[\"right\"], curr[\"right\"][\"bottom\"], i) # FIXIT?\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes_with_gap(curr1, adj1, gap, curr1Left, adj1Right, slice_size)\n",
    "                                                    train_imgs.append(conc)\n",
    "                                                    train_x_delta.append(x_delta)\n",
    "                                                    train_y_delta.append(y_delta)\n",
    "                                                    # mark the image as matched\n",
    "                                                    train_lbls.append([0,1])\n",
    "                                        if (\"left\" in curr.keys()): # enough to check only the curr as the left of the adj is the curr\n",
    "                                            for i in range(5, 101, 5):\n",
    "                                                curr1 = pad_left(curr, curr[\"left\"],i)\n",
    "                                                adj1 = pad_left(adj, adj[\"left\"],i) # essentially the curr\n",
    "                                                if (gap == 0 or (\"left\" in curr.keys() and \"right\" in adj.keys())):\n",
    "                                                    if (gap == 0):\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes(curr1, adj1, slice_size)\n",
    "                                                    else:\n",
    "                                                        curr1Left = pad_left(curr[\"left\"], ZERO_CUBE, i) # FIXIT? + assuming the gap will not be more than 150\n",
    "                                                        adj1Right = pad_left(adj[\"right\"], ZERO_CUBE, i) # FIXIT? + assuming the gap will not be more than 150\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes_with_gap(curr1, adj1, gap, curr1Left, adj1Right, slice_size)\n",
    "                                                    train_imgs.append(conc)\n",
    "                                                    train_x_delta.append(x_delta)\n",
    "                                                    train_y_delta.append(y_delta)\n",
    "                                                    # mark the image as matched\n",
    "                                                    train_lbls.append([0,1])\n",
    "                                        if (\"right\" in adj.keys()): # enough to check only the adj as the right of the curr is the adj\n",
    "                                            for i in range(5, 101, 5):\n",
    "                                                curr1 = pad_right(curr, curr[\"right\"],i) # essentially the adj\n",
    "                                                adj1 = pad_right(adj, adj[\"right\"],i)\n",
    "                                                if (gap == 0 or (\"left\" in curr.keys() and \"right\" in adj.keys())):\n",
    "                                                    if (gap == 0):\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes(curr1, adj1, slice_size)\n",
    "                                                    else:\n",
    "                                                        curr1Left = pad_right(curr[\"left\"], ZERO_CUBE, i) # FIXIT? + assuming the gap will not be more than 150\n",
    "                                                        adj1Right = pad_right(adj[\"right\"], ZERO_CUBE, i) # FIXIT? + assuming the gap will not be more than 150\n",
    "                                                        conc, x_delta, y_delta = concatenate_cubes_with_gap(curr1, adj1, gap, curr1Left, adj1Right, slice_size)\n",
    "                                                    train_imgs.append(conc)\n",
    "                                                    train_x_delta.append(x_delta)\n",
    "                                                    train_y_delta.append(y_delta)\n",
    "                                                    # mark the image as matched\n",
    "                                                    train_lbls.append([0,1])\n",
    "                                    else:\n",
    "                                        # mark the image as not matched\n",
    "                                        train_lbls.append([1,0])\n",
    "                                \n",
    "    warnings.filterwarnings(\"default\")\n",
    "    return train_imgs, train_lbls, train_x_delta, train_y_delta"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# RUN Define model util functions\n",
    "\n",
    "# initialize a shaped matrix of weights with random values\n",
    "def weight_variable(shape):\n",
    "  initial = tf.truncated_normal(shape, stddev=0.1)\n",
    "  return tf.Variable(initial)\n",
    "\n",
    "# initialize a shaped matrix of bias with random values\n",
    "def bias_variable(shape):\n",
    "  initial = tf.constant(0.1, shape=shape)\n",
    "  return tf.Variable(initial)\n",
    "\n",
    "def conv2d(x, W):\n",
    "  return tf.nn.conv2d(x, W, strides=[1, 1, 1, 1], padding='SAME')\n",
    "\n",
    "def max_pool_2x2(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 2, 2, 1],\n",
    "                        strides=[1, 2, 2, 1], padding='SAME')\n",
    "\n",
    "def max_pool_1x2(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 1, 2, 1],\n",
    "                        strides=[1, 1, 2, 1], padding='SAME')\n",
    "\n",
    "def max_pool_2x1(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 2, 1, 1],\n",
    "                        strides=[1, 2, 1, 1], padding='SAME')\n",
    "\n",
    "def max_pool_1x1(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 1, 1, 1],\n",
    "                        strides=[1, 1, 1, 1], padding='SAME')\n",
    "\n",
    "def max_pool_5x5(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 5, 5, 1],\n",
    "                        strides=[1, 5, 5, 1], padding='SAME')\n",
    "\n",
    "def max_pool_5x2(x):\n",
    "  return tf.nn.max_pool(x, ksize=[1, 5, 2, 1],\n",
    "                        strides=[1, 5, 2, 1], padding='SAME')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# RUN Image utility functions (external source)\n",
    "def branchedPoints(skel):\n",
    "    branch1=np.array([[2, 1, 2], [1, 1, 1], [2, 2, 2]])\n",
    "    branch2=np.array([[1, 2, 1], [2, 1, 2], [1, 2, 1]])\n",
    "    branch3=np.array([[1, 2, 1], [2, 1, 2], [1, 2, 2]])\n",
    "    branch4=np.array([[2, 1, 2], [1, 1, 2], [2, 1, 2]])\n",
    "    branch5=np.array([[1, 2, 2], [2, 1, 2], [1, 2, 1]])\n",
    "    branch6=np.array([[2, 2, 2], [1, 1, 1], [2, 1, 2]])\n",
    "    branch7=np.array([[2, 2, 1], [2, 1, 2], [1, 2, 1]])\n",
    "    branch8=np.array([[2, 1, 2], [2, 1, 1], [2, 1, 2]])\n",
    "    branch9=np.array([[1, 2, 1], [2, 1, 2], [2, 2, 1]])\n",
    "    br1=mh.morph.hitmiss(skel,branch1)\n",
    "    br2=mh.morph.hitmiss(skel,branch2)\n",
    "    br3=mh.morph.hitmiss(skel,branch3)\n",
    "    br4=mh.morph.hitmiss(skel,branch4)\n",
    "    br5=mh.morph.hitmiss(skel,branch5)\n",
    "    br6=mh.morph.hitmiss(skel,branch6)\n",
    "    br7=mh.morph.hitmiss(skel,branch7)\n",
    "    br8=mh.morph.hitmiss(skel,branch8)\n",
    "    br9=mh.morph.hitmiss(skel,branch9)\n",
    "    return br1+br2+br3+br4+br5+br6+br7+br8+br9\n",
    "\n",
    "def endPoints(skel):\n",
    "    endpoint1=np.array([[0, 0, 0],\n",
    "                        [0, 1, 0],\n",
    "                        [2, 1, 2]])\n",
    "    \n",
    "    endpoint2=np.array([[0, 0, 0],\n",
    "                        [0, 1, 2],\n",
    "                        [0, 2, 1]])\n",
    "    \n",
    "    endpoint3=np.array([[0, 0, 2],\n",
    "                        [0, 1, 1],\n",
    "                        [0, 0, 2]])\n",
    "    \n",
    "    endpoint4=np.array([[0, 2, 1],\n",
    "                        [0, 1, 2],\n",
    "                        [0, 0, 0]])\n",
    "    \n",
    "    endpoint5=np.array([[2, 1, 2],\n",
    "                        [0, 1, 0],\n",
    "                        [0, 0, 0]])\n",
    "    \n",
    "    endpoint6=np.array([[1, 2, 0],\n",
    "                        [2, 1, 0],\n",
    "                        [0, 0, 0]])\n",
    "    \n",
    "    endpoint7=np.array([[2, 0, 0],\n",
    "                        [1, 1, 0],\n",
    "                        [2, 0, 0]])\n",
    "    \n",
    "    endpoint8=np.array([[0, 0, 0],\n",
    "                        [2, 1, 0],\n",
    "                        [1, 2, 0]])\n",
    "    \n",
    "    ep1=mh.morph.hitmiss(skel,endpoint1)\n",
    "    ep2=mh.morph.hitmiss(skel,endpoint2)\n",
    "    ep3=mh.morph.hitmiss(skel,endpoint3)\n",
    "    ep4=mh.morph.hitmiss(skel,endpoint4)\n",
    "    ep5=mh.morph.hitmiss(skel,endpoint5)\n",
    "    ep6=mh.morph.hitmiss(skel,endpoint6)\n",
    "    ep7=mh.morph.hitmiss(skel,endpoint7)\n",
    "    ep8=mh.morph.hitmiss(skel,endpoint8)\n",
    "    ep = ep1+ep2+ep3+ep4+ep5+ep6+ep7+ep8\n",
    "    return ep\n",
    "\n",
    "def pruning(skeleton, size):\n",
    "    '''remove iteratively end points \"size\" \n",
    "       times from the skeleton\n",
    "    '''\n",
    "    for i in range(0, size):\n",
    "        endpoints = endPoints(skeleton)\n",
    "        endpoints = np.logical_not(endpoints)\n",
    "        skeleton = np.logical_and(skeleton,endpoints)\n",
    "    return skeleton\n",
    "\n",
    "def plot_comparison(original, filtered, filter_name):\n",
    "\n",
    "    fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=(8, 4), sharex=True, sharey=True)\n",
    "    ax1.imshow(original, cmap=plt.cm.gray)\n",
    "    ax1.set_title('original')\n",
    "    ax1.axis('off')\n",
    "    ax1.set_adjustable('box-forced')\n",
    "    ax2.imshow(filtered, cmap=plt.cm.gray)\n",
    "    ax2.set_title(filter_name)\n",
    "    ax2.axis('off')\n",
    "    ax2.set_adjustable('box-forced')\n",
    "    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# RUN Define the model - 250, 125, 62, 25\n",
    "def model_tf(input_width): \n",
    "    global accuracy, correct_prediction, train_step, x, y_, y_conv, keep_prob\n",
    "    \n",
    "    # foundation of the model - the input layer of the image 250 x input_width*2\n",
    "    x = tf.placeholder(tf.float32, [None, 250, input_width*2])\n",
    "    x_image = tf.reshape(x, [-1,250,input_width*2,1]) # 1 is the number of color channels\n",
    "\n",
    "    # the target digits of the model\n",
    "    y_ = tf.placeholder(tf.float32, [None, 2])\n",
    "\n",
    "    # first convolutional layer: one input image and 32 output filters of 5x5\n",
    "    W_conv1 = weight_variable([5, 5, 1, 32])\n",
    "    b_conv1 = bias_variable([32])\n",
    "    h_conv1 = tf.nn.relu(conv2d(x_image, W_conv1) + b_conv1)\n",
    "    if (input_width == 250):\n",
    "        h_pool1 = max_pool_2x2(h_conv1) # size is reduced to 125x250\n",
    "    elif (input_width == 125):\n",
    "        h_pool1 = max_pool_2x1(h_conv1) # size is reduced to 125x250\n",
    "    elif (input_width == 62):\n",
    "        h_pool1 = max_pool_2x1(h_conv1) # size is reduced to 125x125\n",
    "    elif (input_width == 25):\n",
    "        h_pool1 = max_pool_2x1(h_conv1) # size is reduced to 125x50\n",
    "    else:\n",
    "        print(\"ERROR - unsupported slice width\")\n",
    "        return\n",
    "\n",
    "    # second convolutional layer: 32 input (filtered) images and 32 output filters of 5x5\n",
    "    W_conv2 = weight_variable([5, 5, 32, 32])\n",
    "    b_conv2 = bias_variable([32])\n",
    "    h_conv2 = tf.nn.relu(conv2d(h_pool1, W_conv2) + b_conv2)\n",
    "    if (input_width == 62):\n",
    "        h_pool2 = max_pool_1x1(h_conv2) # size is reduced to 125x125\n",
    "    elif (input_width == 25):\n",
    "        h_pool2 = max_pool_1x1(h_conv2) # size is reduced to 125x50\n",
    "    else:\n",
    "        h_pool2 = max_pool_1x2(h_conv2) # size is reduced to 125x125\n",
    "\n",
    "\n",
    "    # third convolutional layer: 32 input (filtered) images and 32 output filters of 5x5\n",
    "    W_conv3 = weight_variable([5, 5, 32, 32])\n",
    "    b_conv3 = bias_variable([32])\n",
    "    h_conv3 = tf.nn.relu(conv2d(h_pool2, W_conv3) + b_conv3)\n",
    "    if (input_width == 25):\n",
    "        h_pool3 = max_pool_5x2(h_conv3) # size is reduced to 25x25\n",
    "    else:\n",
    "        h_pool3 = max_pool_5x5(h_conv3) # size is reduced to 25x25\n",
    "\n",
    "\n",
    "    h_pool3_flat = tf.reshape(h_pool3, [-1, 25*25*32]) # shape as an array \n",
    "\n",
    "    # fourth layer - fully connected with input 25*25*128 and output 1024\n",
    "    W_fc1 = weight_variable([25*25*32, 1024])\n",
    "    b_fc1 = bias_variable([1024])\n",
    "    h_fc1 = tf.nn.relu(tf.matmul(h_pool3_flat, W_fc1) + b_fc1)\n",
    "\n",
    "    # a drop layer with probability \n",
    "    keep_prob = tf.placeholder(tf.float32)\n",
    "    h_fc1_drop = tf.nn.dropout(h_fc1, keep_prob)\n",
    "\n",
    "    # final layer - softmax reduction 2 outputs\n",
    "    W_fc2 = weight_variable([1024, 2])\n",
    "    b_fc2 = bias_variable([2])\n",
    "    y_conv = tf.matmul(h_fc1_drop, W_fc2) + b_fc2\n",
    "\n",
    "    # minimize loss function\n",
    "    cross_entropy = tf.reduce_mean(\n",
    "      tf.nn.softmax_cross_entropy_with_logits(labels=y_, logits=y_conv))\n",
    "\n",
    "    # define train step and rate\n",
    "    train_step = tf.train.AdamOptimizer(1e-4).minimize(cross_entropy)\n",
    "\n",
    "    # evaluate the prediction and the accuracy on the train test - needed only for printing during the training\n",
    "    correct_prediction = tf.equal(tf.argmax(y_conv,1), tf.argmax(y_,1))\n",
    "    accuracy = tf.reduce_mean(tf.cast(correct_prediction, tf.float32))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def train(file, model, train_gap, slice_size):\n",
    "    print(\"#####################################################################\")\n",
    "    print(\"TRAINING:\"+file)\n",
    "    print(\"MODEL:\"+model)\n",
    "    print(\"#####################################################################\")\n",
    "\n",
    "    # TRAIN Read the blue wave-length image and show it\n",
    "    blue = read_and_crop(file) \n",
    "\n",
    "    # Show the relevant image\n",
    "    plt.imshow(blue)\n",
    "\n",
    "    cubes, cols, rows = slice_to_static_slices(blue)\n",
    "\n",
    "    # TRAIN Build a labeled images list (match/no-match) corresponding the the labels list\n",
    "    train_imgs, train_lbls, train_x_delta, train_y_delta = build_train_set(rows, cols, cubes, train_gap, slice_size)\n",
    "    len(train_imgs)\n",
    "\n",
    "    # TRAIN split the original set randomly into a test set\n",
    "    from random import randrange\n",
    "\n",
    "    test_imgs = []\n",
    "    test_lbls = []\n",
    "    len_test = int(len(train_imgs) / 10)\n",
    "    for _ in range(0, len_test):\n",
    "            random_index = randrange(0,len(train_imgs))\n",
    "            test_imgs.append(train_imgs.pop(random_index))\n",
    "            test_lbls.append(train_lbls.pop(random_index))\n",
    "\n",
    "    # train_set = DataSet(np.array(train_imgs), np.array(train_lbls))     \n",
    "    len(test_imgs)\n",
    "\n",
    "    # TRAIN Prepare the session\n",
    "\n",
    "    # create a saver object\n",
    "    saver = tf.train.Saver()\n",
    "\n",
    "    # start session and initialize variables\n",
    "    sess = tf.InteractiveSession()\n",
    "    sess.run(tf.initialize_all_variables())\n",
    "\n",
    "    # TRAIN Train the model\n",
    "    x_batch = []\n",
    "    y_batch = []\n",
    "    # run the train batches\n",
    "    for i in range(400):\n",
    "        x_batch = []\n",
    "        y_batch = []\n",
    "        for _ in range(50):\n",
    "            random_index = randrange(0,len(train_imgs))\n",
    "            x_batch.append(train_imgs[random_index])\n",
    "            y_batch.append(train_lbls[random_index])\n",
    "\n",
    "        # print the accuracy thus far\n",
    "        if (i+1)%50 == 0:\n",
    "            train_accuracy = accuracy.eval(feed_dict={\n",
    "                x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "            print(\"step %d, training accuracy %g\"%(i, train_accuracy))\n",
    "        # train\n",
    "        print(\"step %d\"%(i))\n",
    "        train_step.run(feed_dict={x: x_batch, y_: y_batch, keep_prob: 0.5})\n",
    "\n",
    "    train_accuracy = accuracy.eval(feed_dict={\n",
    "        x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "    print(\"step %d, training accuracy %g\"%(i, train_accuracy))\n",
    "\n",
    "    # Save the variables to disk.\n",
    "    save_path = saver.save(sess, model)\n",
    "    print(\"Model saved in file: %s\" % save_path)\n",
    "\n",
    "    # Close the Session when we're done. If un-commented - need to run next bock of restore...\n",
    "    sess.close()   \n",
    "    print(\"#####################################################################\")\n",
    "    print(\"TRIANING ENDED\")\n",
    "    print(\"#####################################################################\")\n",
    "    print(\" \")\n",
    "    print(\" \")    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def train_all_gaps(file, model, max_train_gap, slice_size):\n",
    "    print(\"#####################################################################\")\n",
    "    print(\"TRAINING:\"+file)\n",
    "    print(\"MODEL:\"+model)\n",
    "    print(\"#####################################################################\")\n",
    "\n",
    "    # TRAIN Read the blue wave-length image and show it\n",
    "    blue = read_and_crop(file) \n",
    "\n",
    "    # Show the relevant image\n",
    "    plt.imshow(blue)\n",
    "\n",
    "    cubes, cols, rows = slice_to_static_slices(blue)\n",
    "\n",
    "    # TRAIN Build a labeled images list (match/no-match) corresponding the the labels list\n",
    "    tmp_train_imgs = []\n",
    "    tmp_train_lbls = []\n",
    "    valid_gaps = [0, 5, 10, 20, 40, 80, 120]\n",
    "    import pdb; pdb.set_trace()\n",
    "    for gap in valid_gaps:\n",
    "        if gap <= max_train_gap:\n",
    "            tmp_train_imgs, tmp_train_lbls, train_x_delta, train_y_delta = build_train_set(rows, cols, cubes, train_gap, slice_size)\n",
    "            train_imgs.append(tmp_train_imgs)\n",
    "            train_lbls.append(tmp_train_lbls)\n",
    "    len(train_imgs)\n",
    "\n",
    "    # TRAIN split the original set randomly into a test set\n",
    "    from random import randrange\n",
    "\n",
    "    test_imgs = []\n",
    "    test_lbls = []\n",
    "    len_test = int(len(train_imgs) / 10)\n",
    "    for _ in range(0, len_test):\n",
    "            random_index = randrange(0,len(train_imgs))\n",
    "            test_imgs.append(train_imgs.pop(random_index))\n",
    "            test_lbls.append(train_lbls.pop(random_index))\n",
    "\n",
    "    # train_set = DataSet(np.array(train_imgs), np.array(train_lbls))     \n",
    "    len(test_imgs)\n",
    "\n",
    "    # TRAIN Prepare the session\n",
    "\n",
    "    # create a saver object\n",
    "    saver = tf.train.Saver()\n",
    "\n",
    "    # start session and initialize variables\n",
    "    sess = tf.InteractiveSession()\n",
    "    sess.run(tf.initialize_all_variables())\n",
    "\n",
    "    # TRAIN Train the model\n",
    "    x_batch = []\n",
    "    y_batch = []\n",
    "    # run the train batches\n",
    "    for i in range(500):\n",
    "        x_batch = []\n",
    "        y_batch = []\n",
    "        for _ in range(100):\n",
    "            random_index = randrange(0,len(train_imgs))\n",
    "            x_batch.append(train_imgs[random_index]) # check: train_imgs.pop(random_index)\n",
    "            y_batch.append(train_lbls[random_index])\n",
    "\n",
    "        # print the accuracy thus far\n",
    "        if (i+1)%50 == 0:\n",
    "            train_accuracy = accuracy.eval(feed_dict={\n",
    "                x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "            print(\"step %d, training accuracy %g\"%(i, train_accuracy))\n",
    "        # train\n",
    "        print(\"step %d\"%(i))\n",
    "        train_step.run(feed_dict={x: x_batch, y_: y_batch, keep_prob: 0.5})\n",
    "\n",
    "    train_accuracy = accuracy.eval(feed_dict={\n",
    "        x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "    print(\"step %d, training accuracy %g\"%(i, train_accuracy))\n",
    "\n",
    "    # Save the variables to disk.\n",
    "    save_path = saver.save(sess, model)\n",
    "    print(\"Model saved in file: %s\" % save_path)\n",
    "\n",
    "    # Close the Session when we're done. If un-commented - need to run next bock of restore...\n",
    "    sess.close()   \n",
    "    print(\"#####################################################################\")\n",
    "    print(\"TRIANING ENDED\")\n",
    "    print(\"#####################################################################\")\n",
    "    print(\" \")\n",
    "    print(\" \")    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "def validate(file, model, gap, slice_size):\n",
    "    print(\"#####################################################################\")\n",
    "    print(\"VALIDATING:\"+file)\n",
    "    print(\"GAP:\"+str(gap))\n",
    "    print(\"#####################################################################\")\n",
    "    \n",
    "    # VALIDATE Read the blue wave-length VALIDATION image and show it\n",
    "    blue = read_and_crop(file) \n",
    "\n",
    "    blue = blue[0:2500,0:3500]\n",
    "\n",
    "    # Show the relevant image\n",
    "    plt.imshow(blue)\n",
    "\n",
    "    cubes, cols, rows = slice_to_static_slices(blue)\n",
    "\n",
    "    # VALIDATE prepare the data sets\n",
    "    test_imgs, test_lbls, test_x_delta, test_y_delta = build_train_set(rows, cols, cubes, gap, slice_size)\n",
    "    len(test_imgs)\n",
    "\n",
    "    # VALIDATE Prepare a test session \n",
    "\n",
    "    # Add ops to save and restore all the variables.\n",
    "    saver = tf.train.Saver()\n",
    "\n",
    "    # start session and initialize variables\n",
    "    sess = tf.InteractiveSession()\n",
    "\n",
    "    # Restore variables from disk.\n",
    "    saver.restore(sess, model)\n",
    "    print(\"Model restored.\")\n",
    "\n",
    "    # VALIDATE Validate the model\n",
    "\n",
    "    # import pdb; pdb.set_trace()\n",
    "    se = 0\n",
    "    st = 0\n",
    "    v1t = []\n",
    "    v2t = []\n",
    "    v1tt = []\n",
    "    v2tt = []\n",
    "    false_x_deltas = []\n",
    "    false_y_deltas = []\n",
    "    count = 0\n",
    "    length = 10\n",
    "    batch = 100\n",
    "    end = (int(len(test_imgs) / batch)) - length\n",
    "    # change the ranges in the loop below - first number is the start point (multiplied by batch size)\n",
    "    # second number is the end point (multiplied by batch size)\n",
    "    # third number is the jump from batch to batch\n",
    "    # use the length about to set the batch length\n",
    "    for start in range(0, end, 20):\n",
    "        # print(\"Validating start at #%d end at %d\"%(start*batch,(start+length)*batch))\n",
    "        for i in range(start, start+length):\n",
    "            x_batch = []\n",
    "            y_batch = []\n",
    "            row_delta_batch = []\n",
    "            col_delta_batch = []\n",
    "            for j in range(batch):\n",
    "                x_batch.append(test_imgs[i*batch+j])\n",
    "                y_batch.append(test_lbls[i*batch+j])\n",
    "                row_delta_batch.append(test_x_delta[i*batch+j])\n",
    "                col_delta_batch.append(test_y_delta[i*batch+j])\n",
    "\n",
    "\n",
    "        #     my_prediction=tf.argmax(y_conv,1)\n",
    "        #     my_target=tf.argmax(y_,1)\n",
    "        #     print(my_prediction.eval(feed_dict={x:x_batch, y_: y_batch, keep_prob: 1.0}))\n",
    "        #     print(my_target.eval(feed_dict={x:x_batch, y_: y_batch, keep_prob: 1.0}))\n",
    "\n",
    "        # # print the accuracy thus far\n",
    "        #     test_accuracy = accuracy.eval(feed_dict={\n",
    "        #         x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "        # #     if (i+1)%10 == 0:\n",
    "        #     print(\"step %d, validation accuracy %f\"%(i, test_accuracy))\n",
    "\n",
    "            my_prediction=tf.argmax(y_conv,1)\n",
    "            my_target=tf.argmax(y_,1)\n",
    "            v1 = my_prediction.eval(feed_dict={x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "            v2 = my_target.eval(feed_dict={x:x_batch, y_: y_batch, keep_prob: 1.0})\n",
    "            v1t = np.concatenate((v1t, v1), axis=0)\n",
    "            v2t = np.concatenate((v2t, v2), axis=0)\n",
    "\n",
    "            c1 = np.sum(np.absolute(np.subtract(v2, v1)))\n",
    "            c2 = np.sum(np.absolute(v2))\n",
    "            se += c1\n",
    "            st += c2\n",
    "\n",
    "            false_guess_indexes = np.nonzero((v1-v2).clip(min=0))\n",
    "            A = np.array(row_delta_batch)\n",
    "            false_x_deltas = np.concatenate((false_x_deltas, A[false_guess_indexes]), axis=0)\n",
    "            B = np.array(col_delta_batch)\n",
    "            false_y_deltas = np.concatenate((false_y_deltas, B[false_guess_indexes]), axis=0)\n",
    "\n",
    "        # print(\"predict=\", v1t)\n",
    "        # print(\"target =\", v2t)\n",
    "        count += ((i+1)*batch - start*batch)\n",
    "        precision, recall, f_score, support = precision_recall_fscore_support(v2t, v1t, average='binary')\n",
    "        print(\"step %d-%d, precision %f, recall %f, f_score %f\"%(start, i, precision, recall, f_score))\n",
    "        # print(\"Accumulated total true = %d\"%(st));\n",
    "        # print(\"Accumulated total error rate = %f\"%(se/count));\n",
    "        v1tt = np.concatenate((v1tt, v1t), axis=0)\n",
    "        v2tt = np.concatenate((v2tt, v2t), axis=0)\n",
    "\n",
    "    precision, recall, f_score, support = precision_recall_fscore_support(v2tt, v1tt, average='binary')\n",
    "    print(\"TOTAL %d, precision %f, recall %f, f_score %f\"%(count, precision, recall, f_score))\n",
    "    print(\"TOTAL true = %d\"%(st));\n",
    "    print(\"TOTAL error rate = %f\"%(se/count));\n",
    "\n",
    "    with open('results.csv', 'a') as csvfile:\n",
    "        csvout = csv.writer(csvfile)\n",
    "        csvout.writerow([file, model, gap, slice_size, count, precision, recall, f_score, st, se/count])\n",
    "\n",
    "    plt.plot(false_x_deltas, false_y_deltas, \"o\")\n",
    "\n",
    "    # Close the Session when we're done.\n",
    "    sess.close()   \n",
    "    print(\"#####################################################################\")\n",
    "    print(\"VALIDATION ENDED\")\n",
    "    print(\"#####################################################################\")\n",
    "    print(\" \")\n",
    "    print(\" \")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "def run_all(file_a, file_b, model, train_gap, slice_size):\n",
    "    model_tf(slice_size)\n",
    "    train_all_gaps(file_a, model, train_gap, slice_size)\n",
    "    validate(file_b, model, 0, slice_size)\n",
    "    validate(file_b, model, 5, slice_size)\n",
    "    validate(file_b, model, 10, slice_size)\n",
    "    validate(file_b, model, 20, slice_size)\n",
    "    validate(file_b, model, 40, slice_size)\n",
    "    validate(file_b, model, 80, slice_size)\n",
    "    validate(file_b, model, 120, slice_size)    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": [],
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "#####################################################################\n",
      "TRAINING:PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\n",
      "MODEL:/Users/il239838/Development/models/ML924__012_160gap_250px_model400x50.ckpt\n",
      "#####################################################################\n",
      "2264 4999\n",
      "> <ipython-input-7-d242125dc537>(20)train_all_gaps()\n",
      "-> for gap in valid_gaps:\n",
      "(Pdb) gap\n",
      "*** NameError: name 'gap' is not defined\n"
     ]
    }
   ],
   "source": [
    "# RUN sequence of files with 250px\n",
    "\n",
    "with open('results.csv', 'a') as csvfile:\n",
    "    csvout = csv.writer(csvfile)\n",
    "    csvout.writerow([])\n",
    "\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_0gap_250px_model400x50.ckpt\", 0, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_5gap_250px_model400x50.ckpt\", 5, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_10gap_250px_model400x50.ckpt\", 10, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_20gap_250px_model400x50.ckpt\", 20, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_40gap_250px_model400x50.ckpt\", 40, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_80gap_250px_model400x50.ckpt\", 80, 250)\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_120gap_250px_model400x50.ckpt\", 120, 250)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_160gap_250px_model400x50.ckpt\", 160, 250)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# RUN sequence of files with 125px\n",
    "\n",
    "with open('results.csv', 'a') as csvfile:\n",
    "    csvout = csv.writer(csvfile)\n",
    "    csvout.writerow([])\n",
    "\n",
    "# run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112459-ML445__001.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120502-ML445__001.jpg\", \"/Users/il239838/Development/models/ML445__001_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112503-ML475__002.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120506-ML475__002.jpg\", \"/Users/il239838/Development/models/ML475__002_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112508-ML499__003.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120510-ML499__003.jpg\", \"/Users/il239838/Development/models/ML499__003_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112512-ML540__004.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120515-ML540__004.jpg\", \"/Users/il239838/Development/models/ML540__004_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112516-ML595__005.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120519-ML595__005.jpg\", \"/Users/il239838/Development/models/ML595__005_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112520-ML638__006.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120522-ML638__006.jpg\", \"/Users/il239838/Development/models/ML638__006_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112526-ML656__007.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120529-ML656__007.jpg\", \"/Users/il239838/Development/models/ML656__007_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112533-ML706__008.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120535-ML706__008.jpg\", \"/Users/il239838/Development/models/ML706__008_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112539-ML728__009.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120541-ML728__009.jpg\", \"/Users/il239838/Development/models/ML728__009_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112545-ML772__010.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120548-ML772__010.jpg\", \"/Users/il239838/Development/models/ML772__010_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112550-ML858__011.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120553-ML858__011.jpg\", \"/Users/il239838/Development/models/ML858__011_125px_model400x50.ckpt\", 125)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_125px_model400x50.ckpt\", 125)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "code_folding": [
     0
    ],
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# RUN sequence of files with 62px\n",
    "\n",
    "with open('results.csv', 'a') as csvfile:\n",
    "    csvout = csv.writer(csvfile)\n",
    "    csvout.writerow([])\n",
    "\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112459-ML445__001.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120502-ML445__001.jpg\", \"/Users/il239838/Development/models/ML445__001_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112503-ML475__002.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120506-ML475__002.jpg\", \"/Users/il239838/Development/models/ML475__002_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112508-ML499__003.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120510-ML499__003.jpg\", \"/Users/il239838/Development/models/ML499__003_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112512-ML540__004.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120515-ML540__004.jpg\", \"/Users/il239838/Development/models/ML540__004_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112516-ML595__005.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120519-ML595__005.jpg\", \"/Users/il239838/Development/models/ML595__005_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112520-ML638__006.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120522-ML638__006.jpg\", \"/Users/il239838/Development/models/ML638__006_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112526-ML656__007.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120529-ML656__007.jpg\", \"/Users/il239838/Development/models/ML656__007_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112533-ML706__008.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120535-ML706__008.jpg\", \"/Users/il239838/Development/models/ML706__008_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112539-ML728__009.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120541-ML728__009.jpg\", \"/Users/il239838/Development/models/ML728__009_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112545-ML772__010.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120548-ML772__010.jpg\", \"/Users/il239838/Development/models/ML772__010_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112550-ML858__011.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120553-ML858__011.jpg\", \"/Users/il239838/Development/models/ML858__011_62px_model400x50.ckpt\", 62)\n",
    "run_all(\"PX303/FG001/PX303-Fg001-V-C01-R01-D05032015-T112602-ML924__012.jpg\", \"PX303/FG006/PX303-Fg006-V-C01-R01-D08032015-T120605-ML924__012.jpg\", \"/Users/il239838/Development/models/ML924__012_62px_model400x50.ckpt\", 62)"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Raw Cell Format",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "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.5.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}