Le portail francophone de la géomatique

hari · Fri 02 December 2016 15:09

Bonjour,

J'aimerias réappliquer les données georeference (bbox) à une image ou un SHP correspondant précisément à cette bbox.

Données georeferences bbox EPSG:4326 : (http://boundingbox.klokantech.com/)
coordonnées TSV : -1.171845, 45.719448, -1.080063, 45.769615)

J'ai exporter la bounding box ( EPSG:4326 ) d'un plan d'une commune en plusieurs tiles (avec python mercantile) puis je les ai assemblé en une image que j'ai ensuite vectorsié puis exporté en un PNG.
J'ai tenté de convertir le PNG en TIFF avec les informations geographiques (geotiff) de la bounding box pour ensuite le convertir en SHP pour l'utiliser dans geopandas.
#https://github.com/dwtkns/gdal-cheat-sheet
#Georeference an unprojected image with known bounding coordinates:
#gdal_translate -of GTiff -a_ullr "<top_left_lon> <top_left_lat> <bottom_right_lon> <bottom_right_lat>" -a_srs EPSG:4269 input.png output.tif

Pas sûr que cette commande soit la bonne pour réappliquer les données georeference:
$ gdal_translate -of GTiff -a_ullr -1.171845 45.719448 -1.080063, 45.769615 -a_srs EPSG:4326 input.png output.tiff

$ gdal_polygonize.py output.tiff -f "ESRI Shapefile" shpfile/out.shp

Seulement avec les commandes ci-dessus le fichier SHP est deformé dans geopandas et notebook-jupyter, ce qui ne me permet pas d'extraire des infos.

Ceci ci-dessous correspond-il au coordonnées bbox (TSV: -1.171845, 45.719448, -1.080063, 45.769615 ) dans GDAL ? :
-a_ullr "<top_left_lon> <top_left_lat> <bottom_right_lon> <bottom_right_lat>"

Dernière modification par hari (Fri 02 December 2016 17:24)

jonathan_g · Sat 03 December 2016 21:16

Bonsoir,

pourquoi utilises-tu -a_srs EPSG:4269 si tes données sont en 4326 ?

tu pourrais essayer les commandes suivantes stp :

Code:

gdaltindex -f "ESRI Shapefile"  out.shp output.tiff

ogrinfo -al out.shp

Jonathan

hari · Sun 04 December 2016 09:46

Bonjour,

Merci pour votre réponse.

J'ai utilisé -a_ullr EPSG:4326,
EPSG:4269 c etait pour l exemple donné sur github. :
https://github.com/dwtkns/gdal-cheat-sheet
Georeference an unprojected image with known bounding coordinates.

gdal_transalte :
-a_ullr ulx uly lrx lry:
Assign/override the georeferenced bounds of the output file. This assigns georeferenced bounds to the output file, ignoring what would have been derived from the source file. So this does not cause reprojection to the specified SRS.

Je vais essayé les commandes que vous m'avez donné dès que possible, je vous tiens au courant.

hari · Sun 04 December 2016 11:55

$ogrinfo -al out.shp

Code:

 ogrinfo -al shptest/out.shp
INFO: Open of `shptest/out.shp'
      using driver `ESRI Shapefile' successful.

Layer name: out
Metadata:
  DBF_DATE_LAST_UPDATE=2016-12-04
Geometry: Polygon
Feature Count: 1
Extent: (-1.171845, 45.719448) - (-1.080063, 45.769615)
Layer SRS WKT:
GEOGCS["GCS_WGS_1984",
    DATUM["WGS_1984",
        SPHEROID["WGS_84",6378137,298.257223563]],
    PRIMEM["Greenwich",0],
    UNIT["Degree",0.017453292519943295]]
location: String (254.0)
OGRFeature(out):0
  location (String) = output.tiff
  POLYGON ((-1.171845 45.719448,-1.171845 45.769615,-1.080063 45.769615,-1.080063 45.719448,-1.171845 45.719448))

je sais pas trop quoi penser du résultat.

Dans geopandas avec notebook-jupyter le rendu plot ne correspond pas à l'image source et le polygon est déformé.C'est peut être juste la prévisualisation qui est faussé.

Avec la commande (gdal_polygonize.py output.tiff -f "ESRI Shapefile" shpfile/out.shp ): l'image tiff représentant une bbox zoom16, le fichiers SHP fait 54M et je m'arrive pas à l'ourvrir avec josm (plugins geotools), Qgis étant pour 64bit je pourrais essayé dans une vm pour voir dans la différence avec geopandas mais cela ne n'arrange pas trop (il me semble avoir déja essayé).

Y aurait-il une autre logiciel pour présualisé le shp ?
Je pense retransformer de nouveau le SHP en image pour comparer, ça doit être possible avec gdal peut-être ?

Dernière modification par hari (Sun 04 December 2016 13:04)

jonathan_g · Sun 04 December 2016 12:37

Hello,

le Shapefile m'a l'air correct, pourrais-tu envoyer les infos de ton GeoTif

Code:

gdalinfo output.tiff

je ne connais pas geopandas pour ma part j'utlise QGIS

hari · Sun 04 December 2016 13:05

$ gdalinfo output.tiff

Code:

Driver: GTiff/GeoTIFF
Files: output.tiff
Size is 12800, 9956
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (-1.171845000000000,45.719448000000000)
Pixel Size = (0.000007170468750,0.000005038871033)
Metadata:
  AREA_OR_POINT=Area
  Software=www.inkscape.org
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  -1.1718450,  45.7194480) (  1d10'18.64"W, 45d43'10.01"N)
Lower Left  (  -1.1718450,  45.7696150) (  1d10'18.64"W, 45d46'10.61"N)
Upper Right (  -1.0800630,  45.7194480) (  1d 4'48.23"W, 45d43'10.01"N)
Lower Right (  -1.0800630,  45.7696150) (  1d 4'48.23"W, 45d46'10.61"N)
Center      (  -1.1259540,  45.7445315) (  1d 7'33.43"W, 45d44'40.31"N)
Band 1 Block=12800x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA 
Band 2 Block=12800x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA 
Band 3 Block=12800x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA 
Band 4 Block=12800x1 Type=Byte, ColorInterp=Alpha

jonathan_g · Sun 04 December 2016 13:35

Hello,

votre image est orientée sud->nord (pixel size y posiitif), est-ce que c'est ce que vous désirez ? sinon il faudrait inverser les positions Y lors de l'application de l'option -a_ullr

hari · Sun 04 December 2016 13:40

Pour info, voici le code utilisé qui n'est peut être pas le bon pour voir le shp dans jupyter-notebook avec geopandas:
testviewshp.ipynb :

Code:

import geopandas as gpd
df = gpd.GeoDataFrame.from_file("shpfile/out.shp") 
import matplotlib.pyplot as plt
ax=df.plot() 
plt.show()

hari · Sun 04 December 2016 14:05

Mince j'ai posté entre temps.

Merci pour votre aide !
je vais tester d' inverser les positions Y :
gdal_translate -of GTiff -a_ullr -1.171845 45.769615 -1.080063 45.719448 -a_srs EPSG:4326 input.png output.tiff

au lieu de:
gdal_translate -of GTiff -a_ullr -1.171845 45.719448 -1.080063 45.769615 -a_srs EPSG:4326 input.png output.tiff

Je vous tiens au courant.

Dernière modification par hari (Sun 04 December 2016 14:05)

hari · Sun 04 December 2016 14:24

Voila le resultat, le pixel size y est passé en négatif.

$gdal_translate -of GTiff -a_ullr -1.171845 45.769615 -1.080063 45.719448 -a_srs EPSG:4326 input.png output2.tiff

$ gdalinfo output2.tiff

Code:

Driver: GTiff/GeoTIFF
Files: output2.tiff
Size is 12800, 9956
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (-1.171845000000000,45.769615000000002)
Pixel Size = (0.000007170468750,-0.000005038871033)
Metadata:
  AREA_OR_POINT=Area
  Software=www.inkscape.org
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  -1.1718450,  45.7696150) (  1d10'18.64"W, 45d46'10.61"N)
Lower Left  (  -1.1718450,  45.7194480) (  1d10'18.64"W, 45d43'10.01"N)
Upper Right (  -1.0800630,  45.7696150) (  1d 4'48.23"W, 45d46'10.61"N)
Lower Right (  -1.0800630,  45.7194480) (  1d 4'48.23"W, 45d43'10.01"N)
Center      (  -1.1259540,  45.7445315) (  1d 7'33.43"W, 45d44'40.31"N)
Band 1 Block=12800x1 Type=Byte, ColorInterp=Red
  Mask Flags: PER_DATASET ALPHA 
Band 2 Block=12800x1 Type=Byte, ColorInterp=Green
  Mask Flags: PER_DATASET ALPHA 
Band 3 Block=12800x1 Type=Byte, ColorInterp=Blue
  Mask Flags: PER_DATASET ALPHA 
Band 4 Block=12800x1 Type=Byte, ColorInterp=Alpha

J'ai lancé la polygonisation, ça va prendre un peu temps pour vérifier:
gdal_polygonize.py output2.tiff -f "ESRI Shapefile" shpfile/out2.shp

Je vous tiens au courant.

Dernière modification par hari (Sun 04 December 2016 14:27)

hari · Sun 04 December 2016 15:12

Le rendu est maintenant bien orienté nord sud mais le shp est toujours déformé:

Image source compressé input.tiff:
[img]http://imageshack.com/a/img921/2642/BwmOSj.jpg[/img]
[img]http://picpaste.com/input-jo20DD9J.jpg[/img]

Image rendu geopandas plot out2.shp
[img]http://imageshack.com/a/img921/9820/3rZZ9z.png[/img]
[img]http://picpaste.com/out2-wg7RNDVC.png[/img]

Dernière modification par hari (Sun 04 December 2016 15:22)

Fichier(s) joint(s) :
Pour accéder aux fichiers vous devez vous inscrire.

hari · Sun 04 December 2016 16:00

Avec le code testviewshp2.ipynb

Code:

import shapefile as shp
import matplotlib.pyplot as plt

sf = shp.Reader("shpfile/out2.shp")

plt.figure()
for shape in sf.shapeRecords():
    x = [i[0] for i in shape.shape.points[:]]
    y = [i[1] for i in shape.shape.points[:]]
    plt.plot(x,y)
plt.show()

j' obtiens cette image qui est toujours déformé mais un peu moins.

[img]http://imageshack.com/a/img921/6522/iKcLyh.png[/img]

hari · Sun 04 December 2016 16:16

Je pense pouvoir dire que ma demande initial de réappliquer les données georeference (bounding box) à une image est résolue.

Données georeferences bbox EPSG:4326
coordonnées TSV : -1.171845, 45.719448, -1.080063, 45.769615

inversez les positions Y
gdal_translate -of GTiff -a_ullr -1.171845 45.769615 -1.080063 45.719448 -a_srs EPSG:4326 input.png output.tiff

Par contre, je dois comprendre un peu mieux geopandas (si vous avez des infos je suis preneur )

Merci Jonathan!

Dernière modification par hari (Sun 04 December 2016 17:17)

Communauté	Ressources	L'emploi	Le marché	L'association	Fils RSS
Forum Blogs Wiki Planet Sites amis	Agenda Biblio Fiches AMO-CNIG Paris PCGIS Téléchargements	Offres d'emploi Banque de CV Formations Démarche Métiers	Geo-Entreprises Geo-Communiqué Appels d'offres	L'association Qui sommes-nous ? Nous aider Bienfaiteurs Partenaires Historique Statuts	Tout le site Job Géomatique Services GeoBlogs ...	Faire un don Mentions légales

Le portail francophone de la géomatique

Annonce

#1 Fri 02 December 2016 15:09

GDAL : Réappliquez les données georeference (bounding box) à une image

#2 Sat 03 December 2016 21:16

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#3 Sun 04 December 2016 09:46

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

#4 Sun 04 December 2016 11:55

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#5 Sun 04 December 2016 12:37

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#6 Sun 04 December 2016 13:05

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#7 Sun 04 December 2016 13:35

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

#8 Sun 04 December 2016 13:40

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#9 Sun 04 December 2016 14:05

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

#10 Sun 04 December 2016 14:24

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#11 Sun 04 December 2016 15:12

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

#12 Sun 04 December 2016 16:00

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Code:

#13 Sun 04 December 2016 16:16

Re: GDAL : Réappliquez les données georeference (bounding box) à une image

Pied de page des forums