Fabio/terrain2
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

first commit

Browse source
This commit is contained in:
Fabio Micheluz 2026-02-12 17:21:04 +01:00
commit 9603735519
17 changed files with 1665 additions and 0 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

188
README.md Normal file
View file

@ -0,0 +1,188 @@
# Creazione MbTiles per tileserver GL utilizzando Copernicus glo30
## Download dei file
🚀 Scaricare i dati da mirror Copernicus GLO‑30
```
aws s3 sync s3://copernicus-dem-30m ./glo30 --no-sign-request
```
Crea glo30 dir con tutto il mondo
ma non servono tutti i dati quindi questi sono i tre programmi python per
scaricare italia europa e mondo in base a quello che vuoi ottenere
```
python3 i.py # italia
python3 e.py # europa
python3 w.py # mondo
```
## Trasformare i DEM geotiff in tiff RGB
con la funzione rgbify_batch.py puoi creare le versioni tif in formato RGB
```
python3 rgbify_batch.py glo30_italia glo30_italia_rgb
```
il comando per la singola tiff sarebbe
```
rio rgbify -b -10000 -i 0.1 a.tif a_rgb.tif
```
## Trasformare gli RGB in tiles
corretto
```
python3 make_tiles_rio_tiler.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif tiles_x
```
## Trasformare i tiles in mbtiles
```
mb-util --scheme=xyz --image_format=png tiles_x/ x.mbtiles
```
## Inserire i metadata in mbtiles
```
python3 set_mbtiles_metadata.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif xxx.mbtiles "Terrain RGB Adamello"
```
## Verificare
```
python3 verify_mbtiles_metadata.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif xxx.mbtiles
python3 verify_and_fix_mbtiles.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif xxx.mbtiles
python3 verify_bounds_tileservergl.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif xxx.mbtiles
```
## multipli
multili con file
```
python3 make_tiles_rio_tiler.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif tiles_m
python3 make_tiles_rio_tiler.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E010_00_DEM_rgb.tif tiles_m
mb-util --scheme=xyz --image_format=png tiles_m/ m.mbtiles
python3 set_mbtiles_metadata_multi.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E010_00_DEM_rgb.tif m.mbtiles
```
multipli con dir
```
python3 make_tiles_rio_tiler_multi.py dir_tif tiles_mm
mb-util --scheme=xyz --image_format=png tiles_mm/ mm.mbtiles
python3 set_mbtiles_metadata_dir.py dir_tif mm.mbtiles
```
ultima idea
```
python3 make_tiles_rio_tiler.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E010_00_DEM_rgb.tif tiles_x
python3 make_tiles_rio_tiler.py glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif tiles_y
rsync -av tiles_x/ tiles_m/
rsync -av --ignore-existing tiles_y/ tiles_m/
mb-util --scheme=xyz --image_format=png tiles_m/ m.mbtiles
python3 set_mbtiles_metadata_dir.py dir_tif m.mbtiles
scp m.mbtiles orangepi@192.168.1.4:/home/nvme/dockerdata/tileserver1/terrain.mbtiles
```
## Controllo
Monte Adamello
Latitude: 46° 9' 19" N
Longitude: 10° 29' 47" E
Lat/Long (dec): 46.1554,10.4966
Köppen climate type: ET : Tundra
Elevation: 3,539m
elevazione=-10000+(R*256^2+G*256+B)*0,1
estrarre l'altezza del monte dal geotiff
```
gdallocationinfo -wgs84 -valonly glo30_italia/Copernicus_DSM_COG_10_N46_00_E010_00_DEM.tif 10.4966 46.1554
```
è da come risultato
3496.90673828125
estrarre l'altezza del monte dal rgb
```
gdallocationinfo -wgs84 glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E010_00_DEM_rgb.tif 10.4966 46.1554
```
e da come report
```
Report:
Location: (1788P,3041L)
Band 1:
Value: 2
Band 2:
Value: 15
Band 3:
Value: 57
```
utilizzando la funzione python che calcola altutudine da RGB si ottiene 3496.90
```
python3 rgb2alt.py 2 15 57
```
verifichiamo con mbtiles
```
gdallocationinfo -wgs84 -valonly \
-b 1 -b 2 -b 3 \
aa.mbtiles 10.4966 46.1554
```
```
gdallocationinfo -wgs84 -valonly \
-b 1 -b 2 -b 3 \
x.mbtiles 10.4966 46.1554
```
Punta Penia
Dolomiti (3.343 m), sono approssimativamente:
Latitudine: 46° 26' 03.84" N (46.43429)
Longitudine: 11° 50' 58.58" E (11.85051)
```
gdallocationinfo -wgs84 -valonly glo30_italia/Copernicus_DSM_COG_10_N46_00_E011_00_DEM.tif 11.85051 46.43429
3236.89526367188
gdallocationinfo -wgs84 glo30_italia_rgb/Copernicus_DSM_COG_10_N46_00_E011_00_DEM_rgb.tif 11.85051 46.43429
Report:
Location: (3062P,2037L)
Band 1:
Value: 2
Band 2:
Value: 5
Band 3:
Value: 16
python3 rgb2alt.py 2 5 16
3236.800000000001
gdallocationinfo -wgs84 -valonly -b 1 -b 2 -b 3 x.mbtiles 11.85051 46.43429
2
6
68
python3 rgb2alt.py 2 6 68
3267.6
```
gdallocationinfo -wgs84 -valonly \
-b 1 -b 2 -b 3 \
aa.mbtiles 10.4966 46.1554
2
15
147
quindi 3505,9

96
bounds.py Normal file
View file

@ -0,0 +1,96 @@
#!/usr/bin/env python3
import argparse, sqlite3, sys
def fetch_one(cur, sql, args=()):
row = cur.execute(sql, args).fetchone()
return row[0] if row else None
def main():
ap = argparse.ArgumentParser(
description="Normalizza bounds (W,S,E,N) e imposta center/minzoom/maxzoom in un MBTiles."
)
ap.add_argument("mbtiles", help="Percorso al file .mbtiles")
ap.add_argument("--zoom", type=int, default=None,
help="Zoom del center da impostare (default: usa center esistente o media min/max dalle tile)")
ap.add_argument("--minzoom", type=int, default=None,
help="Forza minzoom (default: dai livelli in 'tiles')")
ap.add_argument("--maxzoom", type=int, default=None,
help="Forza maxzoom (default: dai livelli in 'tiles')")
args = ap.parse_args()
conn = sqlite3.connect(args.mbtiles)
cur = conn.cursor()
# --- 1) leggi bounds correnti
b_raw = fetch_one(cur, "SELECT value FROM metadata WHERE name='bounds';")
if not b_raw:
sys.exit("ERRORE: 'bounds' mancante nella tabella metadata.")
try:
x1, y1, x2, y2 = [float(t) for t in b_raw.split(",")]
except Exception:
sys.exit(f"ERRORE: bounds non parseable: {b_raw}")
# normalizza ordine W,S,E,N (min/max su lon e lat)
W, E = sorted((x1, x2))
S, N = sorted((y1, y2))
bounds = f"{W:.6f},{S:.6f},{E:.6f},{N:.6f}"
# --- 2) min/max zoom: da tiles salvo override
minz = fetch_one(cur, "SELECT MIN(zoom_level) FROM tiles;")
maxz = fetch_one(cur, "SELECT MAX(zoom_level) FROM tiles;")
# fallback se tiles non esiste o vuota
try:
minz = int(minz) if minz is not None else None
maxz = int(maxz) if maxz is not None else None
except Exception:
minz = maxz = None
if args.minzoom is not None: minz = int(args.minzoom)
if args.maxzoom is not None: maxz = int(args.maxzoom)
# --- 3) center: lon/lat dal bbox; zoom = argomento > zoom esistente > media min/max > 0
center_existing = fetch_one(cur, "SELECT value FROM metadata WHERE name='center';")
center_z_existing = None
if center_existing:
parts = [p.strip() for p in center_existing.split(",")]
if len(parts) >= 3 and parts[2] != "":
center_z_existing = parts[2]
if args.zoom is not None:
z = int(args.zoom)
elif center_z_existing is not None:
try:
z = int(center_z_existing)
except Exception:
z = None
else:
z = int((minz + maxz) // 2) if (minz is not None and maxz is not None) else 0
ctr_lon = (W + E) / 2.0
ctr_lat = (S + N) / 2.0
center = f"{ctr_lon:.6f},{ctr_lat:.6f},{z}"
# --- 4) scrittura: crea metadata se manca, REPLACE dei valori
cur.execute("CREATE TABLE IF NOT EXISTS metadata (name TEXT PRIMARY KEY, value TEXT)")
cur.execute("REPLACE INTO metadata(name,value) VALUES(?,?)", ("bounds", bounds))
cur.execute("REPLACE INTO metadata(name,value) VALUES(?,?)", ("center", center))
if minz is not None:
cur.execute("REPLACE INTO metadata(name,value) VALUES(?,?)", ("minzoom", str(minz)))
if maxz is not None:
cur.execute("REPLACE INTO metadata(name,value) VALUES(?,?)", ("maxzoom", str(maxz)))
conn.commit()
# --- 5) stampa riepilogo
print("OK")
for k in ("bounds","center","minzoom","maxzoom"):
v = fetch_one(cur, "SELECT value FROM metadata WHERE name=?", (k,))
print(f"{k} = {v}")
zminmax = cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;").fetchone()
if zminmax and all(v is not None for v in zminmax):
print(f"zoom (tiles) = {zminmax[0]}..{zminmax[1]}")
conn.close()
if __name__ == "__main__":
main()

86
e.py Normal file
View file

@ -0,0 +1,86 @@
import boto3
import os
import logging
from botocore import UNSIGNED
from botocore.config import Config
BUCKET = "copernicus-dem-30m"
OUTPUT_DIR = "./glo30_europe"
LOG_FILE = "download_europe.log"
# Europa approx:
lat_range = range(34, 73) # N34–N72
lon_east = range(0, 46) # E000–E045
lon_west = range(1, 26) # W001–W025
logging.basicConfig(
filename=LOG_FILE,
level=logging.INFO,
format="%(asctime)s [%(levelname)s] %(message)s",
)
console = logging.getLogger("console")
console.setLevel(logging.INFO)
console.addHandler(logging.StreamHandler())
def log(msg):
logging.info(msg)
console.info(msg)
os.makedirs(OUTPUT_DIR, exist_ok=True)
s3 = boto3.client("s3", config=Config(signature_version=UNSIGNED))
def folder_name(lat, lon, east=True):
if east:
return f"Copernicus_DSM_COG_10_N{lat:02d}_00_E{lon:03d}_00_DEM/"
else:
return f"Copernicus_DSM_COG_10_N{lat:02d}_00_W{lon:03d}_00_DEM/"
log("Inizio scansione directory europee...")
found_folders = []
paginator = s3.get_paginator("list_objects_v2")
for page in paginator.paginate(Bucket=BUCKET, Delimiter="/"):
for prefix in page.get("CommonPrefixes", []):
folder = prefix["Prefix"]
# EASTERN EUROPE
for lat in lat_range:
for lon in lon_east:
if folder == folder_name(lat, lon, east=True):
found_folders.append(folder)
log(f"[FOUND] {folder}")
# WESTERN EUROPE (longitudes W)
for lat in lat_range:
for lon in lon_west:
if folder == folder_name(lat, lon, east=False):
found_folders.append(folder)
log(f"[FOUND] {folder}")
log(f"Trovate {len(found_folders)} cartelle europee")
# -----------------------------
# DOWNLOAD DEM PRINCIPALE
# -----------------------------
for folder in found_folders:
tif_name = folder[:-1] + ".tif"
key = folder + tif_name.split("/")[-1]
out_path = os.path.join(OUTPUT_DIR, tif_name.split("/")[-1])
if os.path.exists(out_path):
log(f"[SKIP] {out_path} già presente")
continue
log(f"[DOWNLOAD] {key}")
try:
s3.download_file(BUCKET, key, out_path)
log(f"[OK] {key} scaricato")
except Exception as e:
log(f"[ERROR] {key} → {e}")
log("Completato.")

74
i.py Normal file
View file

@ -0,0 +1,74 @@
import boto3
import os
import logging
from botocore import UNSIGNED
from botocore.config import Config
BUCKET = "copernicus-dem-30m"
OUTPUT_DIR = "./glo30_italia"
LOG_FILE = "download.log"
lat_range = range(36, 48) # N36–N47
lon_range = range(6, 19) # E006–E018
logging.basicConfig(
filename=LOG_FILE,
level=logging.INFO,
format="%(asctime)s [%(levelname)s] %(message)s",
)
console = logging.getLogger("console")
console.setLevel(logging.INFO)
console.addHandler(logging.StreamHandler())
def log(msg):
logging.info(msg)
console.info(msg)
os.makedirs(OUTPUT_DIR, exist_ok=True)
s3 = boto3.client("s3", config=Config(signature_version=UNSIGNED))
def folder_name(lat, lon):
return f"Copernicus_DSM_COG_10_N{lat:02d}_00_E{lon:03d}_00_DEM/"
log("Inizio scansione delle directory italiane...")
found_folders = []
paginator = s3.get_paginator("list_objects_v2")
for page in paginator.paginate(Bucket=BUCKET, Delimiter="/"):
for prefix in page.get("CommonPrefixes", []):
folder = prefix["Prefix"]
# Controlla se la cartella è italiana
for lat in lat_range:
for lon in lon_range:
if folder == folder_name(lat, lon):
found_folders.append(folder)
log(f"[FOUND] {folder}")
log(f"Trovate {len(found_folders)} cartelle italiane")
# -----------------------------
# DOWNLOAD DEM PRINCIPALE
# -----------------------------
for folder in found_folders:
tif_name = folder[:-1] + ".tif" # rimuove "/" e aggiunge .tif
key = folder + tif_name.split("/")[-1]
out_path = os.path.join(OUTPUT_DIR, tif_name.split("/")[-1])
if os.path.exists(out_path):
log(f"[SKIP] {out_path} già presente")
continue
log(f"[DOWNLOAD] {key}")
try:
s3.download_file(BUCKET, key, out_path)
log(f"[OK] {key} scaricato")
except Exception as e:
log(f"[ERROR] {key} → {e}")
log("Completato.")

96
make_tiles_rio_tiler.py Normal file
View file

@ -0,0 +1,96 @@
#!/usr/bin/env python3
import os
import sys
import mercantile
import numpy as np
from rio_tiler.io import COGReader
from rio_tiler.utils import render
# ---------------------------------------------------------
# USO:
# python3 make_tiles_rio_tiler.py input_rgb.tif tiles_out/
#
# Versione finale:
# - BUFFER per evitare buchi tra quadranti
# - PADDING per generare tile anche se il TIFF copre solo parte della tile
# - Logging avanzato
# - Controllo tile vuote
# ---------------------------------------------------------
ZOOM_MIN = 5
ZOOM_MAX = 14
BUFFER = 0.005 # espansione bounds per evitare tile mancanti
TILESIZE = 256 # tile standard
PADDING = 1 # genera tile anche se parziali
def generate_tiles(input_tif, output_dir):
os.makedirs(output_dir, exist_ok=True)
print(f"\n=== Generazione tile per: {input_tif} ===")
print(f"Zoom: {ZOOM_MIN}..{ZOOM_MAX}")
print(f"Buffer: {BUFFER}°")
print(f"Padding: {PADDING} pixel\n")
with COGReader(input_tif) as cog:
for z in range(ZOOM_MIN, ZOOM_MAX + 1):
print(f"\n== Zoom {z} ==")
# Bounds originali del TIFF
bounds = cog.bounds
# Espansione per evitare buchi tra quadranti
minx = bounds[0] - BUFFER
miny = bounds[1] - BUFFER
maxx = bounds[2] + BUFFER
maxy = bounds[3] + BUFFER
# Tile da generare
tiles = list(mercantile.tiles(minx, miny, maxx, maxy, z))
print(f"Tile da generare: {len(tiles)}")
for t in tiles:
out_dir_zxy = os.path.join(output_dir, str(z), str(t.x))
os.makedirs(out_dir_zxy, exist_ok=True)
out_path = os.path.join(out_dir_zxy, f"{t.y}.png")
try:
tile, mask = cog.tile(
t.x, t.y, z,
tilesize=TILESIZE,
padding=PADDING,
resampling_method="nearest"
)
except Exception as e:
print(f" ⚠️ Errore tile {z}/{t.x}/{t.y}: {e}")
continue
# Controllo tile vuote (tutto zero o tutto trasparente)
if mask is not None and np.all(mask == 0):
print(f" ⚠️ Tile vuota saltata: {z}/{t.x}/{t.y}")
continue
img = render(tile, mask=mask, img_format="PNG")
with open(out_path, "wb") as f:
f.write(img)
print(f" ✓ {z}/{t.x}/{t.y}")
print("\n=== Completato senza errori ===\n")
if __name__ == "__main__":
if len(sys.argv) != 3:
print("Uso: python3 make_tiles_rio_tiler.py <input_rgb.tif> <output_dir>")
sys.exit(1)
INPUT_TIF = sys.argv[1]
OUTPUT_DIR = sys.argv[2]
if not os.path.isfile(INPUT_TIF):
print(f"Errore: file non trovato: {INPUT_TIF}")
sys.exit(1)
generate_tiles(INPUT_TIF, OUTPUT_DIR)

58
make_tiles_rio_tiler.py.old Normal file
View file

@ -0,0 +1,58 @@
#!/usr/bin/env python3
import os
import sys
import mercantile
from rio_tiler.io import COGReader
from rio_tiler.utils import render
# ---------------------------------------------------------
# USO:
# python3 make_tiles_rio_tiler.py input_rgb.tif tiles_out/
#
# Genera tile PNG Terrain-RGB senza alterare i pixel.
# ---------------------------------------------------------
ZOOM_MIN = 5
ZOOM_MAX = 14
def generate_tiles(input_tif, output_dir):
os.makedirs(output_dir, exist_ok=True)
with COGReader(input_tif) as cog:
for z in range(ZOOM_MIN, ZOOM_MAX + 1):
print(f"== Zoom {z} ==")
bounds = cog.bounds
tiles = list(mercantile.tiles(bounds[0], bounds[1], bounds[2], bounds[3], z))
for t in tiles:
out_dir_zxy = os.path.join(output_dir, str(z), str(t.x))
os.makedirs(out_dir_zxy, exist_ok=True)
out_path = os.path.join(out_dir_zxy, f"{t.y}.png")
tile, mask = cog.tile(t.x, t.y, z, resampling_method="nearest")
img = render(tile, mask=mask, img_format="PNG")
with open(out_path, "wb") as f:
f.write(img)
print(f"Tile {z}/{t.x}/{t.y} generata")
print("Completato.")
if __name__ == "__main__":
if len(sys.argv) != 3:
print("Uso: python3 make_tiles_rio_tiler.py <input_rgb.tif> <output_dir>")
sys.exit(1)
INPUT_TIF = sys.argv[1]
OUTPUT_DIR = sys.argv[2]
if not os.path.isfile(INPUT_TIF):
print(f"Errore: file non trovato: {INPUT_TIF}")
sys.exit(1)
generate_tiles(INPUT_TIF, OUTPUT_DIR)

71
make_tiles_rio_tiler.py.old2 Normal file
View file

@ -0,0 +1,71 @@
#!/usr/bin/env python3
import os
import sys
import mercantile
from rio_tiler.io import COGReader
from rio_tiler.utils import render
# ---------------------------------------------------------
# USO:
# python3 make_tiles_rio_tiler.py input_rgb.tif tiles_out/
#
# Genera tile PNG Terrain-RGB senza alterare i pixel.
# Corretto per evitare buchi tra quadranti Copernicus.
# ---------------------------------------------------------
ZOOM_MIN = 5
ZOOM_MAX = 14
BUFFER = 0.001 # espansione per evitare tile mancanti ai bordi
def generate_tiles(input_tif, output_dir):
os.makedirs(output_dir, exist_ok=True)
with COGReader(input_tif) as cog:
for z in range(ZOOM_MIN, ZOOM_MAX + 1):
print(f"== Zoom {z} ==")
# Bounds originali del TIFF
bounds = cog.bounds
# Espansione per evitare buchi tra quadranti
minx = bounds[0] - BUFFER
miny = bounds[1] - BUFFER
maxx = bounds[2] + BUFFER
maxy = bounds[3] + BUFFER
# Tile da generare
tiles = list(mercantile.tiles(minx, miny, maxx, maxy, z))
for t in tiles:
out_dir_zxy = os.path.join(output_dir, str(z), str(t.x))
os.makedirs(out_dir_zxy, exist_ok=True)
out_path = os.path.join(out_dir_zxy, f"{t.y}.png")
# Estrazione tile
tile, mask = cog.tile(t.x, t.y, z, resampling_method="nearest")
# Render PNG
img = render(tile, mask=mask, img_format="PNG")
with open(out_path, "wb") as f:
f.write(img)
print(f"Tile {z}/{t.x}/{t.y} generata")
print("Completato.")
if __name__ == "__main__":
if len(sys.argv) != 3:
print("Uso: python3 make_tiles_rio_tiler.py <input_rgb.tif> <output_dir>")
sys.exit(1)
INPUT_TIF = sys.argv[1]
OUTPUT_DIR = sys.argv[2]
if not os.path.isfile(INPUT_TIF):
print(f"Errore: file non trovato: {INPUT_TIF}")
sys.exit(1)
generate_tiles(INPUT_TIF, OUTPUT_DIR)

139
mbtiles_set_metadata.sh Executable file
View file

@ -0,0 +1,139 @@
#!/usr/bin/env bash
set -euo pipefail
# --------------------------------------------------------------------
# Imposta metadati essenziali in un file .mbtiles (schema XYZ)
# - Calcola bounds (ovest,sud,est,nord) dal raster sorgente (GeoTIFF/COG)
# - Calcola center (lon,lat,zoom)
# - Legge min/max zoom dalla tabella 'tiles' e li scrive nella 'metadata'
# - Imposta name/format/type/attribution
#
# Requisiti: gdalinfo (GDAL), jq, sqlite3
# --------------------------------------------------------------------
# USO:
# ./mbtiles_set_metadata.sh \
# --mbtiles aa.mbtiles \
# --raster aa_rgb.tif \
# [--name "Terrain RGB Adamello"] \
# [--type overlay|baselayer] \
# [--format png|jpg] \
# [--zoom <Zcenter>] \
# [--attribution " Copernicus DEM (GLO-30), ESA"]
#
# NOTE:
# - 'bounds' e 'center' sono calcolati in EPSG:4326 (lon/lat).
# - Se il raster non riassume esattamente lestensione delle tile,
# passa un raster "clip" coerente con larea tilata.
# - Lo script usa REPLACE INTO per compatibilit con SQLite < 3.24.
# --------------------------------------------------------------------
# ---------- parsing argomenti ----------
MBTILES=""
RASTER=""
LAYER_NAME=""
LAYER_TYPE="overlay"
FORMAT="png"
CENTER_Z=""
ATTRIBUTION=""
while [[ $# -gt 0 ]]; do
case "$1" in
--mbtiles) MBTILES="$2"; shift 2;;
--raster) RASTER="$2"; shift 2;;
--name) LAYER_NAME="$2"; shift 2;;
--type) LAYER_TYPE="$2"; shift 2;;
--format) FORMAT="$2"; shift 2;;
--zoom) CENTER_Z="$2"; shift 2;;
--attribution) ATTRIBUTION="$2"; shift 2;;
-h|--help)
sed -n '1,80p' "$0"; exit 0;;
*)
echo "Argomento sconosciuto: $1" >&2; exit 1;;
esac
done
if [[ -z "$MBTILES" || -z "$RASTER" ]]; then
echo "Uso: $0 --mbtiles <file.mbtiles> --raster <raster.tif> [opzioni]" >&2
exit 1
fi
if [[ ! -f "$MBTILES" ]]; then
echo "File MBTiles non trovato: $MBTILES" >&2; exit 1
fi
if [[ ! -f "$RASTER" ]]; then
echo "Raster non trovato: $RASTER" >&2; exit 1
fi
command -v gdalinfo >/dev/null || { echo "gdalinfo non trovato"; exit 1; }
command -v jq >/dev/null || { echo "jq non trovato"; exit 1; }
command -v sqlite3 >/dev/null || { echo "sqlite3 non trovato"; exit 1; }
# ---------- step 1: calcolo bounds dal raster (EPSG:4326) ----------
# Usiamo wgs84Extent se presente, altrimenti cornerCoordinates.
# gdalinfo -json fornisce geometrie e corner coords in lon/lat.
export LC_ALL=C
export LANG=C
GDAL_JSON="$(gdalinfo -json "$RASTER")"
has_wgs84=$(echo "$GDAL_JSON" | jq 'has("wgs84Extent")')
if [[ "$has_wgs84" == "true" ]]; then
# Polygon ring: [ [minLon,minLat], [maxLon,minLat], [maxLon,maxLat], [minLon,maxLat], ... ]
read -r MINLON MINLAT MAXLON MAXLAT <<<"$(echo "$GDAL_JSON" \
| jq -r '.wgs84Extent.coordinates[0] |
[.[0][0], .[0][1], .[2][0], .[2][1]] | @tsv')"
else
# cornerCoordinates: upperLeft, lowerLeft, upperRight, lowerRight (lon,lat)
read -r MINLON MINLAT MAXLON MAXLAT <<<"$(echo "$GDAL_JSON" \
| jq -r '.cornerCoordinates as $c |
[$c.lowerLeft[0], $c.lowerLeft[1], $c.upperRight[0], $c.upperRight[1]] | @tsv')"
fi
# Arrotonda a 6 decimali per compattare la stringa
printf -v BOUNDS "%.6f,%.6f,%.6f,%.6f" "$MINLON" "$MINLAT" "$MAXLON" "$MAXLAT"
# ---------- step 2: calcolo center (lon,lat,zoom) ----------
# lon,lat = met dei bounds; zoom: se non fornito, media tra min e max zoom presenti nelle tile
CTR_LON=$(awk -v a="$MINLON" -v b="$MAXLON" 'BEGIN{printf "%.6f",(a+b)/2.0}')
CTR_LAT=$(awk -v a="$MINLAT" -v b="$MAXLAT" 'BEGIN{printf "%.6f",(a+b)/2.0}')
# Leggi min/max zoom dalla tabella 'tiles'
read -r TMIN TMAX <<<"$(sqlite3 "$MBTILES" "SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")"
if [[ -z "$TMIN" || -z "$TMAX" ]]; then
echo "Attenzione: tabella 'tiles' vuota o mancante in $MBTILES" >&2
TMIN=0; TMAX=0
fi
if [[ -z "${CENTER_Z}" ]]; then
CENTER_Z=$(( (TMIN + TMAX) / 2 ))
fi
CENTER="${CTR_LON},${CTR_LAT},${CENTER_Z}"
# ---------- step 3: definisci name/format/type/attribution ----------
if [[ -z "$LAYER_NAME" ]]; then
# Usa il nome file come default
base="$(basename "$MBTILES")"
LAYER_NAME="${base%.*}"
fi
# ---------- step 4: crea tabella metadata (se non esiste) e scrivi i metadati ----------
sqlite3 "$MBTILES" <<SQL
CREATE TABLE IF NOT EXISTS metadata (name TEXT PRIMARY KEY, value TEXT);
REPLACE INTO metadata(name,value) VALUES('name', '$LAYER_NAME');
REPLACE INTO metadata(name,value) VALUES('format', '$FORMAT');
REPLACE INTO metadata(name,value) VALUES('type', '$LAYER_TYPE');
REPLACE INTO metadata(name,value) VALUES('minzoom', '$TMIN');
REPLACE INTO metadata(name,value) VALUES('maxzoom', '$TMAX');
REPLACE INTO metadata(name,value) VALUES('bounds', '$BOUNDS');
REPLACE INTO metadata(name,value) VALUES('center', '$CENTER');
$( [[ -n "$ATTRIBUTION" ]] && echo "REPLACE INTO metadata(name,value) VALUES('attribution','$ATTRIBUTION');" )
SQL
# ---------- step 5: stampa riepilogo ----------
echo "==> METADATA scritti in: $MBTILES"
sqlite3 "$MBTILES" "SELECT name,value FROM metadata ORDER BY name;"
echo
echo "Zoom effettivi nelle tile: $(sqlite3 "$MBTILES" "SELECT MIN(zoom_level)||'..'||MAX(zoom_level) FROM tiles;")"
echo "Bounds: $BOUNDS"
echo "Center: $CENTER"

19
rgb2alt.py Normal file
View file

@ -0,0 +1,19 @@
#!/usr/bin/env python3
def rgb_to_altitude(R, G, B, scale=0.1, offset=-10000):
value = (R * 65536) + (G * 256) + B
return value * scale + offset
if __name__ == "__main__":
import sys
if len(sys.argv) != 4:
print("Uso: python3 rgb2alt.py R G B")
sys.exit(1)
R = int(sys.argv[1])
G = int(sys.argv[2])
B = int(sys.argv[3])
alt = rgb_to_altitude(R, G, B)
print(alt)

49
rgbify_batch.py Normal file
View file

@ -0,0 +1,49 @@
#!/usr/bin/env python3
import os
import sys
import subprocess
def main():
if len(sys.argv) != 3:
print("Uso: python3 rgbify_batch.py <input_dir> <output_dir>")
sys.exit(1)
input_dir = sys.argv[1]
output_dir = sys.argv[2]
if not os.path.isdir(input_dir):
print(f"Errore: directory di input non trovata: {input_dir}")
sys.exit(1)
os.makedirs(output_dir, exist_ok=True)
# Scansiona tutti i .tif nella directory di input
for filename in os.listdir(input_dir):
if not filename.lower().endswith(".tif"):
continue
in_path = os.path.join(input_dir, filename)
out_name = filename.replace(".tif", "_rgb.tif")
out_path = os.path.join(output_dir, out_name)
print(f"[PROCESS] {filename} → {out_name}")
cmd = [
"rio", "rgbify",
"-b", "-10000",
"-i", "0.1",
in_path,
out_path
]
try:
subprocess.run(cmd, check=True)
print(f"[OK] {out_name}")
except subprocess.CalledProcessError as e:
print(f"[ERROR] {filename}: {e}")
print("Completato.")
if __name__ == "__main__":
main()

130
set_mbtiles_metadata.py Normal file
View file

@ -0,0 +1,130 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
def read_gdal_bounds(raster):
"""Legge i bounds WGS84 dal TIFF usando gdalinfo -json."""
info = subprocess.check_output(["gdalinfo", "-json", raster])
info = json.loads(info)
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
# Normalizzazione latitudine
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
return (minlon, minlat, maxlon, maxlat)
def open_mbtiles(mbtiles):
"""Apre MBTiles e garantisce che la tabella metadata esista."""
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
# Verifica tabella tiles
cur.execute("SELECT COUNT(*) FROM sqlite_master WHERE type='table' AND name='tiles';")
if cur.fetchone()[0] == 0:
raise RuntimeError("La tabella 'tiles' non esiste nell'MBTiles.")
# Conta tile
cur.execute("SELECT COUNT(*) FROM tiles;")
tile_count = cur.fetchone()[0]
if tile_count == 0:
raise RuntimeError("La tabella 'tiles' è vuota.")
# Zoom min/max
cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")
minzoom, maxzoom = cur.fetchone()
# Crea metadata se manca
cur.execute("""
CREATE TABLE IF NOT EXISTS metadata (
name TEXT PRIMARY KEY,
value TEXT
);
""")
# Leggi metadata esistenti
cur.execute("SELECT name, value FROM metadata;")
metadata = dict(cur.fetchall())
return conn, cur, tile_count, minzoom, maxzoom, metadata
def write_metadata(cur, name, value):
cur.execute("REPLACE INTO metadata (name, value) VALUES (?, ?)", (name, value))
def main():
if len(sys.argv) < 3:
print("Uso: python3 set_mbtiles_metadata.py <raster.tif> <file.mbtiles> [name]")
sys.exit(1)
raster = Path(sys.argv[1])
mbtiles = Path(sys.argv[2])
layer_name = sys.argv[3] if len(sys.argv) > 3 else mbtiles.stem
if not raster.exists():
print(f"Raster non trovato: {raster}")
sys.exit(1)
if not mbtiles.exists():
print(f"MBTiles non trovato: {mbtiles}")
sys.exit(1)
print("== Impostazione metadata MBTiles ==")
# 1) Bounds dal TIFF
minlon, minlat, maxlon, maxlat = read_gdal_bounds(raster)
bounds = f"{minlon:.6f},{minlat:.6f},{maxlon:.6f},{maxlat:.6f}"
print(f"Bounds TIFF: {bounds}")
# 2) Lettura MBTiles
conn, cur, tile_count, minzoom, maxzoom, metadata = open_mbtiles(mbtiles)
print(f"Tile count: {tile_count}")
print(f"Zoom effettivi: {minzoom}..{maxzoom}")
# 3) Center
ctr_lon = (minlon + maxlon) / 2
ctr_lat = (minlat + maxlat) / 2
ctr_zoom = (minzoom + maxzoom) // 2
center = f"{ctr_lon:.6f},{ctr_lat:.6f},{ctr_zoom}"
# 4) Scrittura metadata
write_metadata(cur, "name", layer_name)
write_metadata(cur, "format", "png")
write_metadata(cur, "type", "overlay")
write_metadata(cur, "scheme", "xyz")
write_metadata(cur, "bounds", bounds)
write_metadata(cur, "center", center)
write_metadata(cur, "minzoom", str(minzoom))
write_metadata(cur, "maxzoom", str(maxzoom))
conn.commit()
conn.close()
print("\n== Metadata scritti correttamente ==")
print(f"name: {layer_name}")
print(f"format: png")
print(f"type: overlay")
print(f"scheme: xyz")
print(f"bounds: {bounds}")
print(f"center: {center}")
print(f"minzoom: {minzoom}")
print(f"maxzoom: {maxzoom}")
print("\nMBTiles ora è pronto per TileServer GL.")
if __name__ == "__main__":
main()

