Working With COG¶
For this demo we will use the new DigitalGlobe OpenData
dataset https://www.digitalglobe.com/ecosystem/open-data
Requirements¶
- folium
- httpx
pip install folium httpx
In [1]:
Copied!
# Uncomment this line if you need to install the dependencies
# !pip install folium httpx
# Uncomment this line if you need to install the dependencies
# !pip install folium httpx
In [2]:
Copied!
import json
import httpx
from folium import Map, TileLayer
%matplotlib inline
import json
import httpx
from folium import Map, TileLayer
%matplotlib inline
In [3]:
Copied!
titiler_endpoint = "https://titiler.xyz" # Developmentseed Demo endpoint. Please be kind.
url = "https://opendata.digitalglobe.com/events/mauritius-oil-spill/post-event/2020-08-12/105001001F1B5B00/105001001F1B5B00.tif"
titiler_endpoint = "https://titiler.xyz" # Developmentseed Demo endpoint. Please be kind.
url = "https://opendata.digitalglobe.com/events/mauritius-oil-spill/post-event/2020-08-12/105001001F1B5B00/105001001F1B5B00.tif"
Get COG Info¶
In [4]:
Copied!
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/info",
params = {
"url": url,
}
).json()
bounds = r["bounds"]
print(r)
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/info",
params = {
"url": url,
}
).json()
bounds = r["bounds"]
print(r)
{'bounds': [57.664053823239804, -20.55473177712791, 57.84021477996238, -20.25261582755764], 'minzoom': 10, 'maxzoom': 18, 'band_metadata': [['b1', {}], ['b2', {}], ['b3', {}]], 'band_descriptions': [['b1', ''], ['b2', ''], ['b3', '']], 'dtype': 'uint8', 'nodata_type': 'Mask', 'colorinterp': ['red', 'green', 'blue'], 'driver': 'GTiff', 'count': 3, 'width': 38628, 'height': 66247, 'overviews': [2, 4, 8, 16, 32, 64, 128]}
Get COG Metadata¶
In [5]:
Copied!
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/statistics",
params = {
"url": url,
}
).json()
print(json.dumps(r, indent=4))
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/statistics",
params = {
"url": url,
}
).json()
print(json.dumps(r, indent=4))
{ "b1": { "min": 0.0, "max": 255.0, "mean": 36.94901407469342, "count": 574080.0, "sum": 21211690.0, "std": 48.282133573955264, "median": 3.0, "majority": 1.0, "minority": 246.0, "unique": 256.0, "histogram": [ [ 330584.0, 54820.0, 67683.0, 57434.0, 30305.0, 14648.0, 9606.0, 5653.0, 2296.0, 1051.0 ], [ 0.0, 25.5, 51.0, 76.5, 102.0, 127.5, 153.0, 178.5, 204.0, 229.5, 255.0 ] ], "valid_percent": 93.75, "masked_pixels": 38272.0, "valid_pixels": 574080.0, "percentile_2": 0.0, "percentile_98": 171.0 }, "b2": { "min": 0.0, "max": 255.0, "mean": 57.1494356187291, "count": 574080.0, "sum": 32808348.0, "std": 56.300819175100656, "median": 37.0, "majority": 5.0, "minority": 0.0, "unique": 256.0, "histogram": [ [ 271018.0, 34938.0, 54030.0, 69429.0, 70260.0, 32107.0, 29375.0, 9697.0, 2001.0, 1225.0 ], [ 0.0, 25.5, 51.0, 76.5, 102.0, 127.5, 153.0, 178.5, 204.0, 229.5, 255.0 ] ], "valid_percent": 93.75, "masked_pixels": 38272.0, "valid_pixels": 574080.0, "percentile_2": 5.0, "percentile_98": 180.0 }, "b3": { "min": 0.0, "max": 255.0, "mean": 51.251764562430324, "count": 574080.0, "sum": 29422613.0, "std": 39.65505035854822, "median": 36.0, "majority": 16.0, "minority": 252.0, "unique": 254.0, "histogram": [ [ 203263.0, 150865.0, 104882.0, 42645.0, 30652.0, 25382.0, 12434.0, 2397.0, 1097.0, 463.0 ], [ 0.0, 25.5, 51.0, 76.5, 102.0, 127.5, 153.0, 178.5, 204.0, 229.5, 255.0 ] ], "valid_percent": 93.75, "masked_pixels": 38272.0, "valid_pixels": 574080.0, "percentile_2": 14.0, "percentile_98": 158.0 } }
Display Tiles¶
In [6]:
Copied!
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
}
).json()
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=13
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="DigitalGlobe OpenData"
).add_to(m)
m
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
}
).json()
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=13
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="DigitalGlobe OpenData"
).add_to(m)
m
Out[6]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Work with non-byte data¶
In [7]:
Copied!
url = "https://data.geo.admin.ch/ch.swisstopo.swissalti3d/swissalti3d_2019_2573-1085/swissalti3d_2019_2573-1085_0.5_2056_5728.tif"
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/info",
params = {
"url": url,
}
).json()
print(r)
print("Data is of type:", r["dtype"])
# This dataset has statistics metadata
minv, maxv = r["band_metadata"][0][1]["STATISTICS_MINIMUM"], r["band_metadata"][0][1]["STATISTICS_MAXIMUM"]
print("With values from ", minv, "to ", maxv)
url = "https://data.geo.admin.ch/ch.swisstopo.swissalti3d/swissalti3d_2019_2573-1085/swissalti3d_2019_2573-1085_0.5_2056_5728.tif"
# Fetch File Metadata to get min/max rescaling values (because the file is stored as float32)
r = httpx.get(
f"{titiler_endpoint}/cog/info",
params = {
"url": url,
}
).json()
print(r)
print("Data is of type:", r["dtype"])
# This dataset has statistics metadata
minv, maxv = r["band_metadata"][0][1]["STATISTICS_MINIMUM"], r["band_metadata"][0][1]["STATISTICS_MAXIMUM"]
print("With values from ", minv, "to ", maxv)
{'bounds': [7.090624928537461, 45.91605844102821, 7.1035698381384185, 45.92509300025415], 'minzoom': 15, 'maxzoom': 18, 'band_metadata': [['b1', {'STATISTICS_COVARIANCES': '10685.98787505646', 'STATISTICS_EXCLUDEDVALUES': '-9999', 'STATISTICS_MAXIMUM': '2015.0944824219', 'STATISTICS_MEAN': '1754.471184271', 'STATISTICS_MINIMUM': '1615.8128662109', 'STATISTICS_SKIPFACTORX': '1', 'STATISTICS_SKIPFACTORY': '1', 'STATISTICS_STDDEV': '103.37305197708'}]], 'band_descriptions': [['b1', '']], 'dtype': 'float32', 'nodata_type': 'Nodata', 'colorinterp': ['gray'], 'driver': 'GTiff', 'count': 1, 'width': 2000, 'height': 2000, 'overviews': [2, 4, 8], 'nodata_value': -9999.0} Data is of type: float32 With values from 1615.8128662109 to 2015.0944824219
In [9]:
Copied!
# We could get the min/max values using the statistics endpoint
r = httpx.get(
f"{titiler_endpoint}/cog/statistics",
params = {
"url": url,
}
).json()
print(json.dumps(r["b1"], indent=4))
# We could get the min/max values using the statistics endpoint
r = httpx.get(
f"{titiler_endpoint}/cog/statistics",
params = {
"url": url,
}
).json()
print(json.dumps(r["b1"], indent=4))
{ "min": 1615.815185546875, "max": 2015.094482421875, "mean": 1754.49072265625, "count": 1048576.0, "sum": 1839716864.0, "std": 103.40494682173009, "median": 1721.2724609375, "majority": 1957.414794921875, "minority": 1615.815185546875, "unique": 812465.0, "histogram": [ [ 165709.0, 254204.0, 150889.0, 102835.0, 88631.0, 78677.0, 68163.0, 49590.0, 56651.0, 33227.0 ], [ 1615.815185546875, 1655.7431640625, 1695.6710205078125, 1735.5989990234375, 1775.52685546875, 1815.454833984375, 1855.3828125, 1895.3106689453125, 1935.2386474609375, 1975.16650390625, 2015.094482421875 ] ], "valid_percent": 100.0, "masked_pixels": 0.0, "valid_pixels": 1048576.0, "percentile_2": 1626.8169555664062, "percentile_98": 1987.6614990234375 }
Display Tiles¶
Note: By default if the metadata has min/max
statistics, titiler will use those to rescale the data
In [12]:
Copied!
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
).add_to(m)
m
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
).add_to(m)
m
Out[12]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Apply ColorMap
Now that the data is rescaled to byte values (0 -> 255) we can apply a colormap
In [13]:
Copied!
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
"rescale": f"{minv},{maxv}",
"colormap_name": "terrain"
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
).add_to(m)
m
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
"rescale": f"{minv},{maxv}",
"colormap_name": "terrain"
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
).add_to(m)
m
Out[13]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Apply non-linear colormap (intervals)
see https://cogeotiff.github.io/rio-tiler/colormap/#intervals-colormaps
In [14]:
Copied!
import json
cmap = json.dumps(
[
# ([min, max], [r, g, b, a])
([0, 1500], [255,255,204, 255]),
([1500, 1700], [161,218,180, 255]),
([1700, 1900], [65,182,196, 255]),
([1900, 2000], [44,127,184, 255]),
([2000, 3000], [37,52,148, 255]),
]
)
# https://colorbrewer2.org/#type=sequential&scheme=YlGnBu&n=5
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
"colormap": cmap
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
aod_layer = TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
)
aod_layer.add_to(m)
m
import json
cmap = json.dumps(
[
# ([min, max], [r, g, b, a])
([0, 1500], [255,255,204, 255]),
([1500, 1700], [161,218,180, 255]),
([1700, 1900], [65,182,196, 255]),
([1900, 2000], [44,127,184, 255]),
([2000, 3000], [37,52,148, 255]),
]
)
# https://colorbrewer2.org/#type=sequential&scheme=YlGnBu&n=5
r = httpx.get(
f"{titiler_endpoint}/cog/tilejson.json",
params = {
"url": url,
"colormap": cmap
}
).json()
bounds = r["bounds"]
m = Map(
location=((bounds[1] + bounds[3]) / 2,(bounds[0] + bounds[2]) / 2),
zoom_start=r["minzoom"] + 1
)
aod_layer = TileLayer(
tiles=r["tiles"][0],
opacity=1,
attr="Swisstopo"
)
aod_layer.add_to(m)
m
Out[14]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [ ]:
Copied!