Vous n'êtes pas identifié(e).

#1 Re : Général » Changer SRID mais absence de SRID » 27/10/2012 18:13:03

Petite illustration pratique du mécanisme des SRID dans PostGIS et généralisable à tous les SIG.

Prenons une table "transports" dans PostGIS. On peut visualiser ses métadonnées :

formation=# select * from geometry_columns where f_table_name = 'transports';
 f_table_catalog | f_table_schema | f_table_name | f_geometry_column | coord_dimension | srid  |      type       
-----------------+----------------+--------------+-------------------+-----------------+-------+-----------------
                 | public         | transports   | the_geom          |               2 | 27572 | MULTILINESTRING
(1 ligne)

On voit que le SRID est 27572. On peut vérifier que ce SRID est bien présent dans chaque géométrie :

formation=# select gid, st_srid(the_geom) from transports limit 10;
 gid | st_srid 
-----+---------
   1 |   27572
   2 |   27572
   3 |   27572
   4 |   27572
   5 |   27572
   6 |   27572
   7 |   27572
   8 |   27572
   9 |   27572
  10 |   27572
(10 lignes)

C'est bien le même. Si on regarde les objets de la base, on peut voir qu'il y a une contrainte sur la table transports dont la définition est :

ALTER TABLE transports
  ADD CONSTRAINT enforce_srid_the_geom CHECK (st_srid(the_geom) = 27572);

Donc à chaque fois qu'on insère ou modifie de la donnée de cette table, on vérifie que le SRID est bien 27572. Une table ("propre") ne contient donc que des géométries ayant le même SRID.

Si on va voir du côté de spatial_ref_sys :

formation=# select srid, auth_name, auth_srid  from spatial_ref_sys where srid = 27572;
 srid  | auth_name | auth_srid 
-------+-----------+-----------
 27572 | EPSG      |     27572
(1 ligne)

Cela signifie que mon SRID correspond à un système de coordonnées standardisé par l'EPSG, et que le SRID de ce système pour cet organisme est également 27572. Il peut parfois y avoir des différences entre le numéro de l'organisme et celui utilisé dans PostGIS comme SRID, mais c'est rare (actuellement aucun dans la table spatial_ref_sys par défaut).


Il s'agit donc de ce système de référence là :

http://spatialreference.org/ref/epsg/27572/


Allons voir la définition WKT de ce SRID. WKT est un format textuel standardisé par l'OGC de représentation des géométries, et également des système de projection.

formation=# select srid, srtext  from spatial_ref_sys where srid = 27572;

 srid  |                                                                                                                                                                                                                                                                                                                                                          srtext                                                                                                                                                                                                                                                                                                                                                          
-------+--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 27572 | PROJCS["NTF (Paris) / Lambert zone II",GEOGCS["NTF (Paris)",DATUM["Nouvelle_Triangulation_Francaise_Paris",SPHEROID["Clarke 1880 (IGN)",6378249.2,293.4660212936265,AUTHORITY["EPSG","7011"]],TOWGS84[-168,-60,320,0,0,0,0],AUTHORITY["EPSG","6807"]],PRIMEM["Paris",2.33722917,AUTHORITY["EPSG","8903"]],UNIT["grad",0.01570796326794897,AUTHORITY["EPSG","9105"]],AUTHORITY["EPSG","4807"]],UNIT["metre",1,AUTHORITY["EPSG","9001"]],PROJECTION["Lambert_Conformal_Conic_1SP"],PARAMETER["latitude_of_origin",52],PARAMETER["central_meridian",0],PARAMETER["scale_factor",0.99987742],PARAMETER["false_easting",600000],PARAMETER["false_northing",2200000],AUTHORITY["EPSG","27572"],AXIS["X",EAST],AXIS["Y",NORTH]]

C'est un peu dense, mais on voit qu'on retombe bien sur la définition donnée par le site spatialreference.org ( http://spatialreference.org/ref/epsg/27572/prettywkt/ )

En plus lisible :

PROJCS["NTF (Paris) / Lambert zone II",
    GEOGCS["NTF (Paris)",
        DATUM["Nouvelle_Triangulation_Francaise_Paris",
            SPHEROID["Clarke 1880 (IGN)",6378249.2,293.4660212936269,
                AUTHORITY["EPSG","7011"]],
            TOWGS84[-168,-60,320,0,0,0,0],
            AUTHORITY["EPSG","6807"]],
        PRIMEM["Paris",2.33722917,
            AUTHORITY["EPSG","8903"]],
        UNIT["grad",0.01570796326794897,
            AUTHORITY["EPSG","9105"]],
        AUTHORITY["EPSG","4807"]],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    PROJECTION["Lambert_Conformal_Conic_1SP"],
    PARAMETER["latitude_of_origin",52],
    PARAMETER["central_meridian",0],
    PARAMETER["scale_factor",0.99987742],
    PARAMETER["false_easting",600000],
    PARAMETER["false_northing",2200000],
    AUTHORITY["EPSG","27572"],
    AXIS["X",EAST],
    AXIS["Y",NORTH]]

Que peut on comprendre de ces informations. Principalement plusieurs points :

- C'est un système de projection nommé "NTF (Paris) / Lambert zone II"
- Cette projection est de type "Lambert_Conformal_Conic_1SP", donc une projection conique
- Elle est référencée par l'EPSG sous le code 27572
- Elle est en unité métrique
- Elle est basée sur une ellipsoide de référence qui est la "NTF (Paris)". C'est l'objet modélisant la terre sur lequel on effectue la projection conique.
- L'ellipsoide, comme la projection, sont paramétrées. Ces paramètres servent à la bibliothèque PROJ pour pouvoir faire les calculs de projection.


Proj utilise d'ailleurs en interne une autre notation pour ces informations, qui est également dans la table spatial_ref_sys (et sur spatialreference.org) :

formation=# select srid, proj4text from spatial_ref_sys where srid = 27572;
 srid  |                                                                               proj4text                                                                                
-------+------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 27572 | +proj=lcc +lat_1=46.8 +lat_0=46.8 +lon_0=0 +k_0=0.99987742 +x_0=600000 +y_0=2200000 +a=6378249.2 +b=6356515 +towgs84=-168,-60,320,0,0,0,0 +pm=paris +units=m +no_defs 
(1 ligne)

Ce sont les mêmes informations représentées un peu différemment.


Nous avons fait le tour du mécanisme de SRID dans PostGIS. Le reste c'est de la théorie de la géodésie. Les plus courageux pourront se lancer dedans.



Au passage, merci pour ce thread, je comprends mieux pourquoi l'article de comparaison sql server spatial versus PostGIS est rempli d'inepties.

#2 Re : Général » Changer SRID mais absence de SRID » 26/10/2012 16:57:05

Bonjour.


Mettre des gros mots dans les phrases pour que personne comprenne n'autorise pas à dire des bêtises. Je me permets donc de corriger légèrement vos propos.


Commençons par le problème initial.
On parle de PostGIS < 2.0, le fonctionnement est différent (et un peu plus simple) pour du >= 2.0.


A priori Georgie a du faire un chargement de données géographiques sans spécifier le système de coordonnées d'origine. Le chargement a donc été fait en donnant un SRID à -1, qui est la valeur que PostGIS utilise pour dire que le SRID est inconnu.


Comme ce système d'origine est inconnu de PostGIS, on ne peut pas faire de reprojection. La table s'affiche dans QGIS car ce dernier imagine que si le SRID est inconnu on utilise le SRID 4326 (ou alors la projection par défaut du projet en cours). Le SRID 4326 est l'indication de coordonnées en latitude/longitude (donc non projeté). I


Il faut donc remettre tout ça en ordre.


Le SRID est stocké à deux endroits : dans la table geometry_columns (métadonnées sur la base géo), et dans chaque géométrie. Il faut donc modifier les deux. Ce qui est sympa c'est qu'il y a une fonction pour ça :


  UpdateGeometrySRID(varchar table_name, varchar column_name, integer srid)


La doc : http://postgis.refractions.net/document … ySRID.html


Par exemple ici :


  SELECT UpdateGeometrySRID('matable', 'the_geom', 4326);


Sinon il faut le faire à la main :

- Supprimer les contraintes sur la table en question (avec pgadmin par exemple c'est le plus simple). puis :
- update matable set the_geom = st_setsrid(the_geom, 4326);
- remettre les contraintes et les  infos dans geometry_columns (par exemple avec populate_geometry_columns())


Ensuite, corrigeons des idées fausses :


Un SRID est un identifiant de projection (sisi). Il est en fait un numéro qui référence une entrée dans une liste de projections (table spatial_ref_sys pour PostGIS).

Chaque entrée de cette liste contient la définition mathématique de ladite projection : le type de projection, les paramètres d'icelle, et l'ellipsoïde de référence utilisée pour la modélisation de la terre. Les SRID correspondent (quasi toujours) aux codes standardisés par l'EPSG.

À noter que certains SRID déterminent un système de coordonnée non projeté. C'est le cas de EPSG:4326, qui a des coordonnées en latitude, longitude, donc en angles référencés sur une sphéroïde de référence.

Ce SRID est stocké dans la table de métadonnées geometry_columns ET dans la sérialisation de la géométrie qui est stockée dans la colonne géométrique des tables elles-mêmes.

La définition de la projection de la table spatial_ref_sys est utilisée par PROJ pour effectuer les calculs de projection.



Donc en résumé :
- il y a un identifiant SRID de projection dans les métadonnées des tables
- il y a un identifiant SRID de projection dans les géométries d'une colonne de type geometry
- Un SRID correspond à une projection
- Une projection utilise une sphéroïde de référence mais contient plus que ça (sauf cas de coordonnées sphériques).
- Les codes de projections sont (pour la plupart) standardisés, et les données d'un même ensemble partagent généralement la même projection.


Au plaisir,
Vincent

#3 Re : Général » Requête SQL - POSTGIS » 21/09/2012 15:37:36

Bonjour,

Si vous n'avez pas les droits sur le serveur et que votre administrateur ne veut pas installer l'extension, vous ne pourrez pas l'utiliser. Point de salut de ce coté, ce qui est dommage car elle correspond à votre besoin.

Il y a plusieurs autres pistes :
1/ si le langage PL/R est installé et que vous avez le droit de l'utiliser, il y a des fonctions de clustering dans R qui répondront à votre besoin. Peu de chance que ce soit accessible par défaut.
2/ Snapper les points sur une grille, et les regrouper avec un aggrégat sur la géométrie (fonction st_snaptogrid)
3/ faire une grille géométrique, intersecter les points avec cette grille, et regrouper ceux qui sont au dessus d'un certain seuil au centroide de l'élément de grille
4/ Implémenter un algo de clustering en Pl/pgsql.

La 3 vaut mieux que la 2 dans le sens ou les points correspondant aux zones non-denses ne seront pas affectés.
On peut aussi s'en sortir de la meme façon pour la 2,en reprenant tous les points qui sont uniques dans l'agrégats en leur redonnant leur position initiale.

Pas mal d'éléments de réponse, y compris pour la 4, ici :
http://gis.stackexchange.com/questions/ … th-postgis

J'espère que ça aide, ce ne sont que des pistes. Il y a surement d'autres moyens aussi.

#4 Re : Général » Optimisation du calcul du nombre de points dans une surface projetée » 23/02/2012 01:24:01

Jean-Marie a écrit :

merci pour précisions que je vais étudier

les projections ne sont pas pré-définies
la surface peut varier d'une région à la planète entière

dans le cas de la projection sphérique Mercator il est je semble à priori possible
de travailler directement sur l'expression en degré des limites de la surface projetée
ce qui correspond à une requête ayant la syntaxe suivante :

select count(*) from table where coord && ST_
Transform(ST_SetSRID(ST_MakeBox2D(ST_Point(xmin,ymin), ST_Point(xmax, ymax)),900913), 4326) and date ...

ce qui nous permet d'avoir un traitement plus rapide et d'autant plus intéressant que plus la surface
est grande plus il y a de points

Jean-Marie

Heu non, c'est exactement la même chose, faire des calculs de géométrie euclidienne (intersection ici), donc en 2D, avec des coordonnées sphériques (lat/lon, 4326) n'a pas de sens, le résultat sera faux.

#5 Re : Général » Optimisation du calcul du nombre de points dans une surface projetée » 21/02/2012 19:55:37

Jean-Marie a écrit :

Bonjour

effectivement il n'y a pas d'index utilisé ;
mais créer un index sur le st_transform ne me parait pas envisageable car nous sommes amenés à utiliser différentes projections

«différentes» étant combien ? Sont elles nombreuses ? Sont elles prédéfinies ? Quelles sont les zones sur lesquelles vous travaillez, leurs étendues, et leur éventuelle indépendance ? Quelle est la précision requise des calculs d'intersection ? Quel est le ratio lecture/écriture ?

Pour rendre le travail plus simple et plus performant on segmente souvent un domaine en sous domaines, avec des projections spécifiques à chaque fois. Les imports/exports de données vont éventuellement faire automatiquement l'enregistrement dans le bon domaine (= la bonne table = la bonne projection).

Une autre approche dans le cas où la compartimentation des données n'est pas possible, est d'effectivement utiliser seulement des coordonnées géodésiques (sphériques) en latitude / longitude. PostGIS permet de travailler sur de telles géométries avec le type «geography». Cependant ce type de données n'est pas supporté par toutes les fonctions de PostGIS, les calculs géométriques sont plus lents et certains calculs sont faits par projection locale et reconversion sphérique. Tout cela est documenté dans le manuel :

http://postgis.org/documentation/manual … yFunctions

http://postgis.org/documentation/manual … _Geography

Jean-Marie a écrit :

pour revenir à la requête plus rapide mais fausse, je remarque que la conversion en degré des limites projetés
renvoie 5 points (voir ci-dessous la description du polygone) ; il serait intéressant d'obtenir un polygone comportant davantage de points afin d'obtenir une meilleure précision) ; mais je ne sais pas comment faire

SELECT ST_AsText(ST_Transform(ST_SetSRID(ST_MakeBox2D(ST_Point(-111740.67231994, 452481.05159311), ST_Point(1763272.204062, 1722480.3657935)),27572), 4326)) As wgs_geom;

wgs_geom                                                                 
-------------------------------------------------------------------------------------------------------------------------------------------
POLYGON((-4.86421028855627 30.9545265233267,-6.2658910764966 42.1561991488799,16.3040533297558 41.5763321847958,14.0513693532454 30.482942611268,-4.86421028855627 30.9545265233267))

La conversion ne renvoie que 4 points et pas 5 : le dernier point d'un ring de polygone doit être le même que le premier pour le fermer. C'est logique : vous créez un rectangle avec makebox2d, et vous faites une conversion, les 4 coordonnées sont converties.

Le calcul d'intersection sur des géométries (type geometry) définies en EPSG 4326 (lat/lon) se fait dans un espace euclidien (projeté, 2D). Il est donc faux par définition, et d'autant plus faux que les polygones sont grands.

C'est sûrement ce qui se passe ici.

Il y a plusieurs pistes pour avancer, mais il manque des informations pour dire quelle est la meilleure. Globalement :
* on peut utiliser geography, qui fera un meilleur travail sur les données en lat/lon
* on peut compartimenter les données en zones projetées différemment, et indexer
* on peut réduire la taille des objets et augmenter leur nombre, ou les subdiviser (avec st_segmentize par exemple sur le makebox2D)
* on peut définir des indexes sur les transformations généralement utilisées

Mes 2€c,
Vincent

#6 Re : Général » Optimisation du calcul du nombre de points dans une surface projetée » 21/02/2012 00:00:07

Bonjour,
Vous pouvez vérifier avec EXPLAIN si il y a bien des indexes géographiques qui sont utilisés.
À priori non, puisque même si vous en avez défini un sur la géométrie des points, le st_transform empêche de l'utiliser.
La solution de Marc de créer un index sur le st_transform devrait résoudre la situation.

Quant au résultat faux, comment l'évaluez vous ?

#7 Re : Général » Remplir un champ par jointure spatiale » 21/09/2011 11:28:16

Salut,

Pascal24 a écrit :

Bonjour,
J’ai fais une fonction et un trigger qui se créent bien sans erreur, mais ça ne me renvoi pas la donnée du champ num et j’ai un message d’erreur sur le terme « geometry ». J’ai probablement fais une erreur dans ma requête ma là je coince.
Merci de votre aide.

«geometry» est un nom de type. Est ce que votre colonne de type geometry ne s'appelle pas plutôt «the_geom», comme il est coutume de le faire ?

Cdt,
Vincent

#8 Re : Général » [Postgis] Construction de requête SQL » 29/06/2011 15:27:09

Salut,
Quelques réponses ci-dessous.

Routchino a écrit :

Bonjour,
Mes questions :
-    j’ai dans la plupart de mes tables .txt deux colonnes avec des coordonnées X, Y. Dois-je faire une opération particulière pour ajouter une colonne the_geom comme c’est le cas pour les shapes (cela s’est fait automatiquement pour ces tables) ?

Oui il faut effectivement ajouter une colonne de géométrie. Il faut la créer avec la fonction "addgeometrycolumn", puis faire un update pour la remplir avec les x, y avec la fonction "st_makepoint". Cette dernière prend des float en argument, il faudra peut etre faire une conversion si les x,y sont sous forme de chaine de caractère. Dans ce cas un "st_makepoint(x::double precision, y::double precision)". Mettre le tout dans un "st_setsrid" pour spécifier le système de projection.

Routchino a écrit :

-    Dois-je attribuer la même projection à toutes mes tables ?

Dans l'absolu non.

Dans la pratique la réponse est plutôt oui. La plupart des fonctions PostGIS ne peuvent travaillent sur des géométries que si elles sont dans le même système de projection. Reprojeter coute assez cher, donc il vaut mieux le faire une fois pour toute. On peut aussi avoir deux colonnes dans une table avec la meme geometrie dans deux systèmes différents.
Voir : st_transform, addgeometrycolumn

Routchino a écrit :

Ex :
SELECT *
FROM wifi_dep_24, wifi_dep_33……………
WHERE code_insee = '24001'

Problème, la requête ne peut pas se faire car le champ code_insee apparait dans plusieurs tables. Comment puis-je donc réaliser cette requête ?

Là il faut lire la doc SQL.

select t1.*, t2.code_insee as code_insee_co2 from  wifi_dep_24 as t1, wifi_dep_33 as t2 ... WHERE code_insee = '24001'

Mais il faut revoir les jointures aussi, je doute que le cross join soit la vraie requete.

Routchino a écrit :

Serait-il possible de réaliser la requête dans l’autre sens, à savoir rechercher une communauté de commune dans la table de correspondance INSEE-communauté de commune, trouver les code INSEE associé et ensuite réaliser la requête de ces codes INSEE sur toute la table.

Oui c'est possible... Faut faire du SQL, peut etre une requete imbriquée, mais probablement que ce ne serait pas la peine.

Routchino a écrit :

-    pour les table ayant une colonne the_geom, pensez-vous qu’il est possible de faire une requête spatiale ou géométrique, toujours par rapport à la table de correspondance code INSEE-communauté de commune ? Si oui je n’ai aucune idée sur la façon de procéder donc je suis preneur de quelques pistes.

C'est possible, mais pour de la sélection de zones, mieux vaut se baser sur les données attributaires, ça règle un paquet de soucis.
Pour le géométrique il faut faire une jointure spatiale : faire des join on mafonctionspatiale(geom1, geom2)

Sinon voir : st_intersects, st_contains, st_dwithin …

Bon courage,
Vincent

#9 Re : Général » Update de tables Oracle vers PostgreSQL » 07/06/2011 19:02:23

Bonsoir Magalie, salut les gars,


Je suis d'accord avec flo pour dire que le problème principal est de repérer efficacement ce qui a bougé. Une fois ce souci résolu, qui est plutôt fonctionnel que technique, on peut chercher une solution pour la synchro des données.
La solution des trigger pour enregistrer les modifications dans une table spécifique me parait bonne.

Il faut ensuite un programme pour faire le transfert.
* FME
avantages : gère très bien les deux types de base, gère très bien le spatial, permet de faire le transfert en flux sans fichier au milieu, permet de faire les traitements de mise à jour dans la base destination, normalement.
inconvénients : pas très pratique pour automatiser la synchro, licence

* programme maison (python, java, whatever)
avantages : simple et efficace, pas de cout de licence
inconvénient : connaissances minimales en dev requises

Il faut aussi choisir dans le deuxième cas un mécanisme d'import/export
* CSV
C'est encore la solution la plus simple. Attention à dumper dans ce cas les géométries au format WKT, les formats internes binaires PG et Oracles diffèrent. Il faut alors reconstituer les geométries après le chargement dans PG (une ligne de SQL avec st_geomfromtext() )

* OGR2OGR
Permet de faire le transfert directement d'une base à l'autre car il gère Oracle et PG. on a aussi un transfert en flux dans ce cas. Il faut toutefois toujours utiliser un programme pour lancer le SQL qui va mettre à jour les tables de destination à partir des données chargées (qui sont sous forme de diff). Attention en général GDAL/OGR n'est pas compilé avec le support oracle par défaut.

* SHP
les outils oracle et pg peuvent lire et sortir du shapefile. On peut utiliser ce format comme format d'échange, mais on peut alors etre contraint par les limitations du shapefile (nb de caractères dans les noms de champs)

Perso je ferais CSV + Python, mais c'est par gout personnel. En efficacité pure peut etre que FME est à privilégier. Ça dépend du profil de celui qui met tout ça en place.

A retenir : attention aux dumps binaires pas compatibles entre les deux bases, en particulier pour les géométries.

Au passage si je me trompe pas, il me semble que je connais cette base... smile

Pied de page des forums

Propulsé par FluxBB