Merge pull request #70 from developmentseed/develop

version 0.6.1 release

.travis.yml

@@ -1,9 +1,5 @@
 language: python
 
-branches:
- only:
- - /^v0\.5.*$/
-
 python:
  - '2.7'
 
@@ -23,12 +19,10 @@ before_install:
 
 # Install packages
 install:
- - conda install --yes python=$TRAVIS_PYTHON_VERSION numpy scipy matplotlib gdal scikit-image requests six nose dateutil
+ - conda install --yes python=$TRAVIS_PYTHON_VERSION numpy scipy matplotlib scikit-image requests six nose dateutil
+ - conda install --yes -c https://conda.binstar.org/osgeo gdal
  - pip install --install-option="--no-cython-compile" cython
  - pip install -r requirements/travis.txt
 
 script:
  - nosetests
-
-after_success:
- - coveralls

README.rst

@@ -17,41 +17,38 @@ Installation
 
 **On Mac**
 
-.. code-block:: console
-
- $: pip install landsat-util
+ ``$: pip install landsat-util``
 
 **On Ubuntu 14.04**
 
 Use pip to install landsat-util. If you are not using virtualenv, you might have to run ``pip`` as ``sudo``.
 
-.. code-block:: console
-
- $: sudo apt-get update
- $: sudo apt-get install python-pip python-numpy python-scipy libgdal-dev libatlas-base-dev gfortran
- $: pip install landsat-util
+ ``$: sudo apt-get update``
+ ``$: sudo apt-get install python-pip python-numpy python-scipy libgdal-dev libatlas-base-dev gfortran``
+ ``$: pip install landsat-util``
 
 **On Other systems**
 
-.. code-block:: console
-
- $: python setup.py install
+ ``$: python setup.py numpy six``
+ ``$: python setup.py install``
 
 
 **To Upgrade**
 
-.. code-block:: console
-
- $: pip install -U landsat-util
+ ``$: pip install -U landsat-util``
 
 If you have installed previous version of landsat using brew, first run:
 
-.. code-block:: console
+ ``$: brew uninstall landsat-util``
+
+**To Test**
 
- $: brew uninstall landsat-util
+ ``$: pip install -U requirements/dev.txt``
+ ``$: nosetests``
 
 Overview: What can landsat-util do?
-============
+====================================
+
 Landsat-util has three main functions:
 
 - **Search** for landsat tiles based on several search parameters.
@@ -63,7 +60,7 @@ These three functions have to be performed separately.
 **Help**: Type ``landsat -h`` for detailed usage parameters.
 
 Step 1: Search
-============
+===============
 
 Search returns information about all landsat tiles that match your criteria. This includes a link to an unprocessed preview of the tile. The most important result is the tile's *sceneID*, which you will need to download the tile (see step 2 below).
 
@@ -89,7 +86,7 @@ Search by latitude and longitude:
 
 
 Step 2: Download
-============
+=================
 
 You can download tiles using their unique sceneID, which you get from landsat search.
 
@@ -110,7 +107,7 @@ Download multiple sceneIDs:
 ``$: landsat download LC80090452014008LGN00 LC80090452015008LGN00 LC80090452013008LGN00``
 
 Step 3: Image processing
-============
+=========================
 
 You can process your downloaded tiles with our custom image processing algorithms. In addition, you can choose to pansharpen your images and specify which bands to process.
 
@@ -145,7 +142,7 @@ Important Notes
 - Landsat-util requires at least 2GB of Memory (RAM).
 
 Recently Added
-++++++++++
++++++++++++++++
 
 - Add longitude latitude search
 - Improve console output

Vagrantfile

@@ -25,7 +25,8 @@ Vagrant.configure(2) do |config|
  sudo rm -rf /dev/shm
  sudo ln -s /run/shm /dev/shm
 
- conda install --yes numpy scikit-image requests gdal nose dateutil
+ conda install --yes -c https://conda.binstar.org/osgeo gdal
+ conda install --yes numpy scikit-image requests nose dateutil
 
  pip install --install-option="--no-cython-compile" cython
  pip install -r /vagrant/requirements/travis.txt

landsat/__init__.py

@@ -1,7 +1 @@
-import settings
-import sys
-
-if not settings.DEBUG:
- sys.tracebacklimit = 0
-
-__version__ = '0.5.1'
+__version__ = '0.6.1'

landsat/downloader.py

@@ -39,12 +39,19 @@ def download(self, scenes, bands=None):
 
  if isinstance(scenes, list):
  for scene in scenes:
- if bands:
+ # If bands are provided the image is from 2015 or later use Amazon
+ if (bands and int(scene[12]) > 4):
  if isinstance(bands, list):
  # Create a folder to download the specific bands into
  path = check_create_folder(join(self.download_dir, scene))
- for band in bands:
- self.amazon_s3(scene, band, path)
+ try:
+ # Always grab MTL.txt if bands are specified
+ bands_plus = bands
+ bands_plus.append('MTL')
+ for band in bands_plus:
+ self.amazon_s3(scene, band, path)
+ except RemoteFileDoesntExist:
+ self.google_storage(scene, self.download_dir)
  else:
  raise Exception('Expected bands list')
  else:
@@ -63,23 +70,26 @@ def google_storage(self, scene, path):
  url = self.google_storage_url(sat)
 
  if self.remote_file_exists(url):
- self.fetch(url, path, filename)
+ return self.fetch(url, path, filename)
 
  else:
- RemoteFileDoesntExist('%s is not available on Google Storage' % filename)
+ raise RemoteFileDoesntExist('%s is not available on Google Storage' % filename)
 
  def amazon_s3(self, scene, band, path):
  """ Amazon S3 downloader """
  sat = self.scene_interpreter(scene)
 
- filename = '%s_B%s.TIF' % (scene, band)
+ if band != 'MTL':
+ filename = '%s_B%s.TIF' % (scene, band)
+ else:
+ filename = '%s_%s.txt' % (scene, band)
  url = self.amazon_s3_url(sat, filename)
 
  if self.remote_file_exists(url):
- self.fetch(url, path, filename)
+ return self.fetch(url, path, filename)
 
  else:
- RemoteFileDoesntExist('%s is not available on Amazon S3' % filename)
+ raise RemoteFileDoesntExist('%s is not available on Amazon S3' % filename)
 
  def fetch(self, url, path, filename):
 

landsat/image.py

@@ -6,13 +6,16 @@
 import sys
 from os.path import join
 import tarfile
+import glob
+import subprocess
+
 import numpy
 import rasterio
-import glob
 from rasterio.warp import reproject, RESAMPLING, transform
 
-from skimage import img_as_ubyte, exposure
 from skimage import transform as sktransform
+from skimage.util import img_as_ubyte
+from skimage.exposure import rescale_intensity, adjust_gamma
 
 import settings
 from mixins import VerbosityMixin
@@ -66,6 +69,18 @@ def run(self, pansharpen=True):
 
  self.output("* Image processing started for bands %s" % "-".join(map(str, self.bands)), normal=True)
 
+ # Read cloud coverage from mtl file
+ cloud_cover = 0
+ try:
+ with open(self.scene_path + '/' + self.scene + '_MTL.txt', 'rU') as mtl:
+ lines = mtl.readlines()
+ for line in lines:
+ if 'CLOUD_COVER' in line:
+ cloud_cover = float(line.replace('CLOUD_COVER = ', ''))
+ break
+ except IOError:
+ pass
+
  with warnings.catch_warnings():
  warnings.simplefilter("ignore")
  with rasterio.drivers():
@@ -103,10 +118,17 @@ def run(self, pansharpen=True):
  crn = self._get_boundaries(src_data)
 
  dst_shape = src_data['shape']
- y_pixel = (max(crn['ul']['y'][1][0],crn['ur']['y'][1][0])- min(crn['lr']['y'][1][0],crn['ll']['y'][1][0]))/dst_shape[0]
- x_pixel = (max(crn['lr']['x'][1][0],crn['ur']['x'][1][0]) - min(crn['ul']['x'][1][0],crn['ll']['x'][1][0]))/dst_shape[1]
-
- dst_transform = (min(crn['ul']['x'][1][0],crn['ll']['x'][1][0]), x_pixel, 0.0, max(crn['ul']['y'][1][0],crn['ur']['y'][1][0]), 0.0, -y_pixel)
+ y_pixel = (max(crn['ul']['y'][1][0], crn['ur']['y'][1][0]) -
+ min(crn['lr']['y'][1][0], crn['ll']['y'][1][0])) / dst_shape[0]
+ x_pixel = (max(crn['lr']['x'][1][0], crn['ur']['x'][1][0]) -
+ min(crn['ul']['x'][1][0], crn['ll']['x'][1][0])) / dst_shape[1]
+
+ dst_transform = (min(crn['ul']['x'][1][0], crn['ll']['x'][1][0]),
+ x_pixel,
+ 0.0,
+ max(crn['ul']['y'][1][0], crn['ur']['y'][1][0]),
+ 0.0,
+ -y_pixel)
  # Delete crn since no longer needed
  del crn
 
@@ -149,7 +171,7 @@ def run(self, pansharpen=True):
 
  for i, band in enumerate(new_bands):
  # Color Correction
- band = self._color_correction(band, self.bands[i])
+ band = self._color_correction(band, self.bands[i], 0.1, cloud_cover)
 
  # Gamma Correction
  if i == 0:
@@ -158,7 +180,6 @@ def run(self, pansharpen=True):
  if i == 2:
  band = self._gamma_correction(band, 0.9)
 
-
  output.write_band(i+1, img_as_ubyte(band))
 
  new_bands[i] = None
@@ -187,16 +208,19 @@ def _pansharpenning(self, bands):
  return bands
 
  def _gamma_correction(self, band, value):
- return exposure.adjust_gamma(band, value)
+ return adjust_gamma(band, value)
 
- def _color_correction(self, band, band_id):
+ def _color_correction(self, band, band_id, low, cloud_cover):
  band = band.astype(numpy.uint16)
 
  self.output("Color correcting band %s" % band_id, normal=True, color='green', indent=1)
- p2, p98 = self._percent_cut(band)
- band = exposure.rescale_intensity(band, in_range=(p2, p98))
-
- return band
+ p_low, cloud_cut_low = self._percent_cut(band, low, 100 - (cloud_cover * 3 / 4))
+ temp = numpy.zeros(numpy.shape(band), dtype=numpy.uint16)
+ cloud_divide = 65000 - cloud_cover * 100
+ mask = numpy.logical_and(band < cloud_cut_low, band > 0)
+ temp[mask] = rescale_intensity(band[mask], in_range=(p_low, cloud_cut_low), out_range=(256, cloud_divide))
+ temp[band >= cloud_cut_low] = rescale_intensity(band[band >= cloud_cut_low], out_range=(cloud_divide, 65535))
+ return temp
 
  def _read_band(self, band_path):
  """ Reads a band with rasterio """
@@ -248,15 +272,20 @@ def _get_boundaries(self, src):
 
  return output
 
- def _percent_cut(self, color):
- return numpy.percentile(color[numpy.logical_and(color > 0, color < 65535)], (2, 98))
+ def _percent_cut(self, color, low, high):
+ return numpy.percentile(color[numpy.logical_and(color > 0, color < 65535)], (low, high))
 
  def _unzip(self, src, dst, scene):
  """ Unzip tar files """
  self.output("Unzipping %s - It might take some time" % scene, normal=True, arrow=True)
- tar = tarfile.open(src)
- tar.extractall(path=dst)
- tar.close()
+
+ try:
+ tar = tarfile.open(src, 'r')
+ tar.extractall(path=dst)
+ tar.close()
+ except tarfile.ReadError:
+ check_create_folder(dst)
+ subprocess.check_call(['tar', '-xf', src, '-C', dst])
 
  def _get_full_filename(self, band):