Getting familiar with collections¶
The module name of the interactive processing library is inter. Starting from June 2018 all classes and functions of the interactive library can be called also without placing “inter.” before the class or function name.
In this tutorial you can get familiar with the concept of an image collection (group of images sharing some characteristics such as the sensor)
List available collections¶
In order to explore the available collections, the function collectionsSearch() can be used. It takes as first argument the string to search. To have the list of all the available collection simply pass the empty string:
collectionsSearch("Sentinel")
In the results printed by the call to collectionsSearch, the name of the collection is displayed togheter with the path descriptor of the collection (starting with “collections.”):
Sentinel-2 image collection
collections.EarthObservation.Copernicus.Sentinel2
Sentinel-1 mosaic
collections.Products.Mosaics.Global.S1_RGB_2016
Sentinel2 relative orbits published by ESA
collections.BaseData.GeographicalGridSystems.S2orbits
Military Grid Reference System
collections.BaseData.GeographicalGridSystems.MGRS
Ships data in CSV format from Sentinel-1 ships detection
collections.Products.Other.ShipsDetectedFromSentinel1
The function collectionsSearch() takes as second parameter a flag that, if True, displays a code snippet for each of the collections that satisfy the search:
collectionsSearch("Sentinel", True)
Gives as result:
Sentinel-2 image collection
******************************************************************************************************************
coll = Collection(collections.EarthObservation.Copernicus.Sentinel2)
coll = coll.filterOnGeoName("Munchen", 0, "DE")
coll = coll.filterOnDate(2015,1,1, 2016,12,31)
coll = coll.filterOn("cloudcover", "<=", 5)
coll = coll.filterOn("jrc_filepath", "<>", "")
coll.limit(4)
coll.sort("cloudcover",True)
p = coll.process().bands("B04","B03","B02",2)
map.addLayer(p.toLayer())
Sentinel-1 mosaic
******************************************************************************************************************
p = Collection(collections.Products.Mosaics.Global.S1_RGB_2016).process().bands("1","2","3").opacity(255)
map.addLayer(p.toLayer())
...
To use a collection, the function Collection() can be called by passing the path descriptor of the collection:
coll = Collections(collections.EarthObservation.Copernicus.Sentinel2)
Filter Sentinel2 (S2) collection¶
A S2 collection can be filtered based on different criteria:
- geographic location:
ImageCollection.filterOnGeoName()(Name, BufferDistance, Country) that filters all images that intersect a given geoname;- acquisition date:
ImageCollection.filterOnDate()(dateFrom, dateTo) that filters all images acquired between dateFrom and dateTo (date format is: year, month, day);- attributes:
ImageCollection.filterOn()(attributeName, operand, value) that filters all images according to some attribute of the collection. The list of attributes can be printed using the functionImageCollection.printAttributes(). The list of operands are: “<”, “<=”, “>”, “>=”, “=”, “<>” (not equal to), “LIKE”.
Example of geographic filter:
coll = Collections(collections.EarthObservation.Copernicus.Sentinel2)
cityName="Munchen"
enlargeByDegree= 0
countryCode="DE"
coll = coll.filterOnGeoName(cityName, enlargeByDegree, countryCode)
Example of filter on acquisition date:
coll = coll.filterOnDate(2015,1,1, 2016,12,31)
Example of attribute filter:
coll.filterOn("cloudcoverpercentage","<=",10)
Obtain information on image collections¶
Image collections can be sorted (ImageCollection.sort()) and limited (ImageCollection.limit()).
For instance, to sort in ascending mode on the cloud cover percentage of the images and to limit the collection to the first 4 images (those less cloudy), these lines of code can be used:
coll = coll.sort("cloudcoverpercentage")
coll.limit(4)
To list the number of selected images in a collection, use the function ImageCollection.count():
coll.count()
Any attribute of the images in the collection can be printed using the ImageCollection.printImages(), for instance:
coll.printImages("orbitnumber")
will print the orbit numbers of all the images in the collection (not only the distinct ones).
To print the file path to all available Sentinel-2B images acquired in July 2017 over Ispra:
coll = Collections(collections.EarthObservation.Copernicus.Sentinel2)
cityName="Ispra"
enlargeByDegree= 0
countryCode="IT"
coll = coll.filterOnGeoName(cityName, enlargeByDegree, countryCode)
coll = coll.filterOnDate(2017,7,1, 2017,7,31)
coll.printImages('jrc_filepath')
Visualising image collections¶
We first open an interactive map used as basis for the interactive visualisation and processing of a collection:
map = Map()
map
We then need to create a ImageProcess, that will act on our collection.
The most basic process is the selection of a single band (using the function ImageProcess.band()), which can be added as a layer to the map using Map.addLayer(). As an argument we need to convert the process to a layer, using the function ImageProcess.toLayer():
p1 = coll.process().band("B03")
map.clear()
map.zoomToImageExtent(p1)
tlayer1 = map.addLayer(p1.toLayer())
The band “B03” was selected as an example here. Note that you can use the function ImageCollection.printBandNames() to list all the band names of a collection.
To display an RGB composition, use the function ImageProcess.bands(). For instance, by selecting the bands B04, B03 and B02, a “true color” image is created:
p1 = coll.process().bands("B04","B03","B02")
map.clear()
map.zoomToImageExtent(p1)
tlayer1 = map.addLayer(p1.toLayer())
The visibility of a layer can be toggled using the layer’s attribute visible. To make a layer visible use:
tlayer1.visible = True
To make a layre not visible use:
tlayer1.visible = False
We can display NDVI (Normalized Difference Vegetation Index) that is automatically calculated as: (NIR-RED)/(NIR+RED), i.e. (B08-B04)/(B08+B04) for Sentinel 2 products. Values are in the [0,65535] range:
p1 = coll.process().NDVI()
map.clear()
map.zoomToImageExtent(p1)
tlayer1 = map.addLayer(p1.toLayer())
Simple band arithmetic¶
To calculate the average of two S2 bands, B03 and B04, use the ImageProcess.arith() function of a process. The result is a new process, which can thus be re-used for another ImageProcess.arith() function:
p1 = coll.process().band("B04")
p2 = coll.process().band("B03")
p1 = p1.arith("+",p2).arith("*",0.5)
map.clear()
tlayer1 = map.addLayer(p1.toLayer())
Another simple band arithmetic is to mask values. For instance, we can use the result process p1 to mask out all values less than 600, using the function ImageProcess.maskLT(). For better visualisation, we scale the result from value 600 to 1500 to the range [0-255] using the function ImageProcess.scale(). Then we display the result by adding the result to the map, using the function Map.addLayer():
p1.maskLT(600)
p1.scale(600,1500)
map.clear()
tlayer1 = map.addLayer(p1.toLayer())
Using color palettes and transparency¶
Instead of a grey scale image, a single band image can also be mapped to a color palette using a color scheme (ImageProcess.colorScheme()) or even a custom color palette (ImageProcess.colorCustom()).
To get the list of color schemes, use the function ImageCollection.printColorSchemes()
We can get a colored map of the masked layer in the section above as follows:
p1.colorScheme("Oranged_mixed")
map.clear()
tlayer1 = map.addLayer(p1.toLayer())
For a custom color scheme, use the function ImageProcess.colorCustom(), by providing a list of colors. The names of the colors can be printed with the function ImageCollection.printNamedColors(). Other colors can be used with the standard web syntax: “#RRGGBB”. For example: yellow is represented by the code: “#FFFF00”. For the transparency, you can use the function ImageProcess.opacity(). It expects a single argument, ranging from 0 (fully transparent) to 255 (fully opaque).
Please note that when more than one color is provided, the first is mapped to minimum value of the scaling range, the last to maximum value and a linear gradient is calculated in between.
As an example, to display the NDVI values using a custom color palette:
p1 = coll.process().NDVI().maskLT(54000).scale(54000,60000)
colors = ["Red", "Yellow", "Green"]
p1.colorCustom(colors)
map.clear()
map.zoomToImageExtent(p1)
tlayer1 = map.addLayer(p1.toLayer())