extract_timeseries.py

import ee
import logging
import multiprocessing
from retry import retry


"""
This tool downloads data from Earth Engine using parallel requests.

It extracts the timeseries of 8-day max LST from MODIS MOD11A2 per GAUL level-2 region
for all regions in South America, with each time-series written to its own file.
"""

ee.Initialize(opt_url='https://earthengine-highvolume.googleapis.com')


# All the results will be stored here so we can output a single file at the end.

def getRequests():
  """Generates a list of work items to be downloaded.

  Extract the ADM2_CODEs from the GAUL level 2 dataset as work units.
  """
  southAmerica = ee.Geometry.BBox(-84.0, -56.5, -32.0, 12.5)
  gaul2 = (ee.FeatureCollection('FAO/GAUL_SIMPLIFIED_500m/2015/level2')
           .filterBounds(southAmerica))
  return gaul2.aggregate_array('ADM2_CODE').getInfo()


@retry(tries=10, delay=1, backoff=2)
def getResult(index, regionID):
  """Handle the HTTP requests to download one result."""
  region = (ee.FeatureCollection('FAO/GAUL_SIMPLIFIED_500m/2015/level2')
            .filter(ee.Filter.eq('ADM2_CODE', regionID))
            .first())

  def maxLST(image):
    # Mappable function to aggregate the max LST for one image.
    # It builds an output tuple of (max_LST, date, regionID)
    # This function uses -999 to indicate no data.
    date = image.date().format('YYYY-MM-dd')
    image = image.multiply(0.02).subtract(273.15)
    maxValue = (image.reduceRegion(ee.Reducer.max(), region.geometry())
                # set a default in case there's no data in the region
                .combine({'LST_Day_1km': -999}, False)
                .getNumber('LST_Day_1km')
                # format to 2 decimal places.
                .format('%.2f'))

    return image.set('output', [maxValue, date, regionID])

  # Get the max LST for this region, in each image.
  timeSeries = (ee.ImageCollection('MODIS/006/MOD11A2')
                .select('LST_Day_1km')
                .map(maxLST))
  result = timeSeries.aggregate_array('output').getInfo()

  # Write the results to a file.
  filename = 'results_%d.csv' % regionID
  with open(filename, 'w') as out_file:
    for items in result:
      line = ','.join([str(item) for item in items])
      print(line, file=out_file)

  print("Done: ", index)


if __name__ == '__main__':
  logging.basicConfig()
  items = getRequests()

  pool = multiprocessing.Pool(25)
  pool.starmap(getResult, enumerate(items))

  pool.close()
  pool.join()