Downloading images from Google Earth Engine to local drive ( A case study for Timesat image Preparation for Planting Date)

Downloading images from Google Earth Engine to local drive ( A case study for Timesat image Preparation for Planting Date)

Satellite image processing might be very painful, depending upon the volume of data. I don't like to be suggesting the data download and do image analysis on your local device. But there might be a time when we can't follow this suggestion. I passed through a similar situation for analyzing MODIS EVI datasets to determine the planting date for Rice.

I was advised to use Timesat (Welcome to the TIMESAT pages! (lu.se)) for this task. Instead, I had no options to download MODIS EVI data (250 m resolution, 8 Day temporal resolution-terra and Aqua sensors combined). I used Google Earth Engine (GEE) to download the data into my local device and prepare the dataset in the desired format for Timesat (A software package to analyze time-series of satellite sensor data). I would never suggest downloading the Satellite data unless absolutely necessary.

This article will be the perfect solution for those looking for a way to download satellite images from Google Earth Engine (GEE). But I have not yet tested this for large images; that might not be possible due to GEE restriction and is not recommended. 

Let's dive-in directly to the code segment with Google Earth Engine (GEE) python Api.

  1. Import the GEE API and Authenticate

    import ee, geemap
    # import geopandas as gpd
    # Authenticate the earthengine with credentials
    ee.Initialize()
    
    # Import request package later used for Downloading data
    import requests
  2. Import data for declaration of area of interest (AOI)
    # get our Nepal boundary
    # I have taken level 0 data for country data from FAO datasets
    
    # aoi = ee.FeatureCollection("FAO/GAUL/2015/level0") \
    #    .filter(ee.Filter.eq('ADM0_NAME','Nepal')).geometry() # adjust indentation here, May get error
    
    # i am taking shapefile for plain areas of Nepal
    study_area_shp = 'shapefile/terai.shp'
    aoi = geemap.shp_to_ee(study_area_shp).geometry()
    
  3. Define starting and ending year of study (Both year Inclusive)
    # Longer the duration, Longer it will take for Processing and downloading of Images
    StartDate = '2015-01-01'
    EndDate = '2021-12-31'
  4. Filter Image for desired time series
    # add satellite time series: MODIS EVI 250m 16 day -------------
    # terra sensor
    collectionModEvi_terra = ee.ImageCollection('MODIS/006/MOD13Q1').filterDate(StartDate,EndDate) \
        .filterBounds(aoi)\
        .select('EVI')
    
    # Aqua sensor
    collectionModEvi_aqua = ee.ImageCollection('MODIS/061/MYD13Q1').filterDate(StartDate,EndDate) \
        .filterBounds(aoi)\
        .select('EVI');
    
    collectionModEvi = collectionModEvi_terra.merge(collectionModEvi_aqua)
    # THis will provide us 250 m of EVI datasets from MODIS on 8 day interval

    We are combining data from 2 sensors of MODIS (Terra ANd Aqua) so that we will be getting 250 m resolution of MODIS EVI datasets in 8 days interval time.

  5. Loop Over Image in Image Collection and Download individual Image with desired name

    # https://stackoverflow.com/questions/46943061/how-to-iterate-over-and-download-each-image-in-an-image-collection-from-the-goog
    # This is OK for small collections
    collectionList = collectionModEvi.toList(collectionModEvi.size())
    # print(collectionList)
    collectionSize = collectionList.size().getInfo()
    print(collectionSize)
    for i in range(collectionSize):
        image = ee.Image(collectionList.get(i))
        image_name = image.get('system:index').getInfo()
        print(i)
        print(image_name)
        img_name = "Modis_evi_" + str(image_name)
        url = image.getDownloadUrl({
                        'name':img_name,
                        'scale': 250,
                        'crs': 'EPSG:4326',
                        'region': aoi,
                        'format':"GEO_TIFF"
                    })
        print(url)
        r = requests.get(url, stream=True)
        filenameTif = 'modis_data/' + img_name + '.tif'
        with open(filenameTif, "wb") as fd:
                for chunk in r.iter_content(chunk_size=1024):
                    fd.write(chunk)
        

    Up to this step, We have dealt with data download from Google Earth Engine (GEE) directly in our local device. Now we will move on to Data preparation for Timesat software. I have striclty followed the documentation of Timesat for data preparation and this is properly working in my setup.

  6. Now will move on to Image Clipping and converison to Binary Files

    # Importing Rasterio and fiona module for masking rasters 
    import fiona
    import rasterio
    import rasterio.mask
  7. Defining the masking Geometry 
    with fiona.open("shapefile/terai.shp", "r") as shapefile:
        shapes = [feature["geometry"] for feature in shapefile]
  8. Clipping all rasters image with shapefile
    import os
    
    path = "modis_data"
    file_type = '.tif'
    
    for filename in os.listdir(path=path):
        if filename.endswith(file_type):
    #         print(filename)
            print(f"{path}/{filename}")
            filepath = f"{path}/{filename}"
            
            masked_name = filename.replace('.tif', '_masked.tif')
            masked_file_path =  "masked/" + f"{masked_name}"
            print(massked_file_path)
            with rasterio.open(filepath) as src:
                out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
                out_meta = src.meta
    
            out_meta.update({"driver": "GTiff",
                             "height": out_image.shape[1],
                             "width": out_image.shape[2],
                             "transform": out_transform})
    
            with rasterio.open(masked_file_path, "w", **out_meta) as dest:
                dest.write(out_image)
  9. Coversion of Tif image to binary image (Envi format, BIL image)
    path = "masked"
    file_type = '.tif'
    
    for filename in os.listdir(path=path):
        if filename.endswith(file_type):
            input = f"{path}/{filename}"
            new_filename = filename.replace('.tif', '.bil')
            output = f"binary/{new_filename}"
            print(input)
            print(output)
            os.system("gdal_translate -of ENVI " + input + " " +  output)
  10. Making list of BInary (BIL) images and writing into Text files 
    path = "masked"
    file_type = '.tif'
    
    absolute_path = r"C:\projects\planting_date"
    
    absolute_file_path_list = []
    for filename in os.listdir(path=path):
        if filename.endswith(file_type):
            new_filename = filename.replace('.tif', '.bil')
            output = absolute_path + r"\binary\{}".format(new_filename)
            absolute_file_path_list.append(output)
    

    This is required for Timesat image analysis software.

  11. Writing all files path into text files 

    with open('filelist.txt', 'w') as f:
        f.write(str(len(absolute_file_path_list)))
        f.write('\n')
        for line in absolute_file_path_list:
            f.write(line)
            f.write('\n')

Source code: github_link

This code need to be optimized and will be doing it if i get any better ideas. If you have any better idea, please do let me know the better way of optimizing the downloading techinque from Google Earth Engine (GEE) into local Device. 

Please let us know your thoughts on this with following comment section we will definetly reach out to you on your query and feedback.

1 Comments

Shangharsha Thapa | General comments

May 30, 2022, 6:47 p.m.

A very nice and detailed explanation of data preparation for TIMESAT-based vegetation phenology and seasonality parameter extraction. If the study area is the EU region, you can already use TIMESAT derived Plant Phenology Index (PPI) with quality flags as well as the 13 different layers related to the vegetation seasonality (See this: https://www.wekeo.eu/data?view=viewer).

Leave a Reply