118
set_mbtiles_metadata_dir.py Normal file
View file

@ -0,0 +1,118 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
def read_bounds(tif):
"""Legge i bounds W,S,E,N dal TIFF usando gdalinfo -json, arrotondati a 6 decimali."""
info = json.loads(subprocess.check_output(["gdalinfo", "-json", str(tif)]))
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
# Normalizzazione
if minlon > maxlon:
minlon, maxlon = maxlon, minlon
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
return (
round(minlon, 6),
round(minlat, 6),
round(maxlon, 6),
round(maxlat, 6),
)
def main():
if len(sys.argv) != 3:
print("Uso: python3 set_mbtiles_metadata_multi.py <dir_tif> <merged.mbtiles>")
sys.exit(1)
tif_dir = Path(sys.argv[1])
mbtiles = Path(sys.argv[2])
if not tif_dir.exists() or not tif_dir.is_dir():
print(f"Directory TIFF non valida: {tif_dir}")
sys.exit(1)
if not mbtiles.exists():
print(f"MBTiles non trovato: {mbtiles}")
sys.exit(1)
# Trova tutti i TIFF nella directory
tifs = sorted([p for p in tif_dir.iterdir() if p.suffix.lower() in (".tif", ".tiff")])
if not tifs:
print("Nessun TIFF trovato nella directory.")
sys.exit(1)
print(f"Trovati {len(tifs)} TIFF:")
for t in tifs:
print(" -", t.name)
# Calcolo bounds unificati
bounds_list = [read_bounds(t) for t in tifs]
W = min(b[0] for b in bounds_list)
S = min(b[1] for b in bounds_list)
E = max(b[2] for b in bounds_list)
N = max(b[3] for b in bounds_list)
bounds = f"{W:.6f},{S:.6f},{E:.6f},{N:.6f}"
print("\nBounds unificati:", bounds)
# Apri MBTiles
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
cur.execute("""
CREATE TABLE IF NOT EXISTS metadata (
name TEXT PRIMARY KEY,
value TEXT
);
""")
# Zoom reali dal MBTiles
cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")
minzoom, maxzoom = cur.fetchone()
print(f"Zoom reali: {minzoom}..{maxzoom}")
# Center
centerLon = (W + E) / 2
centerLat = (S + N) / 2
centerZoom = (minzoom + maxzoom) // 2
center = f"{centerLon:.6f},{centerLat:.6f},{centerZoom}"
# Scrittura metadata
cur.execute("REPLACE INTO metadata VALUES ('bounds', ?)", (bounds,))
cur.execute("REPLACE INTO metadata VALUES ('center', ?)", (center,))
cur.execute("REPLACE INTO metadata VALUES ('minzoom', ?)", (str(minzoom),))
cur.execute("REPLACE INTO metadata VALUES ('maxzoom', ?)", (str(maxzoom),))
cur.execute("REPLACE INTO metadata VALUES ('format', 'png')")
cur.execute("REPLACE INTO metadata VALUES ('type', 'overlay')")
cur.execute("REPLACE INTO metadata VALUES ('scheme', 'xyz')")
conn.commit()
conn.close()
print("\nMetadata scritti correttamente:")
print("bounds:", bounds)
print("center:", center)
print("minzoom:", minzoom)
print("maxzoom:", maxzoom)
print("format: png")
print("type: overlay")
print("scheme: xyz")
if __name__ == "__main__":
main()

85
set_mbtiles_metadata_multi.py Normal file
View file

@ -0,0 +1,85 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
def read_bounds(tif):
info = json.loads(subprocess.check_output(["gdalinfo", "-json", tif]))
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
if minlon > maxlon:
minlon, maxlon = maxlon, minlon
return (
round(minlon, 6),
round(minlat, 6),
round(maxlon, 6),
round(maxlat, 6),
)
def main():
if len(sys.argv) != 4:
print("Uso: python3 set_mbtiles_metadata_multi.py raster1.tif raster2.tif merged.mbtiles")
sys.exit(1)
tif1 = Path(sys.argv[1])
tif2 = Path(sys.argv[2])
mbtiles = Path(sys.argv[3])
b1 = read_bounds(tif1)
b2 = read_bounds(tif2)
W = min(b1[0], b2[0])
S = min(b1[1], b2[1])
E = max(b1[2], b2[2])
N = max(b1[3], b2[3])
bounds = f"{W:.6f},{S:.6f},{E:.6f},{N:.6f}"
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
cur.execute("CREATE TABLE IF NOT EXISTS metadata (name TEXT PRIMARY KEY, value TEXT);")
cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")
minzoom, maxzoom = cur.fetchone()
centerLon = (W + E) / 2
centerLat = (S + N) / 2
centerZoom = (minzoom + maxzoom) // 2
center = f"{centerLon:.6f},{centerLat:.6f},{centerZoom}"
cur.execute("REPLACE INTO metadata VALUES ('bounds', ?)", (bounds,))
cur.execute("REPLACE INTO metadata VALUES ('center', ?)", (center,))
cur.execute("REPLACE INTO metadata VALUES ('minzoom', ?)", (str(minzoom),))
cur.execute("REPLACE INTO metadata VALUES ('maxzoom', ?)", (str(maxzoom),))
cur.execute("REPLACE INTO metadata VALUES ('format', 'png')")
cur.execute("REPLACE INTO metadata VALUES ('type', 'overlay')")
cur.execute("REPLACE INTO metadata VALUES ('scheme', 'xyz')")
conn.commit()
conn.close()
print("Metadata scritti correttamente:")
print("bounds:", bounds)
print("center:", center)
print("minzoom:", minzoom)
print("maxzoom:", maxzoom)
print("format: png")
print("type: overlay")
print("scheme: xyz")
if __name__ == "__main__":
main()

144
verify_and_fix_mbtiles.py Normal file
View file

@ -0,0 +1,144 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
def read_gdal_bounds(raster):
"""Legge i bounds WGS84 dal TIFF usando gdalinfo -json."""
info = subprocess.check_output(["gdalinfo", "-json", raster])
info = json.loads(info)
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
# Normalizzazione latitudine
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
return (minlon, minlat, maxlon, maxlat)
def read_mbtiles(mbtiles):
"""Legge tile info e metadata da un MBTiles."""
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
# Verifica tabella tiles
cur.execute("SELECT COUNT(*) FROM sqlite_master WHERE type='table' AND name='tiles';")
if cur.fetchone()[0] == 0:
raise RuntimeError("La tabella 'tiles' non esiste nell'MBTiles.")
# Conta tile
cur.execute("SELECT COUNT(*) FROM tiles;")
tile_count = cur.fetchone()[0]
if tile_count == 0:
raise RuntimeError("La tabella 'tiles' è vuota.")
# Zoom min/max
cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")
minzoom, maxzoom = cur.fetchone()
# Metadata
cur.execute("CREATE TABLE IF NOT EXISTS metadata (name TEXT PRIMARY KEY, value TEXT);")
cur.execute("SELECT name, value FROM metadata;")
metadata = dict(cur.fetchall())
return conn, cur, tile_count, minzoom, maxzoom, metadata
def write_metadata(cur, name, value):
cur.execute("REPLACE INTO metadata (name, value) VALUES (?, ?)", (name, value))
def main():
if len(sys.argv) != 3:
print("Uso: python3 verify_and_fix_mbtiles.py <raster.tif> <file.mbtiles>")
sys.exit(1)
raster = Path(sys.argv[1])
mbtiles = Path(sys.argv[2])
if not raster.exists():
print(f"Raster non trovato: {raster}")
sys.exit(1)
if not mbtiles.exists():
print(f"MBTiles non trovato: {mbtiles}")
sys.exit(1)
print("== Verifica e correzione metadata MBTiles ==")
# 1) Bounds dal TIFF
minlon, minlat, maxlon, maxlat = read_gdal_bounds(raster)
expected_bounds = f"{minlon:.6f},{minlat:.6f},{maxlon:.6f},{maxlat:.6f}"
print(f"Bounds TIFF: {expected_bounds}")
# 2) Lettura MBTiles
conn, cur, tile_count, minzoom, maxzoom, metadata = read_mbtiles(mbtiles)
print(f"Tile count: {tile_count}")
print(f"Zoom effettivi: {minzoom}..{maxzoom}")
# 3) Verifica/correzione bounds
mb_bounds = metadata.get("bounds")
if mb_bounds != expected_bounds:
print(f"✘ Bounds NON corretti → FIX")
write_metadata(cur, "bounds", expected_bounds)
else:
print("✔ Bounds OK")
# 4) Verifica/correzione minzoom/maxzoom
mz = int(metadata.get("minzoom", -1))
xz = int(metadata.get("maxzoom", -1))
if mz != minzoom or xz != maxzoom:
print(f"✘ minzoom/maxzoom NON coerenti → FIX")
write_metadata(cur, "minzoom", str(minzoom))
write_metadata(cur, "maxzoom", str(maxzoom))
else:
print("✔ minzoom/maxzoom OK")
# 5) Verifica/correzione center
ctr_lon = (minlon + maxlon) / 2
ctr_lat = (minlat + maxlat) / 2
ctr_zoom = (minzoom + maxzoom) // 2
expected_center = f"{ctr_lon:.6f},{ctr_lat:.6f},{ctr_zoom}"
if metadata.get("center") != expected_center:
print(f"✘ Center NON corretto → FIX")
write_metadata(cur, "center", expected_center)
else:
print("✔ Center OK")
# 6) Verifica/correzione format
if metadata.get("format") != "png":
print("✘ Format NON corretto → FIX")
write_metadata(cur, "format", "png")
else:
print("✔ Format OK")
# 7) Verifica/correzione scheme
if metadata.get("scheme") != "xyz":
print("✘ Scheme NON corretto → FIX")
write_metadata(cur, "scheme", "xyz")
else:
print("✔ Scheme OK")
# 8) Commit
conn.commit()
conn.close()
print("\n== Correzione completata ==")
print("MBTiles ora ha metadata coerenti e validi.")
if __name__ == "__main__":
main()

122
verify_bounds_tileservergl.py Normal file
View file

@ -0,0 +1,122 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
EPS = 1e-6 # tolleranza per confronti float
def almost_equal(a, b):
return abs(a - b) < EPS
def read_gdal_bounds(raster):
"""Legge i bounds W,S,E,N dal TIFF usando gdalinfo -json, arrotondati a 6 decimali."""
info = subprocess.check_output(["gdalinfo", "-json", raster])
info = json.loads(info)
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
# Normalizzazione
if minlon > maxlon:
minlon, maxlon = maxlon, minlon
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
# Arrotondamento a 6 decimali (standard MBTiles)
return (
round(minlon, 6),
round(minlat, 6),
round(maxlon, 6),
round(maxlat, 6),
)
def read_mbtiles_bounds(cur):
cur.execute("SELECT value FROM metadata WHERE name='bounds';")
row = cur.fetchone()
if not row:
return None
try:
w, s, e, n = map(float, row[0].split(","))
return (w, s, e, n)
except:
return None
def write_bounds(cur, bounds):
w, s, e, n = bounds
value = f"{w:.6f},{s:.6f},{e:.6f},{n:.6f}"
cur.execute("REPLACE INTO metadata (name,value) VALUES ('bounds', ?)", (value,))
return value
def main():
if len(sys.argv) != 3:
print("Uso: python3 verify_bounds_tileservergl.py <raster.tif> <file.mbtiles>")
sys.exit(1)
raster = Path(sys.argv[1])
mbtiles = Path(sys.argv[2])
if not raster.exists():
print(f"Raster non trovato: {raster}")
sys.exit(1)
if not mbtiles.exists():
print(f"MBTiles non trovato: {mbtiles}")
sys.exit(1)
print("== Verifica bounds per TileServer GL (W,S,E,N) ==")
# Bounds dal TIFF
tiff_bounds = read_gdal_bounds(raster)
print(f"Bounds TIFF (W,S,E,N): {tiff_bounds}")
# Apri MBTiles
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
# Assicura tabella metadata
cur.execute("""
CREATE TABLE IF NOT EXISTS metadata (
name TEXT PRIMARY KEY,
value TEXT
);
""")
# Bounds MBTiles
mb_bounds = read_mbtiles_bounds(cur)
print(f"Bounds MBTiles: {mb_bounds}")
# Se mancano → scrivili
if mb_bounds is None:
print("✘ Bounds mancanti → FIX")
new_value = write_bounds(cur, tiff_bounds)
print(f"✔ Bounds scritti: {new_value}")
conn.commit()
conn.close()
print("\n== Verifica completata ==")
return
# Confronto con tolleranza
w1, s1, e1, n1 = tiff_bounds
w2, s2, e2, n2 = mb_bounds
if not (almost_equal(w1, w2) and almost_equal(s1, s2) and almost_equal(e1, e2) and almost_equal(n1, n2)):
print("✘ Bounds NON corretti → FIX")
new_value = write_bounds(cur, tiff_bounds)
print(f"✔ Bounds aggiornati a: {new_value}")
conn.commit()
else:
print("✔ Bounds OK (nessuna correzione necessaria)")
conn.close()
print("\n== Verifica completata ==")
if __name__ == "__main__":
main()

122
verify_mbtiles_metadata.py Normal file
View file

@ -0,0 +1,122 @@
#!/usr/bin/env python3
import sqlite3
import json
import subprocess
import sys
from pathlib import Path
def read_gdal_bounds(raster):
"""Legge i bounds WGS84 dal TIFF usando gdalinfo -json."""
info = subprocess.check_output(["gdalinfo", "-json", raster])
info = json.loads(info)
if "wgs84Extent" in info:
coords = info["wgs84Extent"]["coordinates"][0]
minlon, minlat = coords[0]
maxlon, maxlat = coords[2]
else:
cc = info["cornerCoordinates"]
minlon, minlat = cc["lowerLeft"]
maxlon, maxlat = cc["upperRight"]
# Normalizzazione latitudine
if minlat > maxlat:
minlat, maxlat = maxlat, minlat
return (minlon, minlat, maxlon, maxlat)
def read_mbtiles_metadata(mbtiles):
"""Legge metadata e tile info da un MBTiles."""
conn = sqlite3.connect(mbtiles)
cur = conn.cursor()
# Verifica tabella tiles
cur.execute("SELECT COUNT(*) FROM sqlite_master WHERE type='table' AND name='tiles';")
if cur.fetchone()[0] == 0:
raise RuntimeError("La tabella 'tiles' non esiste nell'MBTiles.")
# Conta tile
cur.execute("SELECT COUNT(*) FROM tiles;")
tile_count = cur.fetchone()[0]
if tile_count == 0:
raise RuntimeError("La tabella 'tiles' è vuota.")
# Zoom min/max
cur.execute("SELECT MIN(zoom_level), MAX(zoom_level) FROM tiles;")
minzoom, maxzoom = cur.fetchone()
# Metadata
cur.execute("SELECT name, value FROM metadata;")
metadata = dict(cur.fetchall())
conn.close()
return tile_count, minzoom, maxzoom, metadata
def main():
if len(sys.argv) != 3:
print("Uso: python3 verify_mbtiles_metadata.py <raster.tif> <file.mbtiles>")
sys.exit(1)
raster = Path(sys.argv[1])
mbtiles = Path(sys.argv[2])
if not raster.exists():
print(f"Raster non trovato: {raster}")
sys.exit(1)
if not mbtiles.exists():
print(f"MBTiles non trovato: {mbtiles}")
sys.exit(1)
print("== Verifica metadata MBTiles ==")
# 1) Bounds dal TIFF
minlon, minlat, maxlon, maxlat = read_gdal_bounds(raster)
print(f"Bounds TIFF: {minlon},{minlat},{maxlon},{maxlat}")
# 2) Lettura MBTiles
tile_count, minzoom, maxzoom, metadata = read_mbtiles_metadata(mbtiles)
print(f"Tile count: {tile_count}")
print(f"Zoom effettivi: {minzoom}..{maxzoom}")
# 3) Confronto bounds
mb_bounds = metadata.get("bounds", None)
if mb_bounds:
mb_minlon, mb_minlat, mb_maxlon, mb_maxlat = map(float, mb_bounds.split(","))
print(f"Bounds MBTiles: {mb_bounds}")
if abs(mb_minlon - minlon) < 1e-4 and \
abs(mb_minlat - minlat) < 1e-4 and \
abs(mb_maxlon - maxlon) < 1e-4 and \
abs(mb_maxlat - maxlat) < 1e-4:
print("✔ Bounds OK")
else:
print("✘ Bounds NON corrispondono al TIFF")
else:
print("✘ Metadata 'bounds' mancante")
# 4) Confronto minzoom/maxzoom
if "minzoom" in metadata and "maxzoom" in metadata:
mz = int(metadata["minzoom"])
xz = int(metadata["maxzoom"])
if mz == minzoom and xz == maxzoom:
print("✔ minzoom/maxzoom OK")
else:
print(f"✘ minzoom/maxzoom NON coerenti: metadata={mz}..{xz}, reali={minzoom}..{maxzoom}")
else:
print("✘ Metadata minzoom/maxzoom mancanti")
# 5) Center
if "center" in metadata:
print(f"Center: {metadata['center']}")
else:
print("✘ Metadata 'center' mancante")
print("\n== Verifica completata ==")
if __name__ == "__main__":
main()

68
w.py Normal file
View file

@ -0,0 +1,68 @@
import boto3
import os
import logging
from botocore import UNSIGNED
from botocore.config import Config
BUCKET = "copernicus-dem-30m"
OUTPUT_DIR = "./glo30_world"
LOG_FILE = "download_world.log"
logging.basicConfig(
filename=LOG_FILE,
level=logging.INFO,
format="%(asctime)s [%(levelname)s] %(message)s",
)
console = logging.getLogger("console")
console.setLevel(logging.INFO)
console.addHandler(logging.StreamHandler())
def log(msg):
logging.info(msg)
console.info(msg)
os.makedirs(OUTPUT_DIR, exist_ok=True)
s3 = boto3.client("s3", config=Config(signature_version=UNSIGNED))
log("Inizio scansione di tutte le cartelle globali...")
folders = []
paginator = s3.get_paginator("list_objects_v2")
# Legge SOLO le directory top-level
for page in paginator.paginate(Bucket=BUCKET, Delimiter="/"):
for prefix in page.get("CommonPrefixes", []):
folder = prefix["Prefix"]
# Filtra solo i tile DEM
if folder.startswith("Copernicus_DSM_COG_10_") and
folder.endswith("_DEM/"):
folders.append(folder)
log(f"[FOUND] {folder}")
log(f"Trovate {len(folders)} cartelle globali")
# -----------------------------
# DOWNLOAD SOLO IL DEM PRINCIPALE
# -----------------------------
for folder in folders:
tif_name = folder[:-1] + ".tif" # rimuove "/" e aggiunge .tif
key = folder + tif_name.split("/")[-1]
out_path = os.path.join(OUTPUT_DIR, tif_name.split("/")[-1])
if os.path.exists(out_path):
log(f"[SKIP] {out_path} già presente")
continue
log(f"[DOWNLOAD] {key}")
try:
s3.download_file(BUCKET, key, out_path)
log(f"[OK] {key} scaricato")
except Exception as e:
log(f"[ERROR] {key} → {e}")
log("Completato.")