warp.ipynb

{
 "metadata": {
  "name": "",
  "signature": "sha256:d22620ed99bcbf349dc5bf8f7e825624118898b912a9a96c02d63037a1b862bf"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Warp using random acceleration\n",
      "\n",
      "This notebook demonstrates applying a double integration of an acceleration to get a warping displacement.\n",
      "\n",
      "If you have a random per-pixel displacement, at best you will get some jitter and at worst it will just blur the image. If you have a random velocity, there will be little consistency to the amount of warping in adjacent pixels. If you have random acceleration, the local warping will be consistent.\n",
      "\n",
      "The accumulation is a random walk. This means that by the time you get to the far side of the image, you may have accumulated quite a velocity and there may be too much warp, so you need some hysteresis to damp any long-term acceleration."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from pylab import imread, imshow, rcParams, subplot\n",
      "from IPython.html.widgets import interact\n",
      "from scipy.ndimage import map_coordinates\n",
      "from scipy.ndimage import gaussian_filter\n",
      "import numpy as np\n",
      "\n",
      "rcParams['figure.figsize'] = (18, 9)\n",
      "\n",
      "img = imread(\"/home/leek/coxfs01/data/ac3/images/ac3_input_0000.tif\")\n",
      "def draw_warped(scalex1000, hysteresis_pct, sigma):\n",
      "    np.random.seed(1234)\n",
      "    scale = float(scalex1000) / 1000\n",
      "    hysteresis = float(hysteresis_pct) / 100\n",
      "    displacement_d2x = np.random.normal(size=img.shape[:2], scale=scale)\n",
      "    displacement_d2y = np.random.normal(size=img.shape[:2], scale=scale)\n",
      "    def integrate_with_hysteresis(dx, hysteresis):\n",
      "        x = np.zeros(dx.shape, dx.dtype)\n",
      "        x[0] = dx[0]\n",
      "        x[:, 0] = dx[:, 0]\n",
      "        for i in range(1, dx.shape[0]):\n",
      "            x[i] = x[i-1] * hysteresis + dx[i]\n",
      "        for j in range(1, dx.shape[1]):\n",
      "            x[:, j] = x[:, j-1] * hysteresis + dx[:, j]\n",
      "        return x\n",
      "    displacement_dx = integrate_with_hysteresis(displacement_d2x, hysteresis * hysteresis)\n",
      "    displacement_dy = integrate_with_hysteresis(displacement_d2y, hysteresis)\n",
      "    displacement_x = gaussian_filter(np.cumsum(np.cumsum(displacement_dx, 0), 1), sigma)\n",
      "    displacement_y = gaussian_filter(np.cumsum(np.cumsum(displacement_dy, 0), 1), sigma)\n",
      "    y, x = np.mgrid[0:img.shape[0], 0:img.shape[1]]\n",
      "    imgwarped = map_coordinates(img, ((y+displacement_y, x+displacement_x)), mode='reflect')\n",
      "    subplot(1, 2, 1).imshow(img, cmap='gray')\n",
      "    subplot(1, 2, 1).set_title(\"Original\")\n",
      "    subplot(1, 2, 2).imshow(imgwarped, cmap='gray')\n",
      "    subplot(1, 2, 2).set_title(\"Warped\")\n",
      "interact(draw_warped, scalex1000=(1,20), hysteresis_pct=(1, 100), sigma=(0, 10))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "display_data",