IM_MODIS_processing.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os, glob\n",
    "import xuleta as xu\n",
    "import numpy as np\n",
    "import pylab as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### GENERIC MODIS multi-tile processing and mosaicking\n",
    "* Remember to set the no-data value accordingly.\n",
    "* The boundaries are set no Northeastern Brazil\n",
    "\n",
    "#### MCD43 albedo\n",
    "* White-sky shortwave albedo, band 9\n",
    "* Black-sky shortwave albedo, band 19"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "folder = '/media/denis/seagate/PB/MODIS/mod16a2/'\n",
    "os.chdir(folder)\n",
    "\n",
    "os.system('mkdir 1309 1310 1409 1410')\n",
    "os.system('mv *h13v09* 1309; mv *h13v10* 1310; mv *h14v09* 1409; mv *h14v10* 1410')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "os.chdir('/home/denis/apps/MRT/bin')\n",
    "\n",
    "os.system('java -jar MRTBatch.jar -d /media/denis/seagate/PB/MODIS/mod16a2/1309 -p /media/denis/seagate/PB/MODIS/mod16a2/m1309.prm; chmod +x mrtbatch; ./mrtbatch; rm mrtbatch; java -jar MRTBatch.jar -d /media/denis/seagate/PB/MODIS/mod16a2/1310 -p /media/denis/seagate/PB/MODIS/mod16a2/m1310.prm; chmod +x mrtbatch; ./mrtbatch; rm mrtbatch; java -jar MRTBatch.jar -d /media/denis/seagate/PB/MODIS/mod16a2/1409 -p /media/denis/seagate/PB/MODIS/mod16a2/m1409.prm; chmod +x mrtbatch; ./mrtbatch; rm mrtbatch; java -jar MRTBatch.jar -d /media/denis/seagate/PB/MODIS/mod16a2/1410 -p /media/denis/seagate/PB/MODIS/mod16a2/m1410.prm; chmod +x mrtbatch; ./mrtbatch; rm mrtbatch; ')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "for i,j,k,l in zip(sorted(glob.glob(folder+\"1309/prm/*ET*.tif\")),sorted(glob.glob(folder+\"1310/prm/*ET*.tif\")),\n",
    "                   sorted(glob.glob(folder+\"1409/prm/*ET*.tif\")),sorted(glob.glob(folder+\"1410/prm/*ET*.tif\"))):\n",
    "    g1 = 'gdal_merge.py -ul_lr -48.7551 -1.0445 -32.3925 -18.3494 -n 32767 -a_nodata 32767 -co COMPRESS=DEFLATE -co PREDICTOR=2 -co ZLEVEL=9 -of GTiff -o '\n",
    "\n",
    "    go = folder+i.split('.')[-3]+'.tif '\n",
    "    g2 = i+' '+j+' '+k+' '+l\n",
    "    gf = g1+go+g2\n",
    "    os.system(gf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "xu.renamedate(folder=folder,oldf=\"%Y%j\",newf='%Y-%m-%d',wts=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Addressing further issues with no-data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "import xuleta as xu \n",
    "import numpy as np\n",
    "print(folder)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "for i in glob.glob(folder+\"*.tif\"):\n",
    "    a,m = xu.TifToArray(i)\n",
    "    a = np.where(a>30000,32767,a)\n",
    "    xu.ArrayToTif(a,i.split('/')[-1],Folder=folder,Metadata=m,Type=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "for i in glob.glob(folder+\"*.tif\"):\n",
    "    a,m = xu.TifToArray(i)\n",
    "    b = np.where(a>30000,32767,a)\n",
    "    fn = i.split(\"/\")[-1]\n",
    "    xu.ArrayToTif(b,Folder=folder,Filename=fn,Metadata=m,Type=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "os.chdir(folder)\n",
    "os.system('for f in *.tif; do gdalwarp -overwrite -co COMPRESS=DEFLATE -co PREDICTOR=2 -co ZLEVEL=9 -multi -dstnodata 32767 -of GTiff $f t$f; done')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "folder = '/media/denis/seagate/PB/MODIS/myd17a2h/'\n",
    "os.chdir(folder)\n",
    "os.system('rm A*.tif')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "#CLEANING - MAKE SURE IS EVERYTHING OK BEFORE DOING IT \n",
    "os.chdir(folder)\n",
    "os.system('find . -name \"*.hdf\" -type f -delete')\n",
    "os.system('rm -rf 1309 1310 1409 1410')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "xu.renamedate(folder, oldf='%Y%j',wts=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### REFL, NDVI"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "folder = '/media/denis/seagate/PB/MODIS/myd09a1/'"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "os.chdir(folder)\n",
    "\n",
    "os.system('mkdir 1309 1310 1409 1410')\n",
    "os.system('mv *h13v09* 1309; mv *h13v10* 1310; mv *h14v09* 1409; mv *h14v10* 1410')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**RUN ZHAODELONG NDVI script for each tile**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "for i,j,k,l in zip(sorted(glob.glob(folder+\"t1309/ndvi/*tif\")),sorted(glob.glob(folder+\"t1310/ndvi/*tif\")),\n",
    "                   sorted(glob.glob(folder+\"t1409/ndvi/*tif\")),sorted(glob.glob(folder+\"t1410/ndvi/*tif\"))):\n",
    "    g1 = 'gdal_merge.py  -n -2 -a_nodata -2 -co COMPRESS=DEFLATE -co PREDICTOR=2 -co ZLEVEL=9 -of GTiff -o '\n",
    "\n",
    "    go = folder+i.split('.')[-6]+'.tif '\n",
    "    g2 = i+' '+j+' '+k+' '+l\n",
    "    gf = g1+go+g2\n",
    "    #print(gf)\n",
    "    os.system(gf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "folder = '/media/denis/seagate/PB/MODIS/myd09a1/'\n",
    "for i in glob.glob(folder+\"*.tif\"):\n",
    "    a,m = xu.TifToArray(i)\n",
    "    a = np.where(a<-0.2,-2,a)\n",
    "    a = a*10000\n",
    "    xu.ArrayToTif(Array=a,Folder=folder,Filename=i.split(\"/\")[-1],Metadata=m,Type=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "os.chdir(folder)\n",
    "os.system('for f in *.tif; do gdalwarp -ot Int16 -s_srs \"+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs\" -t_srs EPSG:4326 -dstnodata -20000 -r near -of GTiff -te -48.7551 \" -18.3494\" \" -34.729\" -1.0445 -co COMPRESS=DEFLATE -co PREDICTOR=2 -co ZLEVEL=9 -wo OPTIMIZE_SIZE=TRUE $f r$f; done')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculating montly means and maybe other aggregation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "import xuleta as xu\n",
    "import numpy as np\n",
    "import pylab as plt\n",
    "import bottleneck as bn\n",
    "from glob import glob\n",
    "import os\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": [
    "folder = '/media/denis/seagate/PBv2/mxd16a2_anomalies/out/'\n",
    "ul = 12\n",
    "ll = -12\n",
    "nodata = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/usr/local/lib/python3.5/dist-packages/ipykernel_launcher.py:17: RuntimeWarning: invalid value encountered in greater\n",
      "/usr/local/lib/python3.5/dist-packages/ipykernel_launcher.py:18: RuntimeWarning: invalid value encountered in less\n",
      "/usr/local/lib/python3.5/dist-packages/ipykernel_launcher.py:21: RuntimeWarning: Mean of empty slice\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2002-07.tif saved.\n",
      " \n",
      "2002-08.tif saved.\n",
      " \n",
      "2002-09.tif saved.\n",
      " \n",
      "2002-10.tif saved.\n",
      " \n",
      "2002-11.tif saved.\n",
      " \n",
      "2002-12.tif saved.\n",
      " \n",
      "2003-01.tif saved.\n",
      " \n",
      "2003-02.tif saved.\n",
      " \n",
      "2003-03.tif saved.\n",
      " \n",
      "2003-04.tif saved.\n",
      " \n",
      "2003-05.tif saved.\n",
      " \n",
      "2003-06.tif saved.\n",
      " \n",
      "2003-07.tif saved.\n",
      " \n",
      "2003-08.tif saved.\n",
      " \n",
      "2003-09.tif saved.\n",
      " \n",
      "2003-10.tif saved.\n",
      " \n",
      "2003-11.tif saved.\n",
      " \n",
      "2003-12.tif saved.\n",
      " \n",
      "2004-01.tif saved.\n",
      " \n",
      "2004-02.tif saved.\n",
      " \n",
      "2004-03.tif saved.\n",
      " \n",
      "2004-04.tif saved.\n",
      " \n",
      "2004-05.tif saved.\n",
      " \n",
      "2004-06.tif saved.\n",
      " \n",
      "2004-07.tif saved.\n",
      " \n",
      "2004-08.tif saved.\n",
      " \n",
      "2004-09.tif saved.\n",
      " \n",
      "2004-10.tif saved.\n",
      " \n",
      "2004-11.tif saved.\n",
      " \n",
      "2004-12.tif saved.\n",
      " \n",
      "2005-01.tif saved.\n",
      " \n",
      "2005-02.tif saved.\n",
      " \n",
      "2005-03.tif saved.\n",
      " \n",
      "2005-04.tif saved.\n",
      " \n",
      "2005-05.tif saved.\n",
      " \n",
      "2005-06.tif saved.\n",
      " \n",
      "2005-07.tif saved.\n",
      " \n",
      "2005-08.tif saved.\n",
      " \n",
      "2005-09.tif saved.\n",
      " \n",
      "2005-10.tif saved.\n",
      " \n",
      "2005-11.tif saved.\n",
      " \n",
      "2005-12.tif saved.\n",
      " \n",
      "2006-01.tif saved.\n",
      " \n",
      "2006-02.tif saved.\n",
      " \n",
      "2006-03.tif saved.\n",
      " \n",
      "2006-04.tif saved.\n",
      " \n",
      "2006-05.tif saved.\n",
      " \n",
      "2006-06.tif saved.\n",
      " \n",
      "2006-07.tif saved.\n",
      " \n",
      "2006-08.tif saved.\n",
      " \n",
      "2006-09.tif saved.\n",
      " \n",
      "2006-10.tif saved.\n",
      " \n",
      "2006-11.tif saved.\n",
      " \n",
      "2006-12.tif saved.\n",
      " \n",
      "2007-01.tif saved.\n",
      " \n",
      "2007-02.tif saved.\n",
      " \n",
      "2007-03.tif saved.\n",
      " \n",
      "2007-04.tif saved.\n",
      " \n",
      "2007-05.tif saved.\n",
      " \n",
      "2007-06.tif saved.\n",
      " \n",
      "2007-07.tif saved.\n",
      " \n",
      "2007-08.tif saved.\n",
      " \n",
      "2007-09.tif saved.\n",
      " \n",
      "2007-10.tif saved.\n",
      " \n",
      "2007-11.tif saved.\n",
      " \n",
      "2007-12.tif saved.\n",
      " \n",
      "2008-01.tif saved.\n",
      " \n",
      "2008-02.tif saved.\n",
      " \n",
      "2008-03.tif saved.\n",
      " \n",
      "2008-04.tif saved.\n",
      " \n",
      "2008-05.tif saved.\n",
      " \n",
      "2008-06.tif saved.\n",
      " \n",
      "2008-07.tif saved.\n",
      " \n",
      "2008-08.tif saved.\n",
      " \n",
      "2008-09.tif saved.\n",
      " \n",
      "2008-10.tif saved.\n",
      " \n",
      "2008-11.tif saved.\n",
      " \n",
      "2008-12.tif saved.\n",
      " \n",
      "2009-01.tif saved.\n",
      " \n",
      "2009-02.tif saved.\n",
      " \n",
      "2009-03.tif saved.\n",
      " \n",
      "2009-04.tif saved.\n",
      " \n",
      "2009-05.tif saved.\n",
      " \n",
      "2009-06.tif saved.\n",
      " \n",
      "2009-07.tif saved.\n",
      " \n",
      "2009-08.tif saved.\n",
      " \n",
      "2009-09.tif saved.\n",
      " \n",
      "2009-10.tif saved.\n",
      " \n",
      "2009-11.tif saved.\n",
      " \n",
      "2009-12.tif saved.\n",
      " \n",
      "2010-01.tif saved.\n",
      " \n",
      "2010-02.tif saved.\n",
      " \n",
      "2010-03.tif saved.\n",
      " \n",
      "2010-04.tif saved.\n",
      " \n",
      "2010-05.tif saved.\n",
      " \n",
      "2010-06.tif saved.\n",
      " \n",
      "2010-07.tif saved.\n",
      " \n",
      "2010-08.tif saved.\n",
      " \n",
      "2010-09.tif saved.\n",
      " \n",
      "2010-10.tif saved.\n",
      " \n",
      "2010-11.tif saved.\n",
      " \n",
      "2010-12.tif saved.\n",
      " \n",
      "2011-01.tif saved.\n",
      " \n",
      "2011-02.tif saved.\n",
      " \n",
      "2011-03.tif saved.\n",
      " \n",
      "2011-04.tif saved.\n",
      " \n",
      "2011-05.tif saved.\n",
      " \n",
      "2011-06.tif saved.\n",
      " \n",
      "2011-07.tif saved.\n",
      " \n",
      "2011-08.tif saved.\n",
      " \n",
      "2011-09.tif saved.\n",
      " \n",
      "2011-10.tif saved.\n",
      " \n",
      "2011-11.tif saved.\n",
      " \n",
      "2011-12.tif saved.\n",
      " \n",
      "2012-01.tif saved.\n",
      " \n",
      "2012-02.tif saved.\n",
      " \n",
      "2012-03.tif saved.\n",
      " \n",
      "2012-04.tif saved.\n",
      " \n",
      "2012-05.tif saved.\n",
      " \n",
      "2012-06.tif saved.\n",
      " \n",
      "2012-07.tif saved.\n",
      " \n",
      "2012-08.tif saved.\n",
      " \n",
      "2012-09.tif saved.\n",
      " \n",
      "2012-10.tif saved.\n",
      " \n",
      "2012-11.tif saved.\n",
      " \n",
      "2012-12.tif saved.\n",
      " \n",
      "2013-01.tif saved.\n",
      " \n",
      "2013-02.tif saved.\n",
      " \n",
      "2013-03.tif saved.\n",
      " \n",
      "2013-04.tif saved.\n",
      " \n",
      "2013-05.tif saved.\n",
      " \n",
      "2013-06.tif saved.\n",
      " \n",
      "2013-07.tif saved.\n",
      " \n",
      "2013-08.tif saved.\n",
      " \n",
      "2013-09.tif saved.\n",
      " \n",
      "2013-10.tif saved.\n",
      " \n",
      "2013-11.tif saved.\n",
      " \n",
      "2013-12.tif saved.\n",
      " \n",
      "2014-01.tif saved.\n",
      " \n",
      "2014-02.tif saved.\n",
      " \n",
      "2014-03.tif saved.\n",
      " \n",
      "2014-04.tif saved.\n",
      " \n",
      "2014-05.tif saved.\n",
      " \n",
      "2014-06.tif saved.\n",
      " \n",
      "2014-07.tif saved.\n",
      " \n",
      "2014-08.tif saved.\n",
      " \n",
      "2014-09.tif saved.\n",
      " \n",
      "2014-10.tif saved.\n",
      " \n",
      "2014-11.tif saved.\n",
      " \n",
      "2014-12.tif saved.\n",
      " \n",
      "2015-01.tif saved.\n",
      " \n",
      "2015-02.tif saved.\n",
      " \n",
      "2015-03.tif saved.\n",
      " \n",
      "2015-04.tif saved.\n",
      " \n",
      "2015-05.tif saved.\n",
      " \n",
      "2015-06.tif saved.\n",
      " \n",
      "2015-07.tif saved.\n",
      " \n",
      "2015-08.tif saved.\n",
      " \n",
      "2015-09.tif saved.\n",
      " \n",
      "2015-10.tif saved.\n",
      " \n",
      "2015-11.tif saved.\n",
      " \n",
      "2015-12.tif saved.\n",
      " \n",
      "2016-01.tif saved.\n",
      " \n",
      "2016-02.tif saved.\n",
      " \n",
      "2016-03.tif saved.\n",
      " \n",
      "2016-04.tif saved.\n",
      " \n",
      "2016-05.tif saved.\n",
      " \n",
      "2016-06.tif saved.\n",
      " \n",
      "2016-07.tif saved.\n",
      " \n",
      "2016-08.tif saved.\n",
      " \n",
      "2016-09.tif saved.\n",
      " \n",
      "2016-10.tif saved.\n",
      " \n",
      "2016-11.tif saved.\n",
      " \n",
      "2016-12.tif saved.\n",
      " \n",
      "2017-01.tif saved.\n",
      " \n",
      "2017-02.tif saved.\n",
      " \n"
     ]
    }
   ],
   "source": [
    "if not os.path.exists(folder+\"monthly_mean\"):\n",
    "    os.makedirs(folder+\"monthly_mean\")\n",
    "folderout = folder+\"monthly_mean/\"\n",
    "\n",
    "l = sorted(glob(folder+\"*.tif\"))\n",
    "\n",
    "#list of unique Year-months\n",
    "ym = sorted(list(set([ii[-14:-7] for n,ii in enumerate(l) if ii not in l[:n]])))\n",
    "\n",
    "for m in ym:\n",
    "    #list of files that match the current year-month\n",
    "    toOpen = [x for x in l if m in x]\n",
    "    \n",
    "    arrays = []\n",
    "    for j in toOpen:\n",
    "        imagem, meta = xu.TifToArray(j)\n",
    "        imagem = np.where(imagem > ul, np.nan, imagem)\n",
    "        imagem = np.where(imagem < ll, np.nan, imagem)\n",
    "        arrays.append(imagem)\n",
    "        arrays2 = np.dstack(arrays)\n",
    "    xu.ArrayToTif(np.nanmean(arrays2, axis=2), m+'.tif', folderout, meta, Type=3)  \n",
    "    arrays, media, toOpen, arrays2 = None, None, None, None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "#COMPRESS FILES\n",
    "os.chdir(folderout)\n",
    "#-ot Int16\n",
    "#os.system('for f in *.tif; do gdalwarp -co COMPRESS=DEFLATE -co PREDICTOR=1 -co ZLEVEL=9 -multi -of GTiff $f X$f; done')\n",
    "xu.findreplace(folderout,'X','')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true,
    "jupyter": {
     "outputs_hidden": true
    }
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "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.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}