Discussion :

[Anne-LaureLbl]

Bonjour,
Ceci est mon premier message et j'espère que quelqu'un voudra bien me donner un début de réponse ? J'obtiens de grosses différences de résultats entre la méthode de calcul ST_Distance dans PostGis et la méthode classique Haversine. J'ai sûrement dû faire une erreur, mais où ?

Soit une table t_zon dans pgsql 9.6 (et pas 9.5 mais le préfixe 9.6 n'était pas disponible) avec extension postGIS avec des colonnes zon_id INT, zon_label TEXT, zon_latitude INT, zon_longitude INT, zon_geo geography(point).
On insère deux valeurs (49, 'France', 48.86, 2.35, 'point(48.867 2.35)') et (262, 'Mauritius', -20.2833, 57.55, 'point(-20.2833 57.55)').
Jusque ici tout va bien. Le problème, c'est lors de la recherche de la distance réelle, voici ce que j'obtiens :

1) Haversin
-- > résultat 9 413,520 Km

2) fonction ST_Distance
-- > résultat 8 554,608 Km, soit 900Km de différence !!

Voici le code effectué :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
CREATE EXTENSION postGIS;
CREATE EXTENSION citext;
 
CREATE SCHEMA S_ZON; 	-- > schema geographic zones				
 
 CREATE DOMAIN S_ZON.D_LATITUDE DOUBLE PRECISION					-- latitude coordinates
	CHECK (VALUE BETWEEN -90 AND 90);				
 CREATE DOMAIN S_ZON.D_LONGITUDE DOUBLE PRECISION					-- longitude coordinates
	CHECK (VALUE BETWEEN -180 AND 180);				
 CREATE DOMAIN S_ZON.D_ZOOM SMALLINT					-- zoom level display
	CHECK (VALUE BETWEEN 0 AND 23);						
 
CREATE SCHEMA S_REF;	-- > schema reference		
 
CREATE DOMAIN S_REF.D_NAME as CITEXT					-- max 64 insensitive alphanumeric
	CHECK (length(value) <= 64);		
 
CREATE TABLE S_ZON.T_ZON			-- > geographical zones in schema geographic destinations		
(   ZON_ID				INT 	PRIMARY KEY,							
ZON_label			S_REF.D_NAME		not null,					-- destination name
 ZON_latitude		S_ZON.D_LATITUDE	not null,			-- latitude coordinates
 ZON_longitude		S_ZON.D_LONGITUDE	not null,		-- longitude coordinates
	ZON_geo				geography(point)								-- postGIS type
 );
 
insert into S_ZON.T_ZON values (49, 'France', 48.86, 2.35, 'point(48.867 2.35)’), (262, 'Mauritius', -20.2833, 57.55, 'point(-20.2833 57.55)’);
 
/*POURQUOI LES RESULTATS SONT DIFFERENTS ??*/
-- 1) results 9 413,520 Km
SELECT zon_id, typ_id, zon_label, zon_latitude, zon_longitude,
     111.045* DEGREES(ACOS(COS(RADIANS(latpoint))
                * COS(RADIANS(z.zon_latitude))
                * COS(RADIANS(longpoint) - RADIANS(z.zon_longitude))
                + SIN(RADIANS(latpoint))
                * SIN(RADIANS(z.zon_latitude)))) AS distance_in_km
FROM s_zon.t_zon AS z
JOIN (
    SELECT  48.86  AS latpoint,  2.35 AS longpoint
  ) AS p ON 1=1
WHERE Z.zon_id = 262 ;
 
 
-- 2) results 8 554,608 Km
SELECT (ST_Distance(toZ.zon_geo,fromZ.zon_geo)) / 1000,
  ToZ.zon_label::text || ' - ' || fromZ.zon_label::text as distance_label_km
FROM s_zon.t_zon AS toZ CROSS JOIN s_zon.t_zon AS fromZ              
WHERE  toZ.zon_ID = 262 and fromZ.zon_ID = 49;
 
-- 3) results 8 554,608 Km
select st_distance( st_geographyfromtext ('POINT(48.867 2.35)'),
                   st_geographyfromtext ('point(-20.2833 57.55)'))/1000;
 
-- 4) results 8 554,608 Km
select st_distance( geometry ('0101000020E6100000AE47E17A146E4840CDCCCCCCCCCC0240'),
                   geometry ('0101000020E6100000F1F44A59864834C06666666666C64C40'))/1000;

Merci de votre aide !

[SQLpro]

C'est parfaitement normal... En effet vous avez inversé dans la saisie les coordonnées LATITUDE/LONGITUDE du mode geography :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
insert into T_ZON values (49, 'France', 48.86, 2.35, 'point(2.35 48.867)'), 
                         (262, 'Mauritius', -20.2833, 57.55, 'point(57.55 -20.2833)');

Et tout ira mieux !!!

A +

[SQLpro]

Ha au fait, joli avatar et bienvenue !

A +

[Anne-LaureLbl]

Merci beaucoup !!

[SQLpro]

Attention, le calcul avec la formule d'Haversine ne prend pas en compte l’aplatissement de la terre aux pôles. Il en résulte une approximation de la navigation orthodromique d'environ 0,5%.
Le calcul de PG est donc bien plus précis car il prend en compte l'imperfection de notre globe, à condition d'avoir utiliser le bon SRID !

A +

[Anne-LaureLbl]

Merci SqlPro, tout fonctionne bien comme cela, la distance affichée est de 9 406,769 Km

Code sql :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
select st_distance( st_geographyfromtext ('POINT(2.35 48.867)'),
                    st_geographyfromtext ('point(57.55 -20.2833)'))/1000;

J'ai réintégré toute ma base en indiquant (Long, Lat) au lieu de (Lat, Long) pour le ST_Geography et je pense que ça va rouler.
Merci pour tout !