Tiles¶

In [ ]:

import os
import tempfile
import time
import warnings
from pathlib import Path

import contextily as cx
import geopandas as gpd
import pyproj
import xyzservices.providers as xyz
from dotenv import load_dotenv

import open_geodata as geo

import os import tempfile import time import warnings from pathlib import Path import contextily as cx import geopandas as gpd import pyproj import xyzservices.providers as xyz from dotenv import load_dotenv import open_geodata as geo

---

Credenciais¶

In [ ]:

load_dotenv()

MAPTILER_TOKEN = os.getenv('MAPTILER_TOKEN')
JAWG_TOKEN = os.getenv('JAWG_TOKEN')
MAPBOX_TOKEN = os.getenv('MAPBOX_TOKEN')

load_dotenv() MAPTILER_TOKEN = os.getenv('MAPTILER_TOKEN') JAWG_TOKEN = os.getenv('JAWG_TOKEN') MAPBOX_TOKEN = os.getenv('MAPBOX_TOKEN')

---

Definição de Área de Interesse¶

Obtem uma área para que façamos o download.

In [ ]:

# Instancia classe IBGE
sfb = geo.br.sicar.SICAR(
    tesseract_path='C:\\Program Files\\Tesseract-OCR\\tesseract.exe'
)

# Download Data
with tempfile.TemporaryDirectory() as temp_dir:
    # Cria o caminho temporário em formato Path
    temp_path = Path(temp_dir)

    # Download
    sfb.download_data(
        sigla_estado=geo.br.sicar.State.MA,
        layer=geo.br.sicar.Polygon.AREA_PROPERTY,
        output_path=temp_path,
    )

    # Lista arquivos baixados
    list_files = list(temp_path.glob('*'))
    print(list_files)

    # Lê o shapefile diretamente do zip
    gdf = gpd.read_file(filename=list_files[0])

# Filtra propriedades com área definida em um intervalo de áreas (em hectares)
gdf = gdf[(gdf['num_area'] >= 20) & (gdf['num_area'] <= 21)]

# Obtem Linha Randômica
gdf = gdf.sample(n=1)

# Exibe informações do GeoDataFrame filtrado
gdf.info()
gdf.head(2)

# Instancia classe IBGE sfb = geo.br.sicar.SICAR( tesseract_path='C:\\Program Files\\Tesseract-OCR\\tesseract.exe' ) # Download Data with tempfile.TemporaryDirectory() as temp_dir: # Cria o caminho temporário em formato Path temp_path = Path(temp_dir) # Download sfb.download_data( sigla_estado=geo.br.sicar.State.MA, layer=geo.br.sicar.Polygon.AREA_PROPERTY, output_path=temp_path, ) # Lista arquivos baixados list_files = list(temp_path.glob('*')) print(list_files) # Lê o shapefile diretamente do zip gdf = gpd.read_file(filename=list_files[0]) # Filtra propriedades com área definida em um intervalo de áreas (em hectares) gdf = gdf[(gdf['num_area'] >= 20) & (gdf['num_area'] <= 21)] # Obtem Linha Randômica gdf = gdf.sample(n=1) # Exibe informações do GeoDataFrame filtrado gdf.info() gdf.head(2)

In [ ]:

gdf.explore(column='des_condic', tiles='Esri.WorldImagery')

gdf.explore(column='des_condic', tiles='Esri.WorldImagery')

---

Download Images: Contextily¶

• https://contextily.readthedocs.io/en/latest/

In [ ]:

# Métricas Web Mercator
gdf = gdf.to_crs(epsg=3857)

# Buffer de 100m
gdf = gdf.buffer(distance=100)

# Extract bounding box in WGS84
w, s, e, n = gdf.to_crs(epsg=4326).total_bounds

# Extract bounding box in Web Mercator
# w, s, e, n = gdf.to_crs(epsg=3857).total_bounds
print(w, s, e, n)

# Métricas Web Mercator gdf = gdf.to_crs(epsg=3857) # Buffer de 100m gdf = gdf.buffer(distance=100) # Extract bounding box in WGS84 w, s, e, n = gdf.to_crs(epsg=4326).total_bounds # Extract bounding box in Web Mercator # w, s, e, n = gdf.to_crs(epsg=3857).total_bounds print(w, s, e, n)

Ajusta uma pasta temporária que irá armazenar os dados que serão baixados.

In [ ]:

temp_path = Path(tempfile.gettempdir()) / 'open_geodata' / 'raster' / 'tiles'
temp_path.mkdir(exist_ok=True)
temp_path

temp_path = Path(tempfile.gettempdir()) / 'open_geodata' / 'raster' / 'tiles' temp_path.mkdir(exist_ok=True) temp_path

Provedores¶

Lista alguns provedores de dados. Pode ser usado tiles customizados e também os provedores do pacote xyz.

In [ ]:

# Mapbox
# xyz.MapBox(id='mapbox/satellite-v9', accessToken=MAPBOX_TOKEN)

# Mapbox # xyz.MapBox(id='mapbox/satellite-v9', accessToken=MAPBOX_TOKEN)

In [ ]:

providers = {
    'Google': 'https://mt0.google.com/vt?lyrs=s&x={x}&y={y}&z={z}',
    'ESRI': 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}',
    #'Mapnik': xyz.OpenStreetMap.Mapnik,
    #'BlueMarble': xyz.NASAGIBS.BlueMarble,
    'Mapbox': xyz.MapBox(id='mapbox/satellite-v9', accessToken=MAPBOX_TOKEN),
}

providers = { 'Google': 'https://mt0.google.com/vt?lyrs=s&x={x}&y={y}&z={z}', 'ESRI': 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}', #'Mapnik': xyz.OpenStreetMap.Mapnik, #'BlueMarble': xyz.NASAGIBS.BlueMarble, 'Mapbox': xyz.MapBox(id='mapbox/satellite-v9', accessToken=MAPBOX_TOKEN), }

Faz o download...

In [ ]:

# Download file
for zoom in [
    5,
    10,
    12,
    15,
    16,
    17,
    18,  # Não deu pra ESRI
    19,  # Não deu pra ESRI
    20,  # Não deu pra ESRI
    21,  # Não deu pra ESRI
]:
    # Message
    print(f'Baixando os dados para o zoom {zoom:02d}')

    for key, value in providers.items():
        print(f'> {key}')
        provider = key
        source = value

        # Create Path
        provider_path = temp_path / f'{provider.lower()}'
        provider_path.mkdir(exist_ok=True)

        try:
            # Download ArcGIS
            img, ext = cx.bounds2raster(
                w=w,
                s=s,
                e=e,
                n=n,
                path=provider_path / f'{provider.lower()}_zoom_{zoom:02d}.tif',
                # zoom=f'{zoom}', # TODO: Abrir PR
                zoom=zoom,
                ll=True,
                source=source,
                max_retries=3,
                wait=2,
            )

        except:
            warnings.warn(
                f'> Não foi possível baixar imagem de {provider}, para o zoom {zoom}'
            )

    # Pausa
    time.sleep(2)

# Download file for zoom in [ 5, 10, 12, 15, 16, 17, 18, # Não deu pra ESRI 19, # Não deu pra ESRI 20, # Não deu pra ESRI 21, # Não deu pra ESRI ]: # Message print(f'Baixando os dados para o zoom {zoom:02d}') for key, value in providers.items(): print(f'> {key}') provider = key source = value # Create Path provider_path = temp_path / f'{provider.lower()}' provider_path.mkdir(exist_ok=True) try: # Download ArcGIS img, ext = cx.bounds2raster( w=w, s=s, e=e, n=n, path=provider_path / f'{provider.lower()}_zoom_{zoom:02d}.tif', # zoom=f'{zoom}', # TODO: Abrir PR zoom=zoom, ll=True, source=source, max_retries=3, wait=2, ) except: warnings.warn( f'> Não foi possível baixar imagem de {provider}, para o zoom {zoom}' ) # Pausa time.sleep(2)

In [ ]: