{
 "metadata": {
  "name": "",
  "signature": "sha256:ef6f6afa6b5b8a432f9e74f0c9a65297e115964756970f777728405a6c9adf60"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%matplotlib inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We start by initializing the Earth Engine API."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Load code to setup authorization for an IPython Notebook. Note that this is a temporary step\n",
      "# that is required until the Earth Engine Python API is updated to include this logic.\n",
      "%run 'authorize_earth_engine_in_notebook.ipynb'\n",
      "\n",
      "# Initialize Earth Engine\n",
      "# Note that we are calling a function defined in the previously run IPython Notebook, rather than\n",
      "# the typical call to ee.Initialize()\n",
      "ee_initialize()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Then we import in other Python packages we need."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import datetime\n",
      "from matplotlib import dates\n",
      "import matplotlib.dates as mdates\n",
      "from pylab import *"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next we define the Landsat bands that we would like to plot, along with the starting and ending times."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "xBand = 'time'\n",
      "yBandList = [\n",
      "        '10',\n",
      "        u'20',\n",
      "        u'30',\n",
      "        u'40',\n",
      "        u'50',\n",
      "        u'61',\n",
      "        u'62',\n",
      "        u'70',\n",
      "        u'80',\n",
      "    ]\n",
      "startTime = datetime.datetime(2000, 1, 1)\n",
      "endTime = datetime.datetime(2004, 1, 1)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next we contruct a filtered ImageCollection for the date range, and extract band information for a specified point location."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "collection = ee.ImageCollection('L7_L1T').filterDate(startTime, endTime)\n",
      "point = {'type':'Point', 'coordinates':[ -116.88629,36.56122]};  # death valley (should be stable)\n",
      "info = collection.getRegion(point,500).getInfo()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We separate the information returned into column headers and data."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# extract the header column names\n",
      "header = info[0]\n",
      "# create a Numpy array of the data\n",
      "data = array(info[1:])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next we extract time information and convert it to at Python datatime data type."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# extract the time information\n",
      "iTime = header.index('time')\n",
      "# convert to Python datetime objects\n",
      "time = [datetime.datetime.fromtimestamp(i/1000) for i in (data[0:,iTime].astype(int))]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Extract the data columns what we want to display on the plot."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "iBands = [header.index(b) for b in yBandList]\n",
      "yData = data[0:,iBands].astype(np.float)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Calculate NDVI"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "band3 = yData[:,2]\n",
      "band4 = yData[:,3]\n",
      "ndvi = (band4 - band3) / (band4 + band3)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "And finally, we create a plot of MODIS band values as a function of time."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# matplotlib date format object\n",
      "\n",
      "fig = figure(figsize=(12,8), dpi=80)\n",
      "\n",
      "# plot the band values\n",
      "ax1 = fig.add_subplot(211)\n",
      "ax1.plot(time, yData[:,2], 'o', color=\"red\", label=\"Band 3\")\n",
      "ax1.plot(time, yData[:,3], 'o', color=\"magenta\",  label=\"Band 4\")\n",
      "ax1.legend(loc='best')\n",
      "ax1.grid(True)\n",
      "\n",
      "#plt.title('Band values as a function of time')\n",
      "ax1.set_ylabel('Band Values')\n",
      "\n",
      "# plot NDVI\n",
      "ax2 = fig.add_subplot(212, sharex=ax1)\n",
      "ax2.plot(time, ndvi, 'o', color=\"black\", label=\"NDVI\")\n",
      "ax2.grid(True)\n",
      "start, end = ax2.get_xlim()\n",
      "ax2.xaxis.set_ticks(np.arange(start, end, 64.5))\n",
      "\n",
      "# Format the ticks.\n",
      "years    = mdates.YearLocator()   # every year\n",
      "months   = mdates.MonthLocator()  # every month\n",
      "yearsFmt = mdates.DateFormatter('%Y')\n",
      "\n",
      "ax2.set_ylabel('NDVI')\n",
      "\n",
      "ax2.xaxis.set_major_locator(years)\n",
      "ax2.xaxis.set_major_formatter(yearsFmt)\n",
      "ax2.xaxis.set_minor_locator(months)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# Convert the timestamp to a numpy array\n",
      "t = np.array([i.toordinal() for i in time])\n",
      "t"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "$$\n",
      "\\begin{bmatrix}\n",
      "t_1 & 1 \\\\\\\n",
      "t_2 & 1 \\\\\\\n",
      "\\vdots & \\vdots \\\\\\\n",
      "t_n & 1\n",
      "\\end{bmatrix}\n",
      "\\begin{bmatrix}\n",
      "x_0 \\\\\\\n",
      "x_1 \n",
      "\\end{bmatrix}=\n",
      "\\begin{bmatrix}\n",
      "y_1 \\\\\\\n",
      "y_2 \\\\\\\n",
      "\\vdots \\\\\\\n",
      "y_n\n",
      "\\end{bmatrix}\n",
      "$$\n",
      "\n",
      "$$ \\mathbf{A}\\mathbf{x} = \\mathbf{b} $$\n",
      "\n",
      "$$ \\mathbf{x} = \\mathbf{A} \\backslash \\mathbf{b} $$"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = array([ t, ones(len(t))]).transpose()\n",
      "b = ndvi\n",
      "x = linalg.lstsq(A,b)[0] # obtaining the parameters\n",
      "x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b_hat = A.dot(x)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "fig2 = figure(figsize=(12,4), dpi=80)\n",
      "plot(time,b_hat.transpose(),'r-',t,b,'o')\n",
      "fig2.fmt_xdata = mdates.DateFormatter('%Y-%m-%d')\n",
      "fig2.autofmt_xdate()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": []
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}