Simple Landsat 8 TOA, working with MTLs as of mid-2015. NOT CAREFULLY TESTED! Scale factor of 55k.

l8_toa.py

import click
import rasterio as rio
import os
import numpy as np
import re
from sys import argv, exit
from math import sin, radians

def get_mtl_path(bundle_path):
  mtls = [
    file for file in
    list(os.walk(bundle_path))[0][2]
    if file.endswith('_MTL.txt')
  ]

  if len(mtls) == 0:
    exit('Giving up: could not find an *_MTL.txt in %s' % (bundle_path))
  if len(mtls) > 1:
    exit('Giving up: bundle has more than one *_MTL.txt')

  return os.path.join(bundle_path, mtls[0])

def parse_mtl(mtl):
  # assumes no duplicated key strings between groups
  metadata = {}
  line_number = 0
  for line in mtl:
    line_number += 1
    if (
      'GROUP =' in line # repeated
      or 'END_GROUP =' in line # repeated
      or 'END' == line.strip() # contains no =
      or '' == line # obviously
      or 'ORIGIN =' in line # contains spaces on RHS
    ):
      continue
    try:
      key, equals, value = line.split()
      assert(equals == '=')
      try:
        value = float(value)
      except:
        if value.startswith('"') and value.endswith('"'):
          value = value[1:-1]
      metadata[key] = value
    except:
      exit('Giving up: could not parse MTL (line number %s)' % (line_number))
  return metadata

@click.command()
@click.argument('bundle_path', click.Path(exists=True))
@click.argument('band_index', type=click.INT)
@click.argument('dst_path', click.Path(exists=False)) # existence bug
def toa(bundle_path, band_index, dst_path):
  mtl_path = get_mtl_path(bundle_path)
  print('MTL: %s' % (mtl_path))

  with open(mtl_path) as mtl:
    metadata = parse_mtl(mtl)
  
  src_path = os.path.join(
    bundle_path,
    metadata['FILE_NAME_BAND_%s' % (band_index)]
  )
  
  # http://landsat.usgs.gov/Landsat8_Using_Product.php
  
  Mr = metadata['REFLECTANCE_MULT_BAND_%s' % (band_index)]
  Ar = metadata['REFLECTANCE_ADD_BAND_%s' % (band_index)]
  SE = radians(metadata['SUN_ELEVATION'])
  
  with rio.open(src_path) as src:
    tiff_meta = src.meta
    del tiff_meta['transform']
    #tiff_meta.update({'compress': 'LZW'})
    dn = (src.read()[0]).astype(np.float32)

  toa = ((dn * Mr) + Ar) / sin(SE)
  print np.amin(toa), np.amax(toa)
  scale_factor = 55000 # 16000, 100000, ...
  toa = (toa * scale_factor)
  toa = np.clip(toa, 1, np.iinfo(np.uint16).max).astype(np.uint16)
  print np.amin(toa), np.amax(toa)
  toa[np.where(dn == 0)] = 0

  with rio.open(dst_path, 'w', **tiff_meta) as dst:  
    dst.write(toa, 1)
  

if __name__ == '__main__':
  toa()