better-satcen-00004 data pipeline results (Sentinel-2 DVI, composites and BOA reflectances): loading from pickle file and exploiting¶
This Notebook shows how to: * reload better-satcen-00004 data pipeline results previously stored in a pickle file * Create a new dataframe with data values cropped over the area of interest * Plot data * Create a geotiff * Download data
Import needed modules¶
In [1]:
import pandas as pd
from geopandas import GeoDataFrame
import gdal
import numpy as np
from shapely.geometry import Point
from shapely.geometry import Polygon
import matplotlib
import matplotlib.pyplot as plt
from PIL import Image
%matplotlib inline
from shapely.wkt import loads
from shapely.geometry import box
from urlparse import urlparse
import requests
Read picke data file¶
In [2]:
pickle_filename = 'better-satcen-00004_Jun-Jul_2018_Afghanistan.pkl'
In [3]:
results = GeoDataFrame(pd.read_pickle(pickle_filename))
In [4]:
results.head()
Out[4]:
enclosure | enddate | identifier | self | startdate | title | wkt | |
---|---|---|---|---|---|---|---|
0 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 01be9eec66cb994310e469e012a3a2bc670f92e8 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_2... | POLYGON((64.0648798642054 31.5308017399005,65.... |
1 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 037133bb420b2f5407a223e2f2709e7c6eefbdf0 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_2... | POLYGON((63.9421703624035 30.4373158799178,64.... |
2 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 03aa38a9b430af9f2b2b2cb47cae31788c0fd59c | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... |
3 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 0949c278e99a1b2244ddd0edb3be4e38c7b10e81 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_2... | POLYGON((63.9422209286482 31.5415371135011,64.... |
4 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 0cccedf5d105c98b1bce071eebd3639084ebc43d | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... |
Credentials for ellip platform access¶
- Provide here your ellip username and api key for platform access
In [5]:
import getpass
user = getpass.getpass("Enter your username:")
api_key = getpass.getpass("Enter the API key:")
Define filter parameters¶
Time of interest¶
In [6]:
mydates = ['2018-07-29 06:16:31.024', '2018-07-27 06:26:29.024','2018-07-24 06:16:29.024', '2018-07-22 06:26:31.024','2018-07-19 06:16:31.024']
Area of interest¶
- The user can choose to define an AOI selecting a Point and a buffer size used to build a squared polygon around that point
In [7]:
point_of_interest = Point(64.108, 31.586)
In [8]:
buffer_size = 0.07
In [9]:
aoi_wkt1 = box(*point_of_interest.buffer(buffer_size).bounds)
In [10]:
aoi_wkt1.wkt
Out[10]:
'POLYGON ((64.178 31.516, 64.178 31.656, 64.03800000000001 31.656, 64.03800000000001 31.516, 64.178 31.516))'
- Or define a Polygon (in this case this is exactly the Afghanistan reference AOI) from a WKT geometry
In [11]:
geom = 'MULTIPOLYGON (((-77.24516666666666 0.9424444444444445, -77.24516666666666 2.753833333333333, -79.0253888888889 2.753833333333333, -79.0253888888889 0.9424444444444445, -77.24516666666666 0.9424444444444445)), ((21.29611111111111 39.58638888888889, 21.29611111111111 41.032, 19.89788888888889 41.032, 19.89788888888889 39.58638888888889, 21.29611111111111 39.58638888888889)), ((65.02055555555556 30.43894444444445, 65.02055555555556 33.39566666666666, 63.94222222222222 33.39566666666666, 63.94222222222222 30.43894444444445, 65.02055555555556 30.43894444444445)))'
# Afghanistan area
aoi_wkt2 = loads(geom)[2]
In [12]:
aoi_wkt2.wkt
Out[12]:
'POLYGON ((65.02055555555556 30.43894444444445, 65.02055555555556 33.39566666666666, 63.94222222222222 33.39566666666666, 63.94222222222222 30.43894444444445, 65.02055555555556 30.43894444444445))'
Create a new dataframe with data values¶
- Get the subframe data values related to selected tiles and bands, and create a new GeoDataFrame having the original metadata stack with the specific bands data extended info (data values, size, and geo projection and transformation).
Define auxiliary methods to create new dataframe¶
- The first method gets data url reference, bands list and names, cropping flag, AOI cropping area as inputs and returns the related n bands data array with related geospatial infos
- If crop = False the original bands extent is returned
In [13]:
def get_bands_as_array(url,band_indices,bands_name,crop,bbox):
output = '/vsimem/clip.tif'
bands=[]
ds = gdal.Open('/vsicurl/%s' % url)
bands_count = ds.RasterCount
if bands_count < len(band_indices):
raise('ERROR: Selected bands are more than the data ones')
if crop == True:
ulx, uly, lrx, lry = bbox[0], bbox[3], bbox[2], bbox[1]
ds = gdal.Translate(output, ds, bandList=band_indices, projWin = [ulx, uly, lrx, lry], projWinSRS = 'EPSG:4326')
else:
ds = gdal.Translate(output, ds, bandList=band_indices)
ds = None
ds = gdal.Open(output)
w = ds.GetRasterBand(1).XSize
h = ds.GetRasterBand(1).YSize
geo_transform = ds.GetGeoTransform()
projection = ds.GetProjection()
for i in range(len(band_indices)):
data = ds.GetRasterBand(i+1).ReadAsArray(0, 0, w, h).astype(np.float)
data[data==-9999] = np.nan
bands.append(data)
band = None
data = None
ds = None
return bands,geo_transform,projection,w,h
- The second method gets reference url (for composites data), cropping flag and AOI cropping area as inputs and returns the related n bands data array with geospatial infos
- If crop = False the original bands extent is returned
In [14]:
def get_composites_bands_as_array(url,crop,bbox):
output = '/vsimem/clip.tif'
ds = gdal.Open('/vsicurl/%s' % url)
if crop == True:
ulx, uly, lrx, lry = bbox[0], bbox[3], bbox[2], bbox[1]
ds = gdal.Translate(output, ds, projWin = [ulx, uly, lrx, lry], projWinSRS = 'EPSG:4326')
ds = None
ds = gdal.Open(output)
geo_transform = ds.GetGeoTransform()
projection = ds.GetProjection()
data = ds.ReadAsArray().astype(np.float)
h = data.shape[1]
w = data.shape[2]
data[data==-9999] = np.nan
ds = None
band = None
return data,geo_transform,projection,w,h
- The third one selects the data to be managed (only the actual result related to the selected tile) and returns a new GeoDataFrame containing the requested bands extended info with the original metadata set
In [15]:
def get_GDF_with_datavalues(row, user, api_key, bands_list='all', crop=False, bbox=None):
bands_name = ['B1', 'B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'B8', 'B8A', 'B9', 'B11', 'B12','dvi', 'quality_cloud_confidence', 'quality_snow_confidence', 'quality_scene_classification']
band_indices = range(1,len(bands_name)+1)
if bands_list is not 'all':
band_indices = [bands_name.index(x)+1 for x in bands_list ]
# get data values: everything but the mask tif result
parsed_url = urlparse(row['enclosure'])
url = '%s://%s:%s@%s/api%s' % (list(parsed_url)[0], user, api_key, list(parsed_url)[1], list(parsed_url)[2])
extended_info = list()
if 'DVI_BOA.tif' in row['title']:
print 'Getting bands %s for %s' %(band_indices, row['title'])
bands,geo_transform,projection,w,h = get_bands_as_array(url,band_indices,bands_name,crop,bbox)
for i,bi in enumerate(band_indices):
extended_info.append(dict(row))
extended_info[i]['band'] = bands_name[bi-1]
extended_info[i]['data'] = bands[i]
extended_info[i]['geo_transform'] = geo_transform
extended_info[i]['projection'] = projection
extended_info[i]['xsize'] = w
extended_info[i]['ysize'] = h
else:
print 'Getting composite bands %s' %row['title']
data,geo_transform,projection,w,h = get_composites_bands_as_array(url,crop,bbox)
extended_info.append(dict(row))
extended_info[0]['band'] = '_'.join(row['title'].split('.')[0].split('_')[7:])
extended_info[0]['data'] = data
extended_info[0]['geo_transform'] = geo_transform
extended_info[0]['projection'] = projection
extended_info[0]['xsize'] = w
extended_info[0]['ysize'] = h
#print 'Data values for %s : GOT' %row['title']
return extended_info
- Define the reference tiles
In [16]:
ref_tiles = ['T41SPR', 'T41RNQ', 'T41RNP', 'T41SNS', 'T41SPT', 'T41SPS', 'T41RPQ', 'T41SNR', 'T41RPP', 'T41SNT']
- Define the reference bands list (we want to select)
In [17]:
b_list = ['B4','dvi']
- Select from metadataframe a subframe of results related to more tiles on the different dates:
In [18]:
mydate_results = results[(results.apply(lambda row: row['startdate'] in pd.to_datetime(mydates[0:2]) , axis=1)) &
((results.apply(lambda row: any(ref_tile in row['title'] for ref_tile in ref_tiles[0:2]), axis=1)) |
(results.apply(lambda row: any(b in row['title'] for b in b_list), axis=1)))]
mydate_results.head()
Out[18]:
enclosure | enddate | identifier | self | startdate | title | wkt | |
---|---|---|---|---|---|---|---|
0 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 01be9eec66cb994310e469e012a3a2bc670f92e8 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_2... | POLYGON((64.0648798642054 31.5308017399005,65.... |
2 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 03aa38a9b430af9f2b2b2cb47cae31788c0fd59c | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... |
3 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 0949c278e99a1b2244ddd0edb3be4e38c7b10e81 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_2... | POLYGON((63.9422209286482 31.5415371135011,64.... |
4 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 0cccedf5d105c98b1bce071eebd3639084ebc43d | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... |
5 | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | 13f7b9c57f80b8239602fa8643e850e2810dbe86 | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... |
- Create a new GeoDataFrame with the selected metadata plus the data values of the selected bands (set ‘all’ in order to select all the bands) cropped over the aoi_wkt1 area
In [19]:
myGDF = GeoDataFrame()
aoi_wkt = aoi_wkt1
for row in mydate_results.iterrows():
myGDF = myGDF.append(get_GDF_with_datavalues(row[1], user, api_key, b_list, True, bbox=aoi_wkt.bounds), ignore_index=True)
#myGDF = myGDF.append(get_GDF_with_datavalues(row[1], user, api_key, 'all', True, bbox=aoi_wkt.bounds), ignore_index=True)
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8A_B11_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B8A_B11_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting bands [4, 13] for S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPT_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting bands [4, 13] for S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNP_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNT_20180729T082519_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPQ_20180729T082519_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SNS_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_20180729T082519_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPS_20180729T082519_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2A_MSIL2A_20180729T061631_N0206_R034_T41RPP_20180729T082519_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting bands [4, 13] for S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNR_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPR_20180727T091520_DVI_BOA_B8A_B11_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B8A_B11_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B12_B4_B3.tif
Getting bands [4, 13] for S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T083305_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B11_B8A_B4.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B12_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNQ_20180727T091520_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNT_20180727T083305_DVI_BOA_B8_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPT_20180727T083305_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T083305_DVI_BOA_B8A_B4_B12.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RPQ_20180727T091520_DVI_BOA_B4_B3_B2.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41RNP_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SNS_20180727T091520_DVI_BOA_B8A_B4_B3.tif
Getting composite bands S2B_MSIL2A_20180727T062629_N0206_R077_T41SPS_20180727T091520_DVI_BOA_B8A_B4_B3.tif
*Or*
- Select from metadataframe a subframe of Vegetation Index results related to a single georeferenced tiles of interest and multiple datetimes
mydate_results = results[(results.apply(lambda row: row[‘startdate’] in pd.to_datetime(mydates[0:4]) , axis=1)) & (results.apply(lambda row: ref_tiles[0] in row[‘title’], axis=1))]
mydate_results.head()
- Create a new GeoDataFrame with the selected metadata plus the data values of the selected bands (set ‘all’ in order to select all the bands) cropped over the aoi_wkt2 area
myGDF = GeoDataFrame()
aoi_wkt = aoi_wkt2
for row in mydate_results.iterrows():
myGDF = myGDF.append(get_GDF_with_datavalues(row[1], user, api_key, ['dvi'], True, bbox=aoi_wkt.bounds), ignore_index=True)
#myGDF = myGDF.append(get_GDF_with_datavalues(row[1], user, api_key, 'all', True, bbox=aoi_wkt.bounds), ignore_index=True)
In [20]:
myGDF.head()
Out[20]:
band | data | enclosure | enddate | geo_transform | identifier | projection | self | startdate | title | wkt | xsize | ysize | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | DVI_BOA_B8A_B11_B2 | [[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0... | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | (598410.0, 10.0, 0.0, 3502780.0, 0.0, -10.0) | 01be9eec66cb994310e469e012a3a2bc670f92e8 | PROJCS["WGS 84 / UTM zone 41N",GEOGCS["WGS 84"... | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_2... | POLYGON((64.0648798642054 31.5308017399005,65.... | 1344 | 1538 |
1 | DVI_BOA_B8A_B11_B2 | [[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0... | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | (598410.0, 10.0, 0.0, 3502780.0, 0.0, -10.0) | 03aa38a9b430af9f2b2b2cb47cae31788c0fd59c | PROJCS["WGS 84 / UTM zone 41N",GEOGCS["WGS 84"... | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... | 1344 | 1538 |
2 | DVI_BOA_B8A_B4_B12 | [[[149.0, 165.0, 165.0, 171.0, 171.0, 161.0, 1... | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | (598410.0, 10.0, 0.0, 3502780.0, 0.0, -10.0) | 0949c278e99a1b2244ddd0edb3be4e38c7b10e81 | PROJCS["WGS 84 / UTM zone 41N",GEOGCS["WGS 84"... | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41SNR_2... | POLYGON((63.9422209286482 31.5415371135011,64.... | 1344 | 1538 |
3 | DVI_BOA_B8A_B4_B12 | [[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0... | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | (598410.0, 10.0, 0.0, 3502780.0, 0.0, -10.0) | 0cccedf5d105c98b1bce071eebd3639084ebc43d | PROJCS["WGS 84 / UTM zone 41N",GEOGCS["WGS 84"... | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... | 1344 | 1538 |
4 | DVI_BOA_B11_B8A_B4 | [[[0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0... | https://store.terradue.com/better-satcen-00004... | 2018-07-29 06:16:31.024 | (598410.0, 10.0, 0.0, 3502780.0, 0.0, -10.0) | 13f7b9c57f80b8239602fa8643e850e2810dbe86 | PROJCS["WGS 84 / UTM zone 41N",GEOGCS["WGS 84"... | https://catalog.terradue.com//better-satcen-00... | 2018-07-29 06:16:31.024 | S2A_MSIL2A_20180729T061631_N0206_R034_T41RNQ_2... | POLYGON((63.9421873560077 30.6397709351317,64.... | 1344 | 1538 |
Plotting data and exporting geotif¶
Plot data¶
Choose DVI_BOA bands to plot: we are selecting the DVI_BOA data related to the first two reference tiles
In [21]:
bands = list()
mybands = myGDF[myGDF.apply(lambda row: any(ref_tile in row['title'] for ref_tile in ref_tiles[0:2]), axis=1) &
myGDF.apply(lambda row: 'DVI_BOA.tif' in row['title'], axis=1)]
for i,row in enumerate(mybands.iterrows()):
title = '%s (band %s)'%(row[1]['title'],row[1]['band'])
bmin = row[1]['data'].min()
bmax = row[1]['data'].max()
b = dict()
if bmin != bmax:
b[title] = (row[1]['data'] - bmin)/(bmax - bmin) * 255
else:
b[title] = row[1]['data']
bands.append(b)
numbands = len(bands)
print numbands
8
- Plot the selected bands separately
In [22]:
fig = plt.figure(figsize=(20,20))
m = len(b_list)
n = int(numbands/m)
j=0
for i,x in enumerate(bands):
for k,v in x.iteritems():
j=j+1
a = fig.add_subplot(n, m, j)
imgplot = plt.imshow(v.astype(np.uint8),
cmap='gray')
a.set_title(k)
plt.tight_layout()
fig = plt.gcf()
plt.show()
fig.clf()
plt.close()
Export a Geotiff¶
Export the band4 and dvi band related to the date time mydates[0] and ref_tiles[0,1] in a single Geotif
In [24]:
bands_number = len(b_list)
for i in range(2):
A = mybands[mybands.apply(lambda row: ref_tiles[i] in row['title'], axis=1) &
mybands.apply(lambda row: pd.to_datetime(mydates[0]) == row['startdate'], axis=1)]
cols = A['xsize'].values[0]
rows = A['ysize'].values[0]
geo_transform = A['geo_transform'].values[0]
projection = A['projection'].values[0]
drv = gdal.GetDriverByName('GTiff')
tiff_name = A['title'].values[0].split('.')[0]+'_'+'_'.join(b_list)+'.tif'
ds = drv.Create(tiff_name, cols, rows, bands_number, gdal.GDT_Float32)
ds.SetGeoTransform(geo_transform)
ds.SetProjection(projection)
ds.GetRasterBand(1).WriteArray(A['data'].values[0], 0, 0)
ds.GetRasterBand(2).WriteArray(A['data'].values[1], 0, 0)
ds.FlushCache()
- Select and Plot the composites bands
In [25]:
mybands_composistes = myGDF[myGDF.apply(lambda row: pd.to_datetime(mydates[0]) == row['startdate'], axis=1) &
myGDF.apply(lambda row: ref_tiles[0] in row['title'], axis=1) &
myGDF.apply(lambda row: 'DVI_BOA.tif' not in row['title'], axis=1)]
In [26]:
fig = plt.figure(figsize=(20,20))
title = '_'.join(mybands_composistes['title'][0].split('_')[0:7])
fig.suptitle('%s : DVI composites bands' %title, fontsize=16)
for index,row in enumerate(mybands_composistes.iterrows()):
index=index+1
rgb_uint8 = np.dstack(row[1]['data']).astype(np.uint8)
sub_fig = fig.add_subplot(2, 4, index)
img = Image.fromarray(rgb_uint8)
imgplot = plt.imshow(img)
sub_fig.set_title(row[1]['band'])
Download functionalities¶
Download a product¶
- Define the download function
In [27]:
def get_product(url, dest, api_key):
request_headers = {'X-JFrog-Art-Api': api_key}
r = requests.get(url, headers=request_headers)
open(dest, 'wb').write(r.content)
return r.status_code
- Get the reference download endpoint for the product related to a selected date (second date)
In [28]:
enclosure = myGDF[(myGDF['startdate'] == mydates[0])]['enclosure'].values[0]
enclosure
Out[28]:
'https://store.terradue.com/better-satcen-00004/_results/workflows/ec_better_ewf_satcen_01_03_01_ewf_satcen_01_03_01_0_7/run/26eddb46-20a5-11e9-93f5-0242ac11000f/0002389-181221095105003-oozie-oozi-W/c0f7cfbd-60e6-4842-8c2b-2c8ec3211f61/S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8A_B11_B2.tif'
In [29]:
output_name = myGDF[(myGDF['startdate'] == mydates[0])]['title'].values[0]
output_name
Out[29]:
'S2A_MSIL2A_20180729T061631_N0206_R034_T41SPR_20180729T082519_DVI_BOA_B8A_B11_B2.tif'
In [30]:
get_product(enclosure,
output_name,
api_key)
Out[30]:
200
Bulk Download¶
- Define the bulk download function
In [31]:
def get_product_bulk(row, api_key):
return get_product(row['enclosure'],
row['title'],
api_key)
- Download all the products related to a chosen date (first date) and a chosen reference tile
In [32]:
date_results = myGDF[(myGDF['startdate'] == mydates[0]) &
(myGDF.apply(lambda row: ref_tiles[0] in row['title'], axis=1))]
date_results.apply(lambda row: get_product_bulk(row, api_key), axis=1)
Out[32]:
0 200
12 200
14 200
32 200
38 200
39 200
41 200
45 200
64 200
dtype: int64