Archive

Archive for the ‘GIS’ Category

PostgresSQL / PostGIS CheatSheat

July 1st, 2007 rupert Comments off

This is a quick-command list of Postres. If you want detailed instructions, please visit the Postgres Manual.

1. How do I Show all databases?

postgres=# \l
List of databases
Name | Owner | Encoding
------------------+----------+----------
postgis | postgres | UTF8
postgres | postgres | UTF8
template0 | postgres | UTF8
template1 | postgres | UTF8
template_postgis | postgres | UTF8
(5 rows)

Note: Do not drop template databases if not necessary.

2. How do I run a script from the prompt?
\i
OR
psql -d cybersoftbj -u user < myfile.sql

Its very usuful in reloading user-defined functions.

3. How do I create a user?

CREATE ROLE lbs WITH LOGIN PASSWORD 'tracking' SUPERUSER INHERIT CREATEDB CREATEROLE;

4. How to dump database in a text file?

pg_dump -U lbs -d cybersoftbjv1 -h 127.0.0.1 -W > cybersoftbjv1.sql

4. How to dump database cleanly?

pg_dump -c -d -E UTF8 -h 127.0.0.1 -U lbs -W platform_v1 > platform_v1.sql

5. How to update using two tables?

UPDATE road_for_update u
SET the_geom = r.the_geom
FROM roads r
WHERE r.rd_id = u.rd_id;

6. How to change a column type with Cast?

ALTER TABLE roads ALTER COLUMN class_new TYPE integer USING class_new::integer;

7. How to provide/restrict access privileges to tables?

GRANT SELECT ON TABLE table TO user;
REVOKE SELECT ON TABLE table FROM user;

8. How to add a geometry column to a table?

SYNTAX: AddGeometryColumn(

,
, , , )
EXAMPLE: SELECT AddGeometryColumn(‘public’, ‘poi’, ‘the_geom’, 4326, ‘POINT’, 2)
9. Changing column names with spaces?
ALTER TABLE class_aroundme RENAME "level 1" TO level_1; 10. Setting kernel shmmax for postgres
sysctl -w kernel.shmmax=134217728
Note: For permanent changes see /etc/sysctl.cfg11. How to backup table(s) from pg_dump?
pg_dump poi_beijing -t class -t poi_class -f $BACKUPDIR/test_$MYDATE.sql

Categories: postgis, postgres Tags: ,

Installing QNavigator in QGIS

June 27th, 2007 rupert Comments off
[root@rupert-linux qnavigator]# mkdir build &amp;&amp; cd build
[root@rupert-linux build]# cmake -D QGIS_PREFIX=/usr/local/qgis ..
-- Check for working C compiler: /usr/bin/gcc
-- Check for working C compiler: /usr/bin/gcc -- works
-- Check size of void*
-- Check size of void* - done
-- Check for working CXX compiler: /usr/bin/c++
-- Check for working CXX compiler: /usr/bin/c++ -- works
-- Looking for Q_WS_X11
-- Looking for Q_WS_X11 - found
-- Looking for Q_WS_MAC
-- Looking for Q_WS_MAC - not found.
-- Looking for Q_WS_WIN
-- Looking for Q_WS_WIN - not found.
-- Found Qt-Version 4.2.3
-- Configuring done
-- Generating done
-- Build files have been written to: /home/install/qgis-svn/src/plugins/qnavigator/build
[root@rupert-linux build]# make
Scanning dependencies of target dglib
[  9%] Building CXX object dglib/CMakeFiles/dglib.dir/dgraph.o
[ 18%] Building CXX object dglib/CMakeFiles/dglib.dir/dijkstra.o
[ 27%] Building CXX object dglib/CMakeFiles/dglib.dir/fheap.o
[ 36%] Building CXX object dglib/CMakeFiles/dglib.dir/routing_core.o
Linking CXX static library libdglib.a
Creating makefiles...
[ 36%] Built target dglib
Scanning dependencies of target dgbuild
[ 45%] Building CXX object dglib/CMakeFiles/dgbuild.dir/dgbuild.o
Linking CXX executable dgbuild
[ 45%] Built target dgbuild
Scanning dependencies of target dgpath
[ 54%] Building CXX object dglib/CMakeFiles/dgpath.dir/dgpath.o
Linking CXX executable dgpath
[ 54%] Built target dgpath
[ 54%] Generating ../../python/ui/mainwindow_ui.py
[ 54%] Generating ../../python/ui/findstreetdialog_ui.py
[ 54%] Generating ../../python/ui/aboutdialog_ui.py
[ 54%] Generating ../../python/ui/roadinfodialog_ui.py
[ 54%] Generating ../../python/ui/qnavigator_rc.py
[100%] Built target ui
[root@rupert-linux build]#
Categories: GIS Tags: ,

Mapinfo Utility for gdal_translate.

June 21st, 2007 rupert Comments off

