I would like to read and process Landsat images from Amazon S3 using Python
At this point I haven't implemented this stuff myself, but the approach that I am planning on is combining boto with the information on Amazon's page on Landsat on AWS. At first glance, it appears to be a rather doable task. Another option is landsat-util which can be found on github, which is python-based and can interface with AWS.
You should be able to read Landsat directly out of S3 using
/vsicurl/ file handler. Before coding it up in Python, check if you can access a file from the commandline. I think
/vsicurl/ may work in gdal-2.0 or 1.11, but it's best to use the latest (2.1.2).
There are a couple of environment variables you can set to improve performance. Try this:
env CPL_VSIL_CURL_ALLOWED_EXTENSIONS=tif VSI_CACHE=TRUE GDAL_DISABLE_READDIR_ON_OPEN=TRUE gdalinfo /vsicurl/http://landsat-pds.s3.amazonaws.com/L8/001/003/LC80010032014272LGN00/LC80010032014272LGN00_B1.TIF
Note some of these options only work with gdal >= 2.1.2 which was released in November, 2016. If you can get gdalinfo to display the correct metadata, Python access to the dataset should work. Here's a random example:
from osgeo import gdal import numpy gdal.SetConfigOption('CPL_VSIL_CURL_ALLOWED_EXTENSIONS', 'tif') gdal.SetConfigOption('VSI_CACHE', 'TRUE') gdal.SetConfigOption('GDAL_DISABLE_READDIR_ON_OPEN', 'TRUE') ds = gdal.Open('/vsicurl/http://landsat-pds.s3.amazonaws.com/L8/001/003/LC80010032014272LGN00/LC80010032014272LGN00_B1.TIF') xoff, yoff, xcount, ycount = (4000, 4000, 256, 256) np_array = ds.GetRasterBand(1).ReadAsArray(xoff, yoff, xcount, ycount)