[GIS] How to get data for multiple dates with OWSLib

owslibpythonsentinel-2wms

I am querying sentinel-2 data for multiple dates and want to get all results for these dates.

For example, I queried Arizona crater

from owslib.wms import WebMapService
import matplotlib.pyplot as plt
from pyproj import Proj, transform
import json

instance_id = ...

# in meters
region_size = (2000, 2000)

outProj = Proj(init='epsg:3857') # https://epsg.io/3857

inProj = Proj(init='epsg:4326') # https://epsg.io/4326

arizona_crater=(-111.02278,35.02722)

loc = arizona_crater

x, y = transform(inProj,outProj,loc[0],loc[1])

xupper = int(round(x - region_size[0] / 2))
xlower = int(round(x + region_size[0] / 2))
yupper = int(round(y - region_size[1] / 2))
ylower = int(round(y + region_size[1] / 2))

bbox = (xupper, yupper, xlower, ylower)
#bbox=(-11375363,3757392,-11372226,3759102)
size=region_size
#size=(1333, 716)

#wms_version = '1.1.1'
wms_version = '1.3.0'

wms = WebMapService('http://services.sentinel-hub.com/v1/wms/' + instance_id, version=wms_version)

print 'Type: ', wms.identification.type
print 'Title: ', wms.identification.title
print 'Contents: ', list(wms.contents)
print 'Operations: ',[op.name for op in wms.operations]
#print 'Capabilities: ', wms.getcapabilities()

info = wms.getfeatureinfo(
    layers=['1_NATURAL_COL0R'],
    srs='EPSG:3857',
    bbox = bbox,
    size = size,
    info_format='application/json',
    time='2015-01-01/2017-06-01/P1D',
    xy = (0, 0)
)
info = json.loads(info.read())

but getting only one result in feature list, although there are many

http://sentinel-s2-l1c.s3-website.eu-central-1.amazonaws.com/#tiles/12/S/VD/

How to get all results for all dates or query possible dates with OWSLib?

Best Answer

I don't think it is possible to do what you want with a GetFeatureInfo request. It doesn't matter whether you are constructing the request using OWSLib, or rolling by hand.

When you make a GetMap request for a single layer, then you will get back a single image. With the Sentinel data by default you get a mosaic of the most recent scenes for the Area of Interest (AOI). You can change the default behaviour by specifying a date or date range using the time parameter, you still only get one image in your layer, with the most recent scene from your time/range on top.

When you do a GetFeatureInfo you are requesting the feature information from a point location in your map, so in any one location you are querying the date of the scene used to create that part of the mosaic.

Illustration of GetMap request result

For example using the figure above. The black box is the GetMap request bounding box (AOI). You make a request for scenes that are from the winter months of 2015 using the time parameter. The WMS server finds two scenes (blue and red) and red is the most recent. The WMS generates an image which is composed of all the most recent scene, and most of the earlier scene (where it doesn't overlap the recent scene).

A GetFeatureInfo request in the white are returns no date, a GetFeatureInfo request in the hatched area returns the 'blue' date, and a GetFeatureInfo request in the red area returns the 'red' date, even in the area where the blue scene (geographically) underlies the red. That is the expected behaviour.

I had a look at the Sentinel hub docs on WMS parameters and also this video on the Sentinel Hub - WMS feature info functionality, and can't see any evidence that the WMS behaves any differently from how I have described it above.

EDIT

As pointed out by @grega-m in comments, the suggested way in the Sentinel FAQs to find available dates in a specific area is to use a WFS request with the following format:

http://services.sentinel-hub.com/v1/wfs/%3CINSTANCE_ID%3E?
 REQUEST=GetFeature&
 TYPENAMES=TILE&
 OUTPUTFORMAT=application/json&
 BBOX=1600965,5819367,1610917,5826373&
 TIME=2016-05-01/2017-01-25&
 version=2.0.0&
Related Question