Trasformazioni con grigliati IGM: confronto tra i servizi di Regione Liguria e del Geoportale Nazionale

Ho deciso di resuscitare questo blog che da tempo non utilizzavo, essendomi spostato su altri strumenti, Twitter in primis,  più adatti all'uso che ne facevo.  Mi limitavo infatti a riportare, notizie interessanti sulle tecnologie geospaziali che trovavo in rete, senza commenti o approfondimenti.

In questa nuova serie invece, utilizzerò il blog per approfondire alcuni aspetti tecnici relativi alle tecnologie geospaziali che utilizzo per lavoro.

Per questo primo post ho preso spunto da un interessante articolo di Antonio Falciano pubblicato a febbraio sul blog di TANTO, dove venivano messi a confronto i risultati dei servizi di trasformazione di coordinate della Regione Basilicata con quelli del nuovo servizio messo a disposizione dal Geoportale Nazionale, riscontrando nelle trasformazioni effettuate scostamenti significativi.

Incuriosito (ed anche un poco preoccupato) dalla cosa ho deciso di effettuare un test con i servizi di Regione Liguria, servizi che conosco bene perché realizzati da me.




Ho seguito una metodologia analoga, anche se effettuata con strumenti diversi, a quella utilizzata da Antonio.

I servizi di conversione di Regione Liguria.

Si tratta di una procedura web che permette la trasformazione planimetrica di sistema di riferimento per singoli punti o per shapefile, con un limite massimo nella dimensione dello shape di 30 MB.

Per la trasformazione delle coordinate vengono utilizzati due servizi sviluppati in ambiente ASP.NET, uno per i punti e uno per gli shape.

Per esempio nel caso di punti viene chiamato un servizio REST del tipo:

http://www.cartografiarl.regione.liguria.it/Utilities/web/ConvertPoint.aspx?crsFrom=3003&crsTo=25832&x=1500000&y=4900000&format=json

I parametri sono:

• crsFrom: codice EPSG del sistema di riferimento dei dati di input
• crsTo: codice EPSG del sistema di riferimento desiderato in output
• x: coordinata x, o longitudine nel caso di coordinate geografiche, espresse rispettivamente in metri o gradi sessadecimali
• y: coordinata y, o latitudine nel caso di coordinate geografiche, espresse rispettivamente in metri o gradi sessadecimali
• format: formato di output, può essere json o xml

Un esempio di risposta del servizio è il seguente:

{
    "x": 499977.824,
    "y": 4899971.633,
    "z": 0.0,
    "id": null,
    "retCodeConv": 1,
    "crs": 25832,
    "coordSysRef": {
        "crs": 25832,
        "datum": 3,
        "proiezione": 3,
        "epsgCode": "25832" 
    }
}

Per l'operazione di trasformazione delle coordinate viene utilizzata una particolare versione del software CartLab3 (release 2.02.05) eseguibile tramite linea di comando e facilmente integrabile nel servizio web. 

Per quanto riguarda i grigliati sono stati utilizzati i grigliati IGM di tipo GR1  suddivisi in fogli 1:50000  e che coprono l'intero territorio regionale.

Procedura

Come detto precedentemente ho utilizzato una metodologia simile a quella utilizzata da Falciano e illustrata nell'articolo pubblicato su TANTO:

• Creazione di un livello composto da una semina di punti random situati all'interno del territorio regionale. Il livello è espresso in coordinate Gauss-Boaga.
• Trasformazione delle coordinate in ETRF89/UTM32 utilizzando i vari servizi.
• Confronto tra le coordinate dei livelli ottenuti e calcolo delle statistiche.

Di seguito  il dettaglio delle operazioni effettuate.

Creazione di uno livello puntuale PostGis mediante il comando: 

CREATE TABLE random_points (

    id    SERIAL PRIMARY KEY,

    geom  geometry

);

Creazione di 500 punti random, mediante l'utilizzo della funzione PostGres RandomPointsInPolygon:

INSERT into random_points (geom) SELECT RandomPointsInPolygon(geom, 500) from regione;

Ecco visualizzato in QGIS il risultato ottenuto:

Successivamente ho esportato con OGR il livello postgis in uno shapefile da dare in pasto ai servizi.

Le operazioni di trasformazione sono state effettuate con i seguenti servizi

• CartLab3 versione desktop
• Servizio Regione Liguria
• Servizio Geoportale Nazionale

Il sistema di riferimento di output è stato scelto in tutti i casi come ETRF89/UTM32.

Ho quindi importato con il comando shp2pgsql gli shapefile trasformati su PosGis creando tre livelli distinti, a partire dai quali ho creato una vista contenente le coordinate x e y espresse come campi numerici, i delta calcolati lungo i due assi e la distanza tra le coppie di punti omologhi:

CREATE VIEW random_view AS

SELECT cl3.id id,

ST_X(cl3.the_geom) x_cl3, ST_Y(cl3.the_geom) y_cl3,

ST_X(lig.the_geom) x_lig, ST_Y(lig.the_geom) y_lig,

ST_X(pcn.the_geom) x_pcn, ST_Y(pcn.the_geom) y_pcn,

(ST_X(cl3.the_geom) - ST_X(lig.the_geom)) delta_x_cl3_lig,

(ST_Y(cl3.the_geom) - ST_Y(lig.the_geom)) delta_y_cl3_lig,

(ST_X(pcn.the_geom) - ST_X(lig.the_geom)) delta_x_pcn_lig,

(ST_Y(pcn.the_geom) - ST_Y(lig.the_geom)) delta_y_pcn_lig,
ST_DISTANCE(cl3.the_geom, lig.the_geom) delta_cl3_lig,
ST_DISTANCE(pcn.the_geom, lig.the_geom) delta_pcn_lig
FROM random_cl3_gr cl3

inner join random_pcn_gr pcn on (cl3.id = pcn.id)

inner join random_lig_gr lig on (cl3.id = lig.id);

Ho infine calcolato le statistiche mediante le seguente select sulla vista:

SELECT
round(CAST(avg(abs(delta_x_pcn_lig)) AS numeric),3) abs_avg_delta_x_pcn_lig,
round(CAST(min(abs(delta_x_pcn_lig)) AS numeric),3) abs_min_delta_x_pcn_lig,
round(CAST(max(abs(delta_x_pcn_lig)) AS numeric),3) abs_max_delta_x_pcn_lig,
round(CAST(stddev(abs(delta_x_pcn_lig)) AS numeric),3) abs_stddev_delta_x_pcn_lig,
round(CAST(avg(abs(delta_y_pcn_lig)) AS numeric),3) abs_avg_delta_y_pcn_lig,
round(CAST(min(abs(delta_y_pcn_lig)) AS numeric),3) abs_min_delta_y_pcn_lig,
round(CAST(max(abs(delta_y_pcn_lig)) AS numeric),3) abs_max_delta_y_pcn_lig,
round(CAST(stddev(abs(delta_y_pcn_lig)) AS numeric),3) abs_stddev_delta_y_pcn_lig,
round(CAST(avg(delta_x_pcn_lig) AS numeric),3) avg_delta_x_pcn_lig,
round(CAST(min(delta_x_pcn_lig) AS numeric),3) min_delta_x_pcn_lig,
round(CAST(max(delta_x_pcn_lig) AS numeric),3) max_delta_x_pcn_lig,
round(CAST(stddev(delta_x_pcn_lig) AS numeric),3) stddev_delta_x_pcn_lig,
round(CAST(avg(delta_y_pcn_lig) AS numeric),3) avg_delta_y_pcn_lig,
round(CAST(min(delta_y_pcn_lig) AS numeric),3) min_delta_y_pcn_lig,
round(CAST(max(delta_y_pcn_lig) AS numeric),3) max_delta_y_pcn_lig,
round(CAST(stddev(delta_y_pcn_lig) AS numeric),3) stddev_delta_y_pcn_lig,
round(CAST(avg(delta_pcn_lig) AS numeric),3) avg_delta_pcn_lig,
round(CAST(min(delta_pcn_lig) AS numeric),3) min_delta_pcn_lig,
round(CAST(max(delta_pcn_lig) AS numeric),3) max_delta_pcn_lig,
round(CAST(stddev(delta_pcn_lig) AS numeric),3) stddev_delta_pcn_lig
FROM random_view;

Per le statistiche ho usato sia i valori con segno che i valori assoluti dei vari Δ in modo da prendere in considerazione l'effettiva distanza tra le coppie di punti.



Risultati

Nella seguente tabella sono riportati i risultati ottenuti arrotondati alla terza cifra decimale, non ho riportato il confronto tra CartLab3 Desktop e il servizio di Regione Liguria perché queste trasformazioni sono risultate, come mi aspettavo, identiche.

	ΔX (m)	ΔY (m)	|ΔX| (m)	|ΔY| (m)	Δ (m)
Avg	0.025	0.013	0.029	0.024	0.043
Min	-0.026	-0.022	0.000	0.000	0.027
Max	0.056	0.060	0.056	0.060	0.066
Std	0.024	0.025	0.019	0.016	0.012

Come si può vedere dalla tabella le differenze tra i due servizi sono nell'ordine dei centimetri, con scarti massimi assoluti di circa 6 cm, scarti minori rispetto a quelli evidenziati da Antonio, ma comunque superiori alle aspettative.

Tali scarti potrebbero essere dovute all'utilizzo di grigliati diversi.
I servizi di Regione Liguria utilizzano i grigliati GR1, basati sul modello ITALGEO99.

Dalle prove fatte da Maurizio Foderà in questo commento sembra che anche il GN utilizzi i grigliati GR1, se è così le differenze con il servizio di RL potrebbero essere imputabili al software utilizzato per la trasformazione anche se mi sembrano sinceramente eccessive.

Al riguardo ho inviato una mail a pcn@minambiente.it per avere informazioni sui dettagli tecnici e spero di avere a breve una risposta.

Sarebbe in effetti utile che i servizi di trasformazione riportassero chiaramente i dettagli tecnici sottostanti, in particolare la versione dei grigliati ed il software utilizzato.