I made a small perl utility to automate the gcp’s from Mapinfo Raster TABS to gdal_translate command line. Currently your tabfile would have:

!table
!version 300
!charset WindowsLatin1
 
Definition Table
  File "beijing_6th_1.jpg"
  Type "RASTER"
  (116.522865,40.016316) (347,184) Label "Pt 1",
  (116.681215,40.015286) (7729,243) Label "Pt 2",
  (116.679777,39.777904) (7666,14674) Label "Pt 3",
  (116.523827,39.779108) (397,14606) Label "Pt 4"
  CoordSys Earth Projection 1, 104
  Units "degree"
  RasterStyle 4 1
  RasterStyle 7 1677695

How to use
gdal_mapinfo

1. ls *.TAB > init.sh

2. vi init.sh to reflect the ff:

perl gdal_mapinfo.pl Beijing_6th_1.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_10.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_11.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_12.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_13.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_2.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_3.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_4.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_5.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_6.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_7.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_8.TAB &gt;&gt; final.bat
perl gdal_mapinfo.pl Beijing_6th_9.TAB &gt;&gt; final.bat

3. The resulting final.bat should have the ff:
gdal_translate -gcp 347 184 116.522865 40.016316 -gcp 7729 243 116.681215 40.015286 -gcp 7666 14674 116.679777 39.777904 -gcp 397 14606 116.523827 39.779108 -of GTiff Beijing_6th_1.jpg I:\\satimages\translated\Beijing_6th_1_translated.tif

gdalwarp -s_srs epsg:4326 -t_srs epsg:4326 I:\\satimages\translated\Beijing_6th_1_translated.tif
I:\\satimages\warped\Beijing_6th_1.tif

Categories: GDAL/OGR Tags: , , , ,

Processing Mapinfo Raster JPEG Images using GDAL

June 14th, 2007 rupert Comments off

I have a couple of sat images (raw jpegs) from Google that I want to use with Openlayers/Mapserver. The raw jpegs were registered using Mapinfo via GCP (Ground Control Points).

Mapinfo Raster JPEG Images example:

rupert@rupert-winxp /e/home/map/beijing/new/satimages$ ll
-rw-r–r– 1 rupert None 358 Jan 30 03:23 2NE1.TAB
-rw-r–r– 1 rupert None 7.1M Jan 30 02:38 2NE1.jpg
-rw-r–r– 1 rupert None 356 Feb 1 18:56 2NE2a.TAB
-rw-r–r– 1 rupert None 3.8M Feb 1 08:57 2NE2a.jpg

You cannot fully reference 2NE1.TAB as a Mapserver Layer. I tried to use 2NE1.jpg, but the problem its not georeferenced.

rupert@rupert-winxp /e/home/map/beijing/new/satimages
$ gdalinfo 2NE1.jpg
Driver: JPEG/JPEG JFIF
Size is 8650, 6744
Coordinate System is `'
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 6744.0)
Upper Right ( 8650.0,    0.0)
Lower Right ( 8650.0, 6744.0)
Center      ( 4325.0, 3372.0)
Band 1 Block=8650x1 Type=Byte, ColorInterp=Red
Band 2 Block=8650x1 Type=Byte, ColorInterp=Green
Band 3 Block=8650x1 Type=Byte, ColorInterp=Blue

The georeference coordinates of 2NE1.jpg, just like an ESRI World File, is found in 2NE1.TAB…

rupert@rupert-winxp /e/home/map/beijing/new/satimages
$ cat 2NE1.TAB
!table
!version 300
!charset WindowsLatin1
 
Definition Table
  File "2ne1.jpg"
  Type "RASTER"
  (116.38575,39.906105) (349,6619) Label "Pt 1",
  (116.390072,39.93201) (1160,317) Label "Pt 2",
  (116.42786,39.932296) (8210,253) Label "Pt 3",
  (116.4295878,39.90722318) (8522,6358) Label "Pt 4"
  CoordSys Earth Projection 1, 0
  Units "degree"

I found hurting myself in trying to create an ESRI world file from the current MAPINFO Raster TABS. So, I decided to go for GeoTIFF since its native in Mapserver. Using GDAL utilities my only problem is how to put a coordinate system and reference to the raster.

On windows, you could use Frank’s FWTools. For Linux, compile GDAL by source with python would be extremely helpful later on. For installation of GDAL on Linux, we can use Mapserver’s Verbose Installation in Linux Guide.

GDAL – the saviour!.

GDAL utilities is extremely helpful in reprojection, scaling, image mosaics, etc. For now, we will use gdal_translate and gdal_warp. Please RTFM the utilities.

1. Using gdal_translate to specify the gcp’s registered in Mapinfo.

gdal_translate -gcp pixel line easting northing
Add the indicated ground control point to the output dataset. This option may be provided multiple times to provide a set of GCPs.



$ gdal_translate -gcp 349 6619 116.38575 39.906105 -gcp 1160 317 116.390072 39.93201 -gcp 8210 253 116.42786 39.932296 -gcp 8522 6358 116.4295878 39.90722318 -of GTiff 2NE1.jpg 2NE1translated.tif

Input file size is 8650, 6744
0...10...20...30...40...50...60...70...80...90...100 - done.

Note: even if we specify the gcp’s, gdal_translate would not specify the corner coordinates of the tiff.

rupert@rupert-winxp /e/home/map/beijing/new/satimages
$ gdalinfo 2NE1translated.tif
Driver: GTiff/GeoTIFF
Size is 8650, 6744
Coordinate System is `'
GCP Projection =
GCP[  0]: Id=1, Info=
          (349,6619) -&gt; (116.38575,39.906105,0)
