Created
August 14, 2023 04:51
-
-
Save rchardptrsn/68f06fa54d15b6bd7ec87a87af4d2e13 to your computer and use it in GitHub Desktop.
download NDVI
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import timeit | |
import sys | |
import json | |
import ee | |
# Authenticate to google earth engine | |
service_account = '[email protected]' | |
key_path = 'path_to_your_service_account_key.json' | |
credentials = ee.ServiceAccountCredentials(service_account, key_path) | |
ee.Initialize(credentials) | |
# function to retrieve NDVI data from a single multipolygon | |
def getNDVI(fire_polygon_json): | |
featureCollection = ee.FeatureCollection(json.loads(fire_polygon_json)) #json.loads(fire_bounds_geom)) | |
# Any region of the world | |
# polygon = ee.Geometry.Polygon([112.0, 1.0, | |
# 112.0, 1.5, | |
# 112.5, 1.5, | |
# 112.5, 1.0, | |
# 112.0, 1.0]) | |
startDate = '2021-01-01' | |
endDate = '2023-04-01' | |
# An ImageCollection is a stack or sequence of images | |
modisNDVI = ee.ImageCollection('MODIS/MOD09GA_006_NDVI').select('NDVI').filterDate(startDate, endDate) | |
def getEarthEngineData(n): | |
date = ee.Date(startDate).advance(n,'month') | |
m = date.get("month") | |
y = date.get("year") | |
dic = ee.Dictionary({ # Create the earth engine dictionary | |
'Date':date.format('yyyy-MM') | |
}) | |
tempNDVI = (modisNDVI.filter(ee.Filter.calendarRange(y, y, 'year')) | |
.filter(ee.Filter.calendarRange(m, m, 'month')) | |
.mean() | |
.reduceRegion( | |
reducer = ee.Reducer.mean(), | |
geometry = featureCollection.geometry(), | |
scale = 250)) | |
return dic.combine(tempNDVI) | |
modis_YrMo = ee.List.sequence(0, 12*2-1).map(getEarthEngineData) | |
dataframe = gpd.GeoDataFrame(modis_YrMo.getInfo()) | |
return dataframe | |
################################################# | |
# Function to pass a list of polygons to getNDVI() | |
# and save the results to a file | |
def collectAllNDVI(gdf): | |
for i in range(len(gdf)): | |
# counters | |
record = i | |
next = record + 1 | |
# extract GeoDataFrame geometry to json | |
fire_polygon_json = gdf.iloc[record:next].geometry.to_json() | |
# extract the objectid | |
objectid = gdf.iloc[record:next,0].values[0] | |
# pass the json record for the fire to our getNDVI function | |
# to query Google Earth Engine for NDVI | |
try: | |
start_time = timeit.default_timer() | |
print(f'Getting NDVI data for objectid: {objectid}') | |
# call NDVI function to retrieve NDVI as a pandas dataframe for a given polygon | |
firedf = getNDVI(fire_polygon_json) | |
elapsed = timeit.default_timer() - start_time | |
print(f'Successfully retrieved NDVI data for objectid: {objectid} in {elapsed} seconds') | |
# Add the OBJECTID to each record | |
firedf['OBJECTID'] = objectid | |
# save to file | |
firedf.to_csv(f'ndvi data/{objectid}.csv', index=False) | |
print(f'Saved NDVI data for objectid {objectid} at ndvi data/{objectid}.csv') | |
# write confirmation to log file | |
with open("output.out", 'a') as f: | |
f.write(f'Saved NDVI data for objectid {objectid} at ndvi data/{objectid}.csv in {elapsed} seconds.\n') | |
f.close() | |
except Exception as error: | |
print(f'Error: could not retrieve NDVI data for objectid: {objectid}') | |
print(error) | |
with open("output.out", 'a') as f: | |
f.write(f'Error: could not save NDVI data for objectid {objectid}\n') | |
f.close() | |
################################################# | |
# Execute functions | |
# Define batch download based on row number of GeoPandas dataframe | |
# I found that batches of 250-500 would complete without issues | |
# Some images take longer to download than others based on the size of the polygon | |
ndvi_gdf_run = gdf.iloc[0:250] | |
# Collect all | |
collectAllNDVI(ndvi_gdf_run) | |
print('Complete.') |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment