Tutorials

Topography (topo module)

Presentation

This module provides a set of functions to deal with an EarthDataHub service and to analyze data to provide an array of elevation levels and a base altitude.

Details

The main use case that is detailed here using the following steps as explained. Any of these steps can be used independently for a different use case than the current one.

The steps are:

  • get_topo downloads a terrain elevation data array about the given EPSG:4326 bounding box. The source data raster is resampled to the given width X height.

  • levelize discretizes the raw DEM data to levels of the given thickness and re-level it to the minimal altitude of the data array.

Usage

The inputs are:

  • a bbox (bounding box) that represents the area that the output data will have to cover (the coordinates are expressed as a tuple (X0, Y0, X1, Y1) in EPSG:4326, the two points are the top left and bottom right points),

  • width X height dimension of the final data array,

  • a level height in meters.

Code example

An EARTHDATAHUB_API_KEY env variable must be provided when executing the get_topo function.

from topo import get_topo, levelize
import logging
import asyncio
import numpy

api_key = place here your API key

src_url = 'https://edh:{}@data.earthdatahub.destine.eu/copernicus-dem/GLO-90-v0.zarr'
logger=logging.Logger("test")
loop = asyncio.get_event_loop()
coords=0.8607384750612755,45.82539712786116,1.5793907023568927,45.32235342412723
dt=loop.run_until_complete(get_topo(logger, src_url, coords, 64, 64, api_key))
print(dt)

bl, dt=levelize(dt, 20)
print(bl, dt)
print(numpy.max(dt))

Run

Topography

Climate Zone (cl_zone module)

Presentation

This module provides a climate zone computation implementation of the Köppen-Geiger algorithm.

Details

This implementation of the Köppen-Geiger algorithm outputs one of the following values:

# A     Tropical
# Af    Tropical, Rainforest
# Am    Tropical, Monsoon
# Aw    Tropical, Savannah (with dry-winter)
# B     Arid
# BW    Arid, Desert
# BWh   Arid, Desert (Hot Desert)
# BWk   Arid, Desert (Cold Desert)
# BS    Arid, Steppe
# BSh   Arid, Steppe, Hot (Hot semi-arid)
# BSk   Arid, Steppe, Cold (Cold semi-arid)
# C     Temperate
# Cs    Temperate, Dry Summer
# Csa   Temperate, Dry Summer (Hot-summer Mediterranean)
# Csb   Temperate, Dry Summer (Warm-summer Mediterranean)
# Cw    Temperate, Dry winter
# Cwa   Temperate, Dry winter (Monsoon-influenced humid subtropical)
# Cwb   Temperate, Dry winter (Monsoon-influenced temperate oceanic)
# Cwc   Temperate, Dry winter (Monsoon-influenced subpolar oceanic)
# Cf    Temperate, Without Dry Season
# Cfa   Temperate, Without Dry Season, Hot Summer (Humid subtropical)
# Cfb   Temperate, Without Dry Season, Warm Summer (Temperate Oceanic)
# Cfc   Temperate, Without Dry Season, Cold Summer (Subpolar oceanic)
# D     Cold
# Ds    Cold, Dry summer
# Dsa   Cold, Dry summer (Mediterranean-influenced hot summer humid continental)
# Dsb   Cold, Dry summer (Mediterranean-influenced warm summer humid continental)
# Dsc   Cold, Dry summer (Mediterranean-influenced subartic)
# Dsd   Cold, Dry summer (Mediterranean-influenced extremely cold subartic)
# Dw    Cold, Dry winter
# Dwa   Cold, Dry winter (Monsoon-influenced hot summer humid continental)
# Dwb   Cold, Dry winter (Monsoon-influenced warm summer humid continental)
# Dwc   Cold, Dry winter (Monsoon-influenced subartic)
# Dwd   Cold, Dry winter (Monsoon-influenced extremly cold subartic)
# Df    Cold, Without Dry Season
# Dfa   Cold, Hot Summer (Hot-summer humid continental)
# Dfb   Cold, Warm Summer (Warm-summer humid continental)
# Dfc   Cold, Cold Summer (Subartic)
# Dfd   Cold, Very Cold Winter (Extremely Cold Subartic)
# E     Polar
# ET    Polar, Tundra
# EF    Polar, Frost (Ice Cap)

Usage

The input parameters are:

  • a 12-value data vector of average temperature in Celsius degrees for each month starting with January,

  • a 12-value data vector of cumulative precipitation in mm for each month starting with January.

Code example

import numpy
from cl_zone import compute_from_climate_data

# Define arbitrary value vectors
temps=numpy.asarray([12, 13, 17, 19, 24, 28, 30, 30, 26, 21, 15, 12])
precipitations=numpy.asarray([41, 27, 23, 54, 37, 21, 6, 21, 82, 63, 59, 36])

# Show the result
ret=compute_from_climate_data(temps, precipitations)
print(ret)
print(describe_climate_zone(ret))

Run

ClimateZone

Tiles to WMS (tiles_to_wms module)

Presentation

This module provides a set of functions to deal with a tile server (as per XYZ standard tile server) and provide a WMS like service based on it.

Details

The main use case that is detailed here using the following steps as explained. Any of these steps can be used independently for a different use case than the current one.

The steps are:

  • if Z is not provided, an auto zoom level is computed from the bounding box,

  • compute the covering tiles with margins: this mean select the min and max X and Y values of tiles covering the requested bounding box, covering_tiles(bbox, Z)

  • load all the tiles with the X and Y values selected in the previous step. This task is fully asynchronous: all the download tasks are run in parallel. This task returns the size in pixels (width and height) of each tile. load_tiles(logger, source_url, west_X, north_Y, east_X, south_Y, Z)

  • assemble the downloaded tiles in one full frame big raster. assemble_rasters(dss, nb_tiles_x, nb_tiles_y, size_tiles_x, size_tiles_y)

  • compute the margin in pixel number to remove from the full frame raster to the initial bounding box. compute_removable_margins(west_X, north_Y, east_X, south_Y, size_tiles_x, size_tiles_y)

  • clip the full frame raster for the margins computed in the previous step clip(source_raster, margin_north, margin_south, margin_west, margin_east)

  • return the clipped raster as the final product

Usage

The inputs are:

  • a endpoint to request for source tiles,

  • a bbox (bounding box) that represents the area that the output image will have to cover (the coordinates are expressed as a tuple (X0, Y0, X1, Y1) in EPSG:4326, the two points are the top left and bottom right points),

  • an optional zoom level as Z axis for source tiles - if Z is not provided, an auto zoom level is computed from the bounding box

Code example

import asyncio
import rasterio
import numpy
from tiles_to_wms import tiles_to_wms

async def main():
    # Provide a logger
    import logging
    logger = logging.getLogger(__name__)
    logging.basicConfig(level=logging.INFO)
    # Endpoint for Copernicus land cover tiles
    cop_sc="https://8q6ans41l4.execute-api.eu-west-1.amazonaws.com/cog/tiles/{Z}/{X}/{Y}?url=https%3A%2F%2Fvito-lcv.s3-eu-west-1.amazonaws.com%2Fglobal%2F2019%2Fcog-full_l0-colored-full.tif&resampling_method=mode"
    # Bounding box to request
    bbox_4326=(2.2,49,2.5,48.7)

    raster, width, height = await tiles_to_wms(logger, cop_sc, bbox_4326)

    # Write the result raster data into a TIF file
    with rasterio.open('output.tif', 'w', driver='GTiff', width=width, height=height, count=3, dtype=numpy.uint8) as dst:
        dst.write(raster)
loop = asyncio.get_event_loop()
loop.run_until_complete(main())

Run

TilesToWms

Water (water module)

Presentation

This module provides a set of functions to deal with a geoserver WMS service providing water bodies maps.

Details

The main use case that is detailed here using the following steps as explained. Any of these steps can be used independently for a different use case than the current one.

The steps are:

  • get_water_raster downloads the raster for the source WMS in the requested layer and in width X height dimension. The coordinates must be converted from EPSG:4326 into WMS specific format using bbox_4326_to_wms_bbox,

  • analyze_water_raster analyzes the raster and provides an array of boolean (True when the pixel has water). The water detection is triggered if the transparency pixel value is above a given threshold.

Usage

The inputs are:

  • a geoserver WMS endpoint to request for source image,

  • a geoserver layer name,

  • a bbox (bounding box) that represents the area that the output data will have to cover (the coordinates are expressed as a tuple (X0, Y0, X1, Y1) in EPSG:4326, the two points are the bottom left and top right points),

Code example

import asyncio
from water import WMS_URL, get_water_raster, bbox_4326_to_wms_bbox, analyze_water_raster

water_geoserver_base_url="https://game.demo.aas.atosfordestine.eu/geoserver/destine/"
layer='destine:hydropolys'
bbox=0.2850952148437498,43.32952880859375,2.39556884765625,44.38531494140625

raster=asyncio.run(get_water_raster(water_geoserver_base_url+WMS_URL.format(layer,64,64,bbox_4326_to_wms_bbox(bbox))))
data=analyze_water_raster(raster)
print(data)

Run

Water