GCP[  1]: Id=2, Info=
          (1160,317) -&gt; (116.390072,39.93201,0)
GCP[  2]: Id=3, Info=
          (8210,253) -&gt; (116.42786,39.932296,0)
GCP[  3]: Id=4, Info=
          (8522,6358) -&gt; (116.4295878,39.90722318,0)
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 6744.0)
Upper Right ( 8650.0,    0.0)
Lower Right ( 8650.0, 6744.0)
Center      ( 4325.0, 3372.0)
Band 1 Block=8650x1 Type=Byte, ColorInterp=Red
Band 2 Block=8650x1 Type=Byte, ColorInterp=Green
Band 3 Block=8650x1 Type=Byte, ColorInterp=Blue

2. Use gdalwarp to reproject using the gcp and specify the coordinates.

The gdalwarp utility is an image mosaicing, reprojection and warping utility. The program can reproject to any supported projection, and can also apply GCPs stored with the image if the image is “raw” with control information.


$ gdalwarp -s_srs epsg:4326 -t_srs epsg:4326 2NE1translated.tif warped.tif
Creating output file that is 9422P x 5631L.
Processing input file 2NE1translated.tif.
:0...10...20...30...40...5050...60...70...80...90...

Let’s check after gdalwarp using gdalinfo…

$ gdalinfo warped.tif
Driver: GTiff/GeoTIFF
Size is 9422, 5631
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.2572235630016,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (116.383841930499160,39.933342296341991)
Pixel Size = (0.000004927579869,-0.000004927579869)
Metadata:
  AREA_OR_POINT=Area
Corner Coordinates:
Upper Left  ( 116.3838419,  39.9333423) (116d23'1.83"E, 39d56'0.03"N)
Lower Left  ( 116.3838419,  39.9055951) (116d23'1.83"E, 39d54'20.14"N)
Upper Right ( 116.4302696,  39.9333423) (116d25'48.97"E, 39d56'0.03"N)
Lower Right ( 116.4302696,  39.9055951) (116d25'48.97"E, 39d54'20.14"N)
Center      ( 116.4070558,  39.9194687) (116d24'25.40"E, 39d55'10.09"N)
Band 1 Block=9422x1 Type=Byte, ColorInterp=Red
Band 2 Block=9422x1 Type=Byte, ColorInterp=Green
Band 3 Block=9422x1 Type=Byte, ColorInterp=Blue

Sweet. Now, all we need to do is display the raster images in Mapserver/OpenLayers.

Specifying a raster image in Mapserver

LAYER
NAME “2NE1″
DATA “satimages/2NE1.tif”
TYPE RASTER
STATUS DEFAULT
END

LAYER
NAME “2NE2″
DATA “satimages/2NE2a.tif”
TYPE RASTER
STATUS DEFAULT
END

Here is the end result…
.

Lessons learned, I tried to specify coordinate extents using gdal_translate -a_ullr ulx uly lrx lry. Specifying the coordinates was subjective by just looking at the cursor location of the registered image in Mapinfo. It is still best to use the GCP’s. Simply put, we need to be accurate in specifying corner coordinates in raster images to project them accurately.

.

Spatial Databases

June 2nd, 2007 rupert Comments off

Enterprise SDMS (Spatial Database Management Systems) such as Oracle and IBM (Oracle XE and IBM DB2 Express C) moved their flagship database products with free versions just like MS SQL Server Express (MS SQL Server 2008 Spatial to be released on 2008).

GIM International made a product survey comparing Spatial Databases IBM DB2, MySQL, Oracle and Postgres here.. Although subconciously, postgis would fit most of my purposes such as integration with Mapserver and Routing (pgRouting), I am still hopeful that Oracle Spatial or IBM DB2 Spatial Extender might give other advantages that I might overlooked.

IBM made significant efforts to tie up DB2 Spatial Extender with other opensource products such as Geoserver. And DB2 is not as limiting as Oracle XE. Currently, DB2 can be installed upto dual (2) CPUS, 4 GB of Memory and unlimited number of databases which quickly topple Oracle XE’s limit of single (1) CPU and upto 1 GB of Memory.

Categories: GIS Tags: