Ubuntu-16.04/Mint-18 GIS hackery

howto_gis_ubuntu.md

Add ubuntugis PPA:

$ sudo add-apt-repository ppa:ubuntugis/ubuntugis-unstable

Add POSTGIS repository:

$ wget --quiet -O - https://www.postgresql.org/media/keys/ACCC4CF8.asc | sudo apt-key add -
$ sudo sh -c 'echo "deb http://apt.postgresql.org/pub/repos/apt/ xenial-pgdg main" >> /etc/apt/sources.list.d/postgresql.list'

Update APT:

$ sudo apt-get update

Install GDAL:

$ sudo apt-get install gdal-bin python-gdal

Install Postgres+Postgis:

$ sudo apt-get install postgis postgresql-9.6-postgis-2.4 pgadmin3 postgresql-contrib-9.6 

Install QGIS:

$ sudo apt-get install qgis python-qgis

Install GDAL+Postgis dev packages:

$ sudo apt-get install libgdal-dev postgresql-server-dev-9.6

Install OGR Foreign Data Wrappers from https://github.com/pramsey/pgsql-ogr-fdw:

$ cd ~/gitrepos
$ git clone [email protected]:pramsey/pgsql-ogr-fdw.git
$ cd pgsql-ogr-fdw/
$ make 
$ sudo make install

TODO: Integrate WWW Foreign Data Wrappers examples from https://github.com/cyga/www_fdw/wiki/Examples:

$ cd ~/gitrepos
$ git clone [email protected]:cyga/www_fdw.git
$ cd www-fdw/
$ make 
$ sudo make install

Create user in Postgres:

$ sudo -u postgres createuser --pwprompt --superuser gisuser

Create template in Postgres for Postgis:

$ createdb -h 127.0.0.1 -U gisuser -O gisuser template_postgis2.4
$ psql -h 127.0.0.1 -U gisuser -d template_postgis2.4 -c "CREATE EXTENSION postgis;"
$ psql -h 127.0.0.1 -U gisuser -d template_postgis2.4 -c "CREATE EXTENSION postgis_topology;"
$ psql -h 127.0.0.1 -U gisuser -d template_postgis2.4 -c "CREATE EXTENSION ogr_fdw;"
$ psql -h 127.0.0.1 -U gisuser -d template_postgis2.4 -c "CREATE EXTENSION www_fdw;"

Install PostGISAddons:

$ cd /tmp/
$ wget https://raw.githubusercontent.com/pedrogit/postgisaddons/master/postgis_addons.sql
$ psql -h 127.0.0.1 -U gisuser -d template_postgis2.4 -a -f postgis_addons.sql
$ rm postgis_addons.sql

Create a database based on the template:

$ createdb -h 127.0.0.1 -U gisuser -T template_postgis2.4 postgis_demo

Set GDAL raster output drivers to enabled as they ship disabled (doing this on the template doesn't propagate it through so need to do it each time):

$ psql -h 127.0.0.1 -U gisuser -d postgis_demo -c "ALTER DATABASE \"postgis_demo\" SET postgis.gdal_enabled_drivers TO 'ENABLE_ALL';"

Test access:

$ psql -h 127.0.0.1 -U gisuser -d postgis_demo -c "SELECT postgis_version();"
$ psql -h 127.0.0.1 -U gisuser -d postgis_demo -c "SELECT * FROM ST_GDALDrivers();"

Useful tutorial:

http://docs.qgis.org/2.8/en/docs/training_manual/database_concepts/data_model.html

Create QGIS database connection:

• Open QGIS
• In the Browser box, right click on PostGIS and choose New Connection...
• Make up a name for the connection and enter it into the Name box
• Enter localhost, or the IP address of the PostGIS database in Host
• Enter the database name, for example, postgis_demo from the instructions above
• Enter the username and password you created for that database, for example, gisuser from the instructions above
• Click on Save Username and Save Password if you want to avoid typing them in again
• Click on Test Connect to check that the details were correct
• Click on OK when done to save the connection

Importing shapefile using QGIS:

• Open QGIS
• Click on the Database menu and choose DB Manager->DB Manager
• Click on the arrow next to PostGIS and open the list
• Click on the array next to your database connection
• Click on the array next to public
• Click on the Table menu and choose Import layer/file
• Click on "..." to choose the file to import
• Type your desired table name into the Table field
• Check the box next to "Create spatial index"
• Click OK to import the file
• Wait for the dialog saying the import was successful to come up and click OK on it
• Click on public
• Click on the Database menu and choose Refresh to show the new table
• Close the DB Manager window
• Click on the arrow next to PostGIS and open the list
• Right click on the connection and choose Refresh
• Click on the arrow next to the connection name and open the list if now open
• Right click on public and choose Refresh to show the new table

Export from Postgres to a CSV file:

$ psql -d postgis_demo -U gisuer -c "\copy \"2011Census_I01A_AUST_IARE_short\" TO 'temp/test.csv' WITH CSV HEADER"

Remove a fixed substring from a Postgres field and assign to another field:

$ psql -d postgis_demo -U gisuser -c "UPDATE \"2011Census_I01A_AUST_IARE_short\" SET "region_id_number" = SUBSTRING("region_id", 5, LENGTH("region_id"));"

Merge two tables and export from Postgres to a CSV file:

$ psql -h 127.0.0.1 -U gisuser -d postgis_demo -c "\copy (SELECT * FROM \"2011Census_I01A_AUST_IARE_short\" AS c, \"IARE_2011_AUST\" AS s WHERE c.region_id_number = s.ia_code11) TO 'temp/test.csv' WITH CSV HEADER"

Merge two tables and export from Postgres to a SHP file:

$ pgsql2shp -f temp/test.shp -h localhost -u gisuser -P something postgis_demo "SELECT * FROM \"2011Census_I01A_AUST_IARE_short\" AS c, \"IARE_2011_AUST\" AS s WHERE c.region_id_number = s.ia_code11"

Note the above terminates with "buffer overflow detected" on my machine, so may need to try the following way.

$ ogr2ogr -f "ESRI Shapefile" temp/test.shp PG:"host=localhost user=gisuser dbname=postgis_demo password=something" -sql "SELECT * FROM \"2011Census_I01A_AUST_IARE_short\" AS c, \"IARE_2011_AUST\" AS s WHERE c.region_id_number = s.ia_code11"

Join vector table with statistics:

CREATE TABLE "2011Census_I01A_AUST_IARE_short_with_geom" AS
SELECT *
FROM "IARE_2011_AUST" AS g, "2011Census_I01A_AUST_IARE_short" AS v
WHERE g.ia_code11 = v.region_id_number;

CREATE INDEX "2011Census_I01A_AUST_IARE_short_with_geom_gist" ON "2011Census_I01A_AUST_IARE_short_with_geom" USING gist (geom);

Rasterising a vector:

$ psql "${POSTGRES_PSQL_CONNECTION}" -c "CREATE TABLE \"dummy_rast\" ( rid integer PRIMARY KEY , rast raster );"
$ psql "${POSTGRES_PSQL_CONNECTION}" -c "INSERT INTO \"dummy_rast\" ( rid, rast ) VALUES (1, ST_MakeEmptyRaster(1800, 1800, 0, 0, 0.1, -0.1, 0, 0, 4283));"
$ time psql "${POSTGRES_PSQL_CONNECTION}" -c "CREATE TABLE \"${2}_with_raster_${1}\" AS SELECT ST_ExtractToRaster(ST_AddBand(rast, '32BF'::text, -9999, -9999), 'public',  '${2}_with_geom', 'geom', '${1}', 'AREA_WEIGHTED_MEAN_OF_VALUES') rast FROM \"dummy_rast\";"

Possible way to export raster to file:

$ psql -h 127.0.0.1 -U gisuser -d postgis_demo -c "COPY (SELECT encode(ST_AsGDALRaster(rast, 'GTiff'), 'hex') AS png FROM \"2011Census_I01A_AUST_IARE_short_with_geom_raster_Tot_p_Indig_P\") TO '/tmp/test.hex';"
$ xxd -p -r /tmp/test.hex > /tmp/test.tiff
$ gdal_translate -of EHdr "/tmp/test.tiff" "temp/raster_export.bil"
$ (cd temp/ && zip -r "raster_export.zip" "raster_export".*)

Possible way to convert raster to vector:

SELECT  (ST_DumpASPolygons(rast)).* FROM rasttable

$ ogr2ogr -f "ESRI Shapefile" temp/test.shp PG:"host=localhost user=gisuser dbname=postgis_demo password=something" -sql "SELECT (ST_DumpASPolygons(rast)).* FROM \"2011Census_I01A_AUST_IARE_short_with_geom_raster_Density_p_Indi\""

Query WMSServer with gdalinfo

$ gdalinfo "WMS:https://spatial.industry.nsw.gov.au/arcgis/services/CrownLands/TravellingStockRoutesforNSW/MapServer/WMSServer"

Then query subdataset using a URL found by the above command

$ gdalinfo "WMS:https://spatial.industry.nsw.gov.au/arcgis/services/CrownLands/TravellingStockRoutesforNSW/MapServer/WmsServer?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&LAYERS=0&SRS=EPSG:4326&BBOX=140.999831,-37.263162,153.545131,-28.290220"

Generate an XML file (which may require edits to include specific information) from WMSServer:

gdal_translate "WMS:https://spatial.industry.nsw.gov.au/arcgis/services/CrownLands/TravellingStockRoutesforNSW/MapServer/WmsServer?SERVICE=WMS&VERSION=1.1.1&REQUEST=GetMap&LAYERS=0&SRS=EPSG:4326&BBOX=140.999831,-37.263162,153.545131,-28.290220" TravellingStockRoutes-WMS.xml -of WMS

Generate a JPEG file using the WMS XML file created above:

gdal_translate -of JPEG -outsize 500 250 TravellingStockRoutes-WMS.xml ./temp/TravellingStockRoutes-WMS.jpg