Difference between revisions of "Suporte Postgresql para rasters"

From Geomaster, Lda
Jump to navigation Jump to search
(Created page with "[https://github.com/lcalisto/workshop-postgis-raster Suporte Postgis para Raster], workshop realizada no SASIG 2017, no Porto.")
 
 
Line 1: Line 1:
 
[https://github.com/lcalisto/workshop-postgis-raster Suporte Postgis para Raster], workshop realizada no SASIG 2017, no Porto.
 
[https://github.com/lcalisto/workshop-postgis-raster Suporte Postgis para Raster], workshop realizada no SASIG 2017, no Porto.
 +
 +
http://suite.opengeo.org/docs/latest/dataadmin/pgGettingStarted/raster2pgsql.html
 +
 +
https://duncanjg.wordpress.com/2012/11/20/the-basics-of-postgis-raster/
 +
 +
=== Ortos 2012 ===
 +
 +
<syntaxhighlight>
 +
jgr@dusseldorf:/media/bonn/cme/ortos2012/Finais$ gdalinfo 16004200.tif
 +
Driver: GTiff/GeoTIFF
 +
Files: 16004200.tif
 +
      16004200.tfw
 +
Size is 17000, 17000
 +
Coordinate System is `'
 +
Origin = (-40050.000000000000000,125049.899999999994179)
 +
Pixel Size = (0.300000000000000,-0.300000000000000)
 +
Metadata:
 +
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
 +
  TIFFTAG_XRESOLUTION=2116
 +
  TIFFTAG_YRESOLUTION=2116
 +
Image Structure Metadata:
 +
  COMPRESSION=JPEG
 +
  INTERLEAVE=BAND
 +
Corner Coordinates:
 +
Upper Left  (  -40050.000,  125049.900)
 +
Lower Left  (  -40050.000,  119949.900)
 +
Upper Right (  -34950.000,  125049.900)
 +
Lower Right (  -34950.000,  119949.900)
 +
Center      (  -37500.000,  122499.900)
 +
Band 1 Block=128x128 Type=Byte, ColorInterp=Red
 +
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
 +
Band 2 Block=128x128 Type=Byte, ColorInterp=Green
 +
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
 +
Band 3 Block=128x128 Type=Byte, ColorInterp=Blue
 +
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
 +
Band 4 Block=128x128 Type=Byte, ColorInterp=Undefined
 +
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
 +
</syntaxhighlight>
 +
 +
=== Preparação ===
 +
 +
==== Atribuir uma projeção ====
 +
 +
<syntaxhighlight>
 +
#!/bin/bash
 +
for f in *.tif
 +
do
 +
echo "Processing $f"
 +
  # assign projection to $f
 +
gdal_edit.py -a_srs "EPSG:3763" $f
 +
done
 +
</syntaxhighlight>
 +
 +
[http://qgis.pt/apresentacoes_qgis2014/jag.pdf Notas do Prof. José Alberto Gonçalves]
 +
 +
[https://www.mapbox.com/help/processing-satellite-imagery/ Notas da MapBox], sobre processamento de imagens de satélite.
 +
 +
Exemplo "Ortos do Porto"
 +
 +
<syntaxhighlight>
 +
gdalbuildvrt estarreja.vrt *.tif
 +
 +
gdal_translate -b 1 -scale 31 255 0 255 estarreja.vrt R.tif
 +
gdal_translate -b 2 -scale 36 255 0 255 estarreja.vrt G.tif
 +
gdal_translate -b 3 -scale 40 255 0 255 estarreja.vrt B.tif
 +
 +
gdalbuildvrt -separate estarreja_all.vrt R.tif G.tif B.tif
 +
 +
# gdalwarp -t_srs "EPSG:3763" -dstalpha -co TILED=YES -co COMPRESS=JPEG estarreja_all.vrt estarreja.tif
 +
# gdalwarp -t_srs "EPSG:3763" -dstalpha -co TILED=YES -co COMPRESS=JPEG -multi estarreja_all.vrt estarreja.tif
 +
gdal_translate estarreja_all.vrt estarreja.tif -b 1 -b 2 -b 3 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES --config GDAL_TIFF_INTERNAL_MASK YES
 +
 +
gdaladdo -r average estarreja.tif 4 8 16 32 64 128 256
 +
</syntaxhighlight>
 +
 +
Alternativa (a funcionar com o QGIS e com o Geoserver)
 +
 +
<syntaxhighlight>
 +
gdalbuildvrt estarreja.vrt *.tif
 +
gdal_translate estarreja.vrt translate/estarreja.tif -b 1 -b 2 -b 3 -mask 4 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES --config GDAL_TIFF_INTERNAL_MASK YES
 +
gdaladdo -r average translate/estarreja.tif 4 8 16 32 64 128 256
 +
</syntaxhighlight>
 +
 +
==== Recriar os tiffs (o geotools estoirava nalgumas imagens) ====
 +
 +
http://gis.stackexchange.com/questions/89444/file-size-inflation-normal-with-gdalwarp
 +
 +
<syntaxhighlight>
 +
#!/bin/bash
 +
mkdir warp
 +
for f in *.tif
 +
do
 +
echo "Processing $f"
 +
  # gdal_translate -a_srs "EPSG:3763" $f warp/$f
 +
gdalwarp -t_srs "EPSG:3763" $f warp/$f
 +
done
 +
</syntaxhighlight>
 +
 +
==== Raster support PostgreSQL ====
 +
 +
Nota: Para adicionar ao QGIS, usar o DBManager, em vez do habitual ''Adicionar camada PostGIS''.
 +
 +
Variáveis de ambiente. http://postgis.net/docs/postgis_installation.html#install_short_version
 +
 +
<syntaxhighlight lang="bash">
 +
sudo vi /etc/postgresql/9.5/main/environment
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
POSTGIS_ENABLE_OUTDB_RASTERS=1
 +
POSTGIS_GDAL_ENABLED_DRIVERS=ENABLE_ALL
 +
</syntaxhighlight>
 +
 +
Ortos registados na base de dados (mas guardados fora da base de dados):
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -M -R -l 4 /media/bonn/cme/ortos2012/Finais/15504150.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
-a append
 +
-C enforce constraints
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -M -R -c -l 4 /media/bonn/cme/ortos2012/Finais/16004200.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -M -R -a -l 4 /media/bonn/cme/ortos2012/Finais/15504150.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -R -M /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
Dá... Mas fica muito lento!
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
Muito mais rápido, mas ainda não está perfeito.
 +
 +
[[File:Piramides com -l 2 4 8 16.png]]
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos
 +
</syntaxhighlight>
 +
 +
Muito mais rápido! Quase perfeito.
 +
 +
[[File:Piramides com -l 2,4,8,16,32,64,128,256.png]]
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256,512 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U sig ortos
 +
 +
</syntaxhighlight>
 +
 +
===== No servidor 192.168.100.45 =====
 +
 +
<syntaxhighlight lang="bash">
 +
root@cme-sig01:~# adduser geoserver
 +
20170223
 +
root@cme-sig01:~# mkdir ~geoserver/data
 +
root@cme-sig01:~# chmod a+rwx ~geoserver/
 +
root@cme-sig01:~# chmod a+rwx ~geoserver/data
 +
root@cme-sig01:~# mkdir ~geoserver/data/ortos
 +
root@cme-sig01:~# chmod a+rwx ~geoserver/data/ortos
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
scp estarreja.tif admsig@192.168.100.45:/home/geoserver/data/ortos
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
apt install postgis
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
sudo vi /etc/postgresql/9.6/main/environment
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
POSTGIS_ENABLE_OUTDB_RASTERS=1
 +
POSTGIS_GDAL_ENABLED_DRIVERS=ENABLE_ALL
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
sudo service postgresql restart
 +
</syntaxhighlight>
 +
 +
<syntaxhighlight lang="bash">
 +
raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256,512 /home/geoserver/data/ortos/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U sig ortos
 +
</syntaxhighlight>
 +
 +
=== Virtual layer ===
 +
 +
Funciona muito bem (no QGIS).
 +
 +
<syntaxhighlight lang="bash">
 +
cd /media/bonn/cme/ortos2012/Finais
 +
gdalbuildvrt all_tif.vrt *.tif -a_srs "EPSG:3763"
 +
</syntaxhighlight>
 +
 +
Gerar uma shapefile com os rectangulos onde encaixam os ortos:
 +
 +
<syntaxhighlight lang="bash">
 +
cd /media/bonn/cme/ortos2012/Finais
 +
gdaltindex -t_srs "EPSG:3763" mosaico.shp *.tif
 +
</syntaxhighlight>
 +
 +
==== Geoserver: ImageMosaic ====
 +
 +
Ao criar um store ImageMosaic no Geoserver, é gerada uma shapefile Finais.* (de índice, como a criada pelo gdaltindex) e um documento de configuração Finais.properties.

Latest revision as of 06:52, 20 February 2018

Suporte Postgis para Raster, workshop realizada no SASIG 2017, no Porto.

http://suite.opengeo.org/docs/latest/dataadmin/pgGettingStarted/raster2pgsql.html

https://duncanjg.wordpress.com/2012/11/20/the-basics-of-postgis-raster/

Ortos 2012

jgr@dusseldorf:/media/bonn/cme/ortos2012/Finais$ gdalinfo 16004200.tif
Driver: GTiff/GeoTIFF
Files: 16004200.tif
       16004200.tfw
Size is 17000, 17000
Coordinate System is `'
Origin = (-40050.000000000000000,125049.899999999994179)
Pixel Size = (0.300000000000000,-0.300000000000000)
Metadata:
  TIFFTAG_RESOLUTIONUNIT=2 (pixels/inch)
  TIFFTAG_XRESOLUTION=2116
  TIFFTAG_YRESOLUTION=2116
Image Structure Metadata:
  COMPRESSION=JPEG
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  -40050.000,  125049.900) 
Lower Left  (  -40050.000,  119949.900) 
Upper Right (  -34950.000,  125049.900) 
Lower Right (  -34950.000,  119949.900) 
Center      (  -37500.000,  122499.900) 
Band 1 Block=128x128 Type=Byte, ColorInterp=Red
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
Band 2 Block=128x128 Type=Byte, ColorInterp=Green
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
Band 3 Block=128x128 Type=Byte, ColorInterp=Blue
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16
Band 4 Block=128x128 Type=Byte, ColorInterp=Undefined
  Overviews: 8500x8500, 4250x4250, 2125x2125, 1062x1062, 531x531, 265x265, 132x132, 66x66, 33x33, 16x16

Preparação

Atribuir uma projeção

#!/bin/bash
for f in *.tif
do
	echo "Processing $f"
  	# assign projection to $f
	gdal_edit.py -a_srs "EPSG:3763" $f
done

Notas do Prof. José Alberto Gonçalves

Notas da MapBox, sobre processamento de imagens de satélite.

Exemplo "Ortos do Porto" 

gdalbuildvrt estarreja.vrt *.tif

gdal_translate -b 1 -scale 31 255 0 255 estarreja.vrt R.tif
gdal_translate -b 2 -scale 36 255 0 255 estarreja.vrt G.tif
gdal_translate -b 3 -scale 40 255 0 255 estarreja.vrt B.tif

gdalbuildvrt -separate estarreja_all.vrt R.tif G.tif B.tif

# gdalwarp -t_srs "EPSG:3763" -dstalpha -co TILED=YES -co COMPRESS=JPEG estarreja_all.vrt estarreja.tif
# gdalwarp -t_srs "EPSG:3763" -dstalpha -co TILED=YES -co COMPRESS=JPEG -multi estarreja_all.vrt estarreja.tif
gdal_translate estarreja_all.vrt estarreja.tif -b 1 -b 2 -b 3 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES --config GDAL_TIFF_INTERNAL_MASK YES

gdaladdo -r average estarreja.tif 4 8 16 32 64 128 256

Alternativa (a funcionar com o QGIS e com o Geoserver) 

gdalbuildvrt estarreja.vrt *.tif
gdal_translate estarreja.vrt translate/estarreja.tif -b 1 -b 2 -b 3 -mask 4 -co COMPRESS=JPEG -co PHOTOMETRIC=YCBCR -co TILED=YES --config GDAL_TIFF_INTERNAL_MASK YES
gdaladdo -r average translate/estarreja.tif 4 8 16 32 64 128 256

Recriar os tiffs (o geotools estoirava nalgumas imagens)

http://gis.stackexchange.com/questions/89444/file-size-inflation-normal-with-gdalwarp 

#!/bin/bash
mkdir warp
for f in *.tif
do
	echo "Processing $f"
  	# gdal_translate -a_srs "EPSG:3763" $f warp/$f
	gdalwarp -t_srs "EPSG:3763" $f warp/$f
done

Raster support PostgreSQL

Nota: Para adicionar ao QGIS, usar o DBManager, em vez do habitual Adicionar camada PostGIS.

Variáveis de ambiente. http://postgis.net/docs/postgis_installation.html#install_short_version 

sudo vi /etc/postgresql/9.5/main/environment

POSTGIS_ENABLE_OUTDB_RASTERS=1
POSTGIS_GDAL_ENABLED_DRIVERS=ENABLE_ALL

Ortos registados na base de dados (mas guardados fora da base de dados): 

raster2pgsql -s 3763 -d -I -C -M -R -l 4 /media/bonn/cme/ortos2012/Finais/15504150.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos

-a append -C enforce constraints 

raster2pgsql -s 3763 -d -I -M -R -c -l 4 /media/bonn/cme/ortos2012/Finais/16004200.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos

raster2pgsql -s 3763 -d -I -M -R -a -l 4 /media/bonn/cme/ortos2012/Finais/15504150.tif -F -t 1000x1000 ortos2012x1000x1000 | psql -h localhost -U geobox -d ortos

raster2pgsql -s 3763 -d -I -C -R -M /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos

Dá... Mas fica muito lento! 

raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos

Muito mais rápido, mas ainda não está perfeito.

File:Piramides com -l 2 4 8 16.png 

raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U geobox -d ortos

Muito mais rápido! Quase perfeito.

File:Piramides com -l 2,4,8,16,32,64,128,256.png 

raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256,512 /media/bonn/cme/ortos2012/Finais/translate/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U sig ortos
No servidor 192.168.100.45
root@cme-sig01:~# adduser geoserver
20170223
root@cme-sig01:~# mkdir ~geoserver/data
root@cme-sig01:~# chmod a+rwx ~geoserver/
root@cme-sig01:~# chmod a+rwx ~geoserver/data
root@cme-sig01:~# mkdir ~geoserver/data/ortos
root@cme-sig01:~# chmod a+rwx ~geoserver/data/ortos

scp estarreja.tif admsig@192.168.100.45:/home/geoserver/data/ortos

apt install postgis

sudo vi /etc/postgresql/9.6/main/environment

POSTGIS_ENABLE_OUTDB_RASTERS=1
POSTGIS_GDAL_ENABLED_DRIVERS=ENABLE_ALL

sudo service postgresql restart

raster2pgsql -s 3763 -d -I -C -R -M -l 2,4,8,16,32,64,128,256,512 /home/geoserver/data/ortos/estarreja.tif -F -t 100x100 public.ortos2012 | psql -h localhost -U sig ortos

Virtual layer

Funciona muito bem (no QGIS). 

cd /media/bonn/cme/ortos2012/Finais
gdalbuildvrt all_tif.vrt *.tif -a_srs "EPSG:3763"

Gerar uma shapefile com os rectangulos onde encaixam os ortos: 

cd /media/bonn/cme/ortos2012/Finais
gdaltindex -t_srs "EPSG:3763" mosaico.shp *.tif

Geoserver: ImageMosaic

Ao criar um store ImageMosaic no Geoserver, é gerada uma shapefile Finais.* (de índice, como a criada pelo gdaltindex) e um documento de configuração Finais.properties.

Retrieved from "https://sig.geomaster.pt/index.php?title=Suporte_Postgresql_para_rasters&oldid=147"