Rischio Subsidenza Italia: Come Abbiamo Costruito il Dataset
Se lavori con dati territoriali in Italia — protezione civile, urbanistica, assicurazioni, infrastrutture — probabilmente hai incrociato il problema della subsidenza del suolo almeno una volta. Il terreno si muove, lentamente ma inesorabilmente, e sapere dove e quanto è un'informazione che vale oro per qualsiasi analisi di rischio strutturale.
Il problema è che i dati esistono, ma non sono facili da usare. Il servizio EGMS (European Ground Motion Service) di Copernicus pubblica misurazioni InSAR derivate da Sentinel-1 con precisione millimetrica, ma il formato è pensato per chi fa ricerca, non per chi deve fare una query spaziale su un database.
Abbiamo preso quei dati e li abbiamo trasformati in qualcosa di immediatamente operativo: un layer poligonale a 100 m di risoluzione con classi di rischio precompilate, scaricabile in tre formati e pronto per l'importazione in PostGIS, QGIS o qualsiasi stack GIS.
Cosa c'è dentro il dataset
Il dataset copre l'intera Italia con una griglia regolare di celle 100×100 m in proiezione EPSG:3035. Ogni cella corrisponde a un punto di misura EGMS L3 Ortho Verticale e contiene:
| Campo | Tipo | Descrizione |
|---|---|---|
mean_velocity |
float | Velocità verticale media 2019–2023 (mm/anno) |
acceleration |
float | Accelerazione del movimento (mm/anno²) |
risk_class |
int | Classe di rischio da 1 a 5 |
risk_class_label |
string | Etichetta: trascurabile, basso, moderato, alto, molto_alto |
risk_index |
float | Indice di rischio continuo 0–1 |
direction |
string | Direzione del moto: subsidenza, sollevamento, stabile |
velocity_class |
int | Classe di velocità da 1 a 5 |
trend |
string | Tendenza: accelerazione, decelerazione, stabile |
signal_quality |
string | Qualità del segnale: buona, moderata, scarsa |
Ogni poligono è un quadrato perfetto di 100 m di lato. Non c'è aggregazione, non c'è interpolazione: è una mappatura 1:1 dal punto EGMS alla cella che lo contiene.
I numeri
- 8.486.494 celle poligonali
- Risoluzione: 100 m (EPSG:3035, griglia regolare)
- Copertura temporale: 2019–2023
- Range velocità: da −59,8 a +130,7 mm/anno
- 5 classi di rischio con soglie a 2, 5, 10 e 20 mm/anno
La pipeline: dal download alla cella
Il processo è interamente automatizzato e riproducibile. Ecco cosa succede, passo per passo.
1. Download delle tile EGMS
Il servizio EGMS pubblica i dati in tile da 100×100 km. Per l'Italia servono 118 tile della componente verticale (Ortho U), periodo 2019–2023. Il download passa dall'API ufficiale di Copernicus con autenticazione via token.
Non è banale come sembra. L'API ha rate limiting aggressivo (errori 429), alcune tile esistono in due versioni e 4 tile su 118 sono risultate corrotte al primo tentativo. Lo script gestisce tutto: retry con backoff esponenziale, validazione ZIP post-download, re-download automatico dei file corrotti.
Tempo totale di download: circa 8 ore per ~50 GB di dati compressi.
2. Lettura e pulizia
Ogni tile è un file ZIP che contiene un CSV con le coordinate in EPSG:3035 e la serie temporale completa delle misurazioni. Dalla serie, estraiamo le colonne rilevanti:
mean_velocity: velocità verticale media sull'intero periodomean_velocity_std: deviazione standard della velocitàacceleration: trend di accelerazione/decelerazionermse: errore quadratico medio del fitseasonality: componente stagionale
I duplicati (punti nelle zone di sovrapposizione tra tile) vengono rimossi tenendo la prima occorrenza. I punti con velocità nulla o mancante vengono scartati.
3. Trasformazione di coordinate e clipping
Le coordinate native EGMS sono in EPSG:3035 (LAEA Europe). Per il clipping spaziale, trasformiamo temporaneamente in WGS84 (EPSG:4326) e applichiamo un buffer di 50 km attorno al confine italiano, scaricato da Natural Earth 10m.
Questo rimuove circa il 20% dei punti originali — le tile EGMS coprono un'area più ampia dell'Italia, includendo pezzi di Francia, Austria, Slovenia e mare.
4. Classificazione del rischio
La classificazione combina velocità e accelerazione:
Classi di velocità (basate su |velocità|):
| Classe | Etichetta | Soglia |
|---|---|---|
| 1 | trascurabile | < 2 mm/a |
| 2 | basso | 2–5 mm/a |
| 3 | moderato | 5–10 mm/a |
| 4 | alto | 10–20 mm/a |
| 5 | molto_alto | ≥ 20 mm/a |
Bump di accelerazione: se l'accelerazione è > 1 mm/a² e ha lo stesso segno della velocità (il movimento sta peggiorando), la classe di rischio viene incrementata di 1.
Indice di rischio continuo: R = min(|vel| / 20, 1) × (1 + modificatore), dove il modificatore è +0,2 per accelerazione concordante e −0,1 per decelerazione.
5. Da punto a poligono
Qui sta il passaggio chiave. I punti EGMS L3 siedono su una griglia perfettamente regolare: ogni coordinata è un multiplo di 100 + 50 (es. 4.000.350, 2.500.150). Ogni punto rappresenta il centroide della cella 100×100 m che lo contiene.
Lo script verifica questa regolarità (assertion sul modulo 100) e poi genera i poligoni con shapely.box(e-50, n-50, e+50, n+50). Nessuna interpolazione, nessun Voronoi, nessun buffer arbitrario: solo la geometria nativa della griglia EGMS.
Il risultato è un GeoPackage dove ogni riga è un poligono quadrato con i suoi attributi di rischio. La query PostGIS diventa banale:
SELECT risk_class, risk_class_label, mean_velocity, risk_index
FROM rischio_subsidenza_100m
WHERE ST_Intersects(
geom,
ST_Transform(ST_SetSRID(ST_MakePoint(:lon, :lat), 4326), 3035)
)
LIMIT 1;
Con un indice GIST, i tempi di risposta sono sotto i 5 ms.
Cosa mostrano i dati
Abbiamo prodotto una serie di visualizzazioni per esplorare il dataset. Ecco i risultati principali.
La mappa nazionale
La mappa mostra tutti gli 8,5 milioni di punti colorati per velocità verticale, dal rosso (subsidenza) al teal (sollevamento). I tre hotspot principali emergono con chiarezza:
- Pianura Padana: la grande macchia rossa nel nord, con subsidenza diffusa legata all'estrazione di fluidi sotterranei e alla compattazione naturale dei sedimenti alluvionali
- Campi Flegrei: l'unica area con sollevamento significativo (fino a +48 mm/anno), causato dal bradisismo vulcanico
- Etna: cedimento del fianco orientale del vulcano
La radiografia latitudinale
Questo ridge plot stile Joy Division mostra la distribuzione della velocità verticale per fascia di latitudine, dal sud (in basso) al nord (in alto). Ogni cresta rappresenta una banda di ~0,25° di latitudine.
La Pianura Padana (attorno a Bologna, 44,5°N) mostra una coda rossa pronunciata verso la subsidenza. I Campi Flegrei (40,8°N) spingono una coda teal/cyan verso il sollevamento. Il resto d'Italia è sorprendentemente simmetrico attorno allo zero.
Lo zoom sulla Pianura Padana
La Pianura Padana è il caso studio più critico. Con oltre 2,2 milioni di punti EGMS nella regione, la subsidenza media è di −1,37 mm/anno con un massimo di −59,8 mm/anno. Le aree più critiche si concentrano lungo la fascia Bologna–Modena–Ferrara e sulla costa ravennate, dove la subsidenza si combina con l'innalzamento del livello del mare.
Campi Flegrei: il vulcano che respira
Il bradisismo dei Campi Flegrei è il fenomeno più estremo nel dataset: sollevamento fino a +130 mm/anno nell'epicentro presso la Solfatara di Pozzuoli. I cerchi concentrici a 2, 5 e 10 km mostrano come il gradiente di deformazione decresca rapidamente con la distanza. A 10 km dal centro, il moto è già quasi nullo.
L'anatomia del rischio per regione
Il ranking regionale rivela che l'Emilia-Romagna è la regione con la subsidenza media più alta (−1,91 mm/anno), seguita da Toscana (−1,53) e Veneto (−0,35). La Campania è l'unica regione con velocità media positiva (+0,74 mm/anno) per effetto dei Campi Flegrei — un caso unico in Europa.
Punti deboli e limiti
Essere trasparenti sui limiti è importante quanto mostrare i risultati. Questo dataset ha diversi punti di attenzione:
Copertura non uniforme
I punti EGMS non sono distribuiti uniformemente. Le aree urbane e le infrastrutture hanno densità di punti molto più alta rispetto alle aree rurali, montane o boschive. Questo significa che una cella "dato non disponibile" non implica assenza di rischio — implica assenza di misurazioni radar affidabili in quella zona.
Risoluzione temporale limitata
Il dataset fornisce un valore medio sull'intero periodo 2019–2023. Non è possibile distinguere eventi puntuali (un cedimento improvviso nel 2021) da trend graduali. L'accelerazione cattura parzialmente la tendenza temporale, ma non sostituisce una serie storica completa.
Sensibilità alla vegetazione
La tecnica InSAR funziona bene su superfici rigide (edifici, strade, roccia) ma perde coerenza su terreno vegetato. Le aree agricole della Pianura Padana, ad esempio, hanno meno punti di misura rispetto alle aree urbanizzate circostanti — proprio là dove la subsidenza potrebbe essere più intensa.
Soglie di classificazione arbitrarie
Le soglie a 2, 5, 10 e 20 mm/anno sono convenzioni ragionevoli ma non universali. La letteratura scientifica proposta varia: l'ARPAE Emilia-Romagna usa soglie diverse, l'EGMS stesso non fornisce classificazioni ufficiali. La nostra scelta è coerente con i valori usati in contesti di protezione civile, ma non è l'unica possibile.
Non è una perizia strutturale
Il dataset è uno strumento di screening, non una diagnosi. Un edificio che cade in classe "alto" non è necessariamente a rischio crollo — può essere su un terreno che si muove uniformemente senza conseguenze strutturali. Viceversa, un edificio in classe "trascurabile" ma con fondazioni degradate potrebbe avere problemi seri. L'InSAR misura il moto del suolo, non la vulnerabilità dell'edificio.
Stack tecnico
Tutto il codice è Python puro con dipendenze standard:
- pandas + numpy: manipolazione dati e calcoli vettorializzati
- geopandas + shapely: geometrie e I/O spaziale
- pyproj: trasformazione di coordinate EPSG:3035 ↔ EPSG:4326
- matplotlib: visualizzazioni (con tiles Carto Dark Matter via
contextily)
La pipeline gira su un laptop con 32 GB di RAM. Il passaggio più pesante è la conversione punto→poligono degli 8,5 milioni di record, che impiega circa 15 minuti. La scrittura del GeoPackage finale richiede altri 20 minuti.
Come usare il dataset
In QGIS
- Scarica il file GeoPackage dalla pagina download
- Trascina il file in QGIS → seleziona il layer
rischio_subsidenza_100m - Il CRS è EPSG:3035 — QGIS lo riproietta automaticamente
- Stilizza per
risk_class_labelcon una scala categoriale
In PostGIS
ogr2ogr -f "PostgreSQL" PG:"host=localhost dbname=geo" \
rischio_subsidenza_poligoni.gpkg \
-nln rischio_subsidenza_100m \
-lco GEOMETRY_NAME=geom \
-lco SPATIAL_INDEX=GIST
In Python
import geopandas as gpd
gdf = gpd.read_file("rischio_subsidenza.gpkg",
layer="rischio_subsidenza_100m")
# Filtra celle ad alto rischio nella Pianura Padana
padana = gdf[
(gdf["risk_class"] >= 3) &
(gdf.geometry.centroid.y.between(4900000, 5100000))
]
Download
Il dataset è disponibile in tre formati dalla pagina download dataset:
| Formato | Dimensione | Note |
|---|---|---|
| GeoPackage (.gpkg) | ~2,6 GB | Layer: rischio_subsidenza_100m, EPSG:3035 |
| Parquet | ~800 MB | Compressione zstd, senza geometria |
| GeoJSON | ~4 GB | EPSG:4326, per interoperabilità web |
Licenza: ODbL (Open Database License) — i dati sorgente EGMS sono rilasciati sotto la licenza Copernicus, che consente la redistribuzione con attribuzione.
Fonti e attribuzione
- EGMS — European Ground Motion Service, Copernicus Land Monitoring Service. Dati InSAR Sentinel-1, release 2019–2023, componente Ortho Verticale L3.
- Natural Earth — Confini amministrativi 10m per il clipping spaziale.
- Le soglie di classificazione sono allineate con la letteratura di protezione civile e le linee guida ARPAE per il monitoraggio della subsidenza in Emilia-Romagna.
Risorse correlate
Download dataset territoriali Italia — formati GIS gratuiti Rischio Idrogeologico Italia: Dataset ISPRA Gratuiti per Alluvioni e Frane Rischio Idrogeologico e Sismico in Italia: Analisi Completa Dataset GIS Italia: guida completa ai download gratuiti