change_detection_tests.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import os \n",
    "os.environ[\"ISISROOT\"] = \"/usgs/cpkgs/anaconda3_linux/envs/isis3.9.0\"\n",
    "\n",
    "from pysis import isis\n",
    "from pysis.exceptions import ProcessError\n",
    "\n",
    "import matplotlib\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "from plio.io.io_gdal import GeoDataset \n",
    "import numpy as np\n",
    "\n",
    "from autocnet.examples import get_path\n",
    "from autocnet.cg import change_detection as cd\n",
    "from autocnet.utils.utils import bytescale\n",
    "\n",
    "import matplotlib"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def cartesian_product(*arrays):\n",
    "    la = len(arrays)\n",
    "    dtype = \"object\"\n",
    "    arr = np.empty([len(a) for a in arrays] + [la], dtype=dtype)\n",
    "    for i, a in enumerate(np.ix_(*arrays)):\n",
    "        arr[...,i] = a\n",
    "    return arr.reshape(-1, la)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Set up images for now change \n",
    "nochange_before_geo = GeoDataset(get_path(\"no_change.before.cub\"))\n",
    "nochange_after_geo = GeoDataset(get_path(\"no_change.after.cub\"))\n",
    "\n",
    "nochange_before_arr = nochange_before_geo.read_array()\n",
    "nochange_after_arr = nochange_after_geo.read_array()\n",
    "nochange_before_arr[nochange_before_arr==nochange_before_arr.min()] = np.nan\n",
    "nochange_after_arr[nochange_after_arr==nochange_after_arr.min()] = np.nan\n",
    "\n",
    "\n",
    "boulder_shift_before_geo = GeoDataset(get_path(\"boulder_shift.before.cub\"))\n",
    "boulder_shift_after_geo = GeoDataset(get_path(\"boulder_shift.after.cub\"))\n",
    "\n",
    "boulder_shift_before_arr = boulder_shift_before_geo.read_array()\n",
    "boulder_shift_after_arr = boulder_shift_after_geo.read_array()\n",
    "boulder_shift_before_arr[boulder_shift_before_arr==boulder_shift_before_arr.min()] = np.nan\n",
    "boulder_shift_after_arr[boulder_shift_after_arr==boulder_shift_after_arr.min()] = np.nan"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "9X3 Matrix\n"
     ]
    }
   ],
   "source": [
    "#{\"nfeatures\": 2000, \"scaleFactor\": 1.1, \"nlevels\": 1}\n",
    "# cluster_params={\"min_samples\": 10, \"max_eps\": 10, \"eps\": None, \"xi\":.5}\n",
    "\n",
    "nfeatures = [500, 1000, 2000]\n",
    "scalefactor = [1.1]\n",
    "nlevels = [1] \n",
    "extractor_params = [dict(zip([\"nfeatures\", \"scaleFactor\", \"nlevels\"], ps)) for ps in cartesian_product(nfeatures, scalefactor, nlevels)] \n",
    "\n",
    "\n",
    "min_samples = [10, 15, 20]\n",
    "max_eps = [10, 20, 30]\n",
    "xi = [.5]\n",
    "eps = [.5]\n",
    "cluster_params =  [dict(zip([\"min_samples\", \"max_eps\", \"eps\", \"xi\"], ps)) for ps in cartesian_product(min_samples, max_eps, eps, xi)] \n",
    "\n",
    "print(f\"{len(cluster_params)}X{len(extractor_params)} Matrix\")\n",
    "parameters = cartesian_product(extractor_params, cluster_params)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "boulder_shift_results = [cd.okubogar_detector(boulder_shift_before_geo, boulder_shift_after_geo, nbins=100, extractor_kwargs=params[0]) for params in parameters]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Known Boulder Shift change "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 68,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "nx = len(extractor_params)\n",
    "ny = len(cluster_params)\n",
    "\n",
    "fig, axs = plt.subplots(ny, nx, sharex='col', sharey='row')\n",
    "fig.set_figheight(50)\n",
    "fig.set_figwidth(20)\n",
    "for ax in axs.flat:\n",
    "    ax.set(xlabel='x-label', ylabel='y-label')\n",
    "    \n",
    "xs, ys = np.mgrid[:ny, :nx]\n",
    "cmap = matplotlib.cm.get_cmap('plasma')\n",
    "\n",
    "for i, coord in enumerate(zip(xs.flatten(),ys.flatten())):\n",
    "    x,y = coord\n",
    "    polys, weights, diff = boulder_shift_results[i]\n",
    "    \n",
    "    indices = [i for (i,w) in enumerate(weights) if w > 0]\n",
    "    filtered_polys = []\n",
    "    filtered_polys = [polys[i] for i in indices]\n",
    "    \n",
    "    axs[x,y].set(xlabel=f\"ex {parameters[i][1]}\", ylabel=f\"cl {parameters[i][0]}\")\n",
    "    axs[x,y].imshow(diff, cmap=\"Greys_r\")\n",
    "    axs[x,y].imshow(boulder_shift_after_arr, cmap=\"Greys_r\", alpha=.5)\n",
    "\n",
    "    if weights:\n",
    "        colors = bytescale(np.asarray(weights))\n",
    "        ax = [axs[x,y].fill(*p.exterior.xy, c=\"r\", alpha=.6) for i,p in enumerate(filtered_polys)]\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {
    "collapsed": true,
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "no_change_results = [cd.okubogar_detector(nochange_before_geo, nochange_after_geo, nbins=100, extractor_kwargs=params[0]) for params in parameters]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## No Change Graphs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "nx = len(extractor_params)\n",
    "ny = len(cluster_params)\n",
    "\n",
    "fig, axs = plt.subplots(ny, nx, sharex='col', sharey='row')\n",
    "fig.set_figheight(50)\n",
    "fig.set_figwidth(20)\n",
    "for ax in axs.flat:\n",
    "    ax.set(xlabel='x-label', ylabel='y-label')\n",
    "    \n",
    "xs, ys = np.mgrid[:ny, :nx]\n",
    "cmap = matplotlib.cm.get_cmap('plasma')\n",
    "\n",
    "for i, coord in enumerate(zip(xs.flatten(),ys.flatten())):\n",
    "    x,y = coord\n",
    "    polys, weights, diff = no_change_results[i]\n",
    "    \n",
    "    indices = [i for (i,w) in enumerate(weights) if w > 6]\n",
    "    filtered_polys = []\n",
    "    filtered_polys = [polys[i] for i in indices]\n",
    "    \n",
    "    axs[x,y].set(xlabel=f\"ex {parameters[i][1]}\", ylabel=f\"cl {parameters[i][0]}\")\n",
    "    axs[x,y].imshow(diff, cmap=\"Greys_r\")\n",
    "    axs[x,y].imshow(nochange_after_arr, cmap=\"Greys_r\", alpha=.5)\n",
    "\n",
    "    if weights:\n",
    "        colors = bytescale(np.asarray(weights))\n",
    "        ax = [axs[x,y].fill(*p.exterior.xy, c=\"r\", alpha=.6) for i,p in enumerate(filtered_polys)]\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "nc_weights = []\n",
    "for r in no_change_results:\n",
    "    nc_weights.extend(r[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-4.7621055, 0.2109375, -2.3901317)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.min(nc_weights), np.max(nc_weights), np.mean(nc_weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "bs_weights = []\n",
    "for r in boulder_shift_results:\n",
    "    bs_weights.extend(r[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(1, figsize=(10,5))\n",
    "_ = plt.hist(nc_weights, color=\"orange\", alpha=.5, label=\"no change\")\n",
    "plt.legend(loc='upper right')\n",
    "plt.figure(2, figsize=(10,5))\n",
    "_ = plt.hist(bs_weights, alpha=.5, label=\"change\")\n",
    "plt.legend(loc='upper right')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(1, figsize=(10,5))\n",
    "_ = plt.hist(boulder_shift_results[3][2].flatten(), bins=80, alpha=.6, label=\"change\")\n",
    "_ = plt.hist(no_change_results[3][2].flatten(), bins=40, alpha=.6, label=\"no change\")\n",
    "plt.legend(loc='upper right')\n",
    "plt.show()\n",
    "\n",
    "plt.figure(1, figsize=(10,5))\n",
    "_ = plt.hist(np.append(no_change_results[3][2].flatten(), boulder_shift_results[3][2].flatten()), bins=90, alpha=.5, label=\"change\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [],
   "source": [
    "# np.append(nc_weights, bs_weights)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-67.0, 43.0, -9.827137)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.min(boulder_shift_results[0][2]), np.max(boulder_shift_results[0][2]), np.mean(boulder_shift_results[0][2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-25.0, 16.0, -2.4594834)"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.min(no_change_results[0][2]), np.max(no_change_results[0][2]), np.mean(no_change_results[0][2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "AutocnetDev",
   "language": "python",
   "name": "autocnetdev"
  },
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}