After the upgrading to etch I had to recompile PostGIS (version 1.0.4) for the new version of PostgreSQL:
pic@stakhanov:~$ /usr/lib/postgresql/7.4/bin/postmaster --version
postmaster (PostgreSQL) 7.4.8
I followed the hints in my previous post but had to change something:
$ ./debian/rules config isn't valid any more, one must instead use $ ./debian/rules build and is forced to compile all PostgreSQL's sources.
After this, I tested the new installation with some free GIS data downloaded from Mapping Hacks:
pic@stakhanov:~/devel$ mkdir geodata
pic@stakhanov:~/devel$ cd geodata
pic@stakhanov:~/devel/geodata$ wget http://mappinghacks.com/data/world_borders.zip
--23:03:33-- http://mappinghacks.com/data/world_borders.zip
=> `world_borders.zip'
Resolving mappinghacks.com... 216.218.203.219
Connecting to mappinghacks.com|216.218.203.219|:80... connected.
HTTP request sent, awaiting response... 200 OK
Length: 3,310,260 (3.2M) [application/zip]
100%[====================================>] 3,310,260 293.94K/s ETA 00:00
23:03:47 (248.37 KB/s) - `world_borders.zip' saved [3310260/3310260]
pic@stakhanov:~/devel/geodata$ unzip world_borders.zip
Archive: world_borders.zip
inflating: world_borders.dbf
inflating: world_borders.shp
inflating: world_borders.shx
pic@stakhanov:~/devel/geodata$ /usr/lib/postgresql/7.4/bin/shp2pgsql -s 4326 world_borders.shp country > country.sql
Shapefile type: Polygon
Postgis type: MULTIPOLYGON[2]
pic@stakhanov:~/devel/geodata$ createdb testgis
CREATE DATABASE
pic@stakhanov:~/devel/geodata$ createlang plpgsql testgis
pic@stakhanov:~/devel/geodata$ # if the next command works,
pic@stakhanov:~/devel/geodata$ # PostGIS has been installed correctly
pic@stakhanov:~/devel/geodata$ psql -d testgis -f /usr/share/postgresql/7.4/contrib/lwpostgis.sql
BEGIN
[...]
COMMIT
pic@stakhanov:~/devel/geodata$ psql -d testgis -f /usr/share/postgresql/7.4/contrib/spatial_ref_sys.sql
BEGIN
[...]
COMMIT
VACUUM
pic@stakhanov:~/devel/geodata$ psql -d testgis -f country.sql
BEGIN
[...]
COMMIT
pic@stakhanov:~/devel/geodata$ psql testgis
Welcome to psql 7.4.8, the PostgreSQL interactive terminal.
Type: \copyright for distribution terms
\h for help with SQL commands
\? for help on internal slash commands
\g or terminate with semicolon to execute query
\q to quit
testgis=# select cntry_name, area, sum(area(transform(the_geom, 26591))) as calculated_area from country where cntry_name = 'Italy' group by cntry_name, area ;
cntry_name | area | calculated_area
------------+--------+-----------------
Italy | 301230 | 301199476062.85
(1 row)
testgis=# /* note that the "area" column (in square kilometers, included in
testgis*# shapefile) and the "calculate_area" column (in square meters) are
testgis*# practically equal */
Other examples, indipendent of the imported data:
testgis=# select length2d(geometryFromText('LINESTRING(12 34, 13 35)')) ;
length2d
-----------------
1.4142135623731
(1 row)
testgis=# /* that is the usual distance between that two points in a cartesian
testgis*# plan, but: */
testgis=# select length2d_spheroid(geometryFromText('LINESTRING(12 34, 13 35)'), 'SPHEROID["GRS_1980",6378137,298.257222101]') ; length2d_spheroid
-------------------
144013.139256991
(1 row)
testgis=# /* the result is much bigger because in this case the points'
testgis*# coordinates are interpreted as longitude, latitude on a big
testgis*# spheroid.
testgis*# Changing the longitude doesn't affect the result: */
testgis=# select length2d_spheroid(geometryFromText('LINESTRING(42 34, 43 35)'), 'SPHEROID["GRS_1980",6378137,298.257222101]') ;
length2d_spheroid
-------------------
144013.139256991
(1 row)
testgis=# /* but changing the latitude does: */
testgis=# select length2d_spheroid(geometryFromText('LINESTRING(12 64, 13 65)'), 'SPHEROID["GRS_1980",6378137,298.257222101]') ;
length2d_spheroid
-------------------
121397.522102769
(1 row)
testgis=# /* the same length using the distance function and (I hope)
testgis*# an appropriate SRID: */
testgis=# select distance(transform(setSRID(geometryFromText('POINT(12 34)'), 4326), 26591), transform(setSRID(geometryFromText('POINT(13 35)'), 4326), 26591)) ;
distance
------------------
144144.751580593
(1 row)
Sincerely, I have a lot of doubts regard of the usage of spatial reference identifier (SRID). I chose 4326 because It's suggested in many documents :-P , then, for calculating the Italy's area, I chose 26591 as I saw in spatial_ref_sys table that it is related to that country and has meters as unit of measurement.
Calculating every country's area with this last SRID yields big gaps respect to the values in the shapefile.
Otherwise, using directly 4326 results, for Italy's area, in a small 33.1005615782892:
testgis=# select cntry_name, area, sum(area(the_geom)) as calculated_area from country where cntry_name = 'Italy' group by cntry_name, area ;
cntry_name | area | calculated_area
------------+--------+------------------
Italy | 301230 | 33.1005615782892
(1 row)
What are these? Degrees?
An experiment shows that, in these conditions, the area is calculated as if the poligons were on a cartesian plan:
testgis=# select area(setSRID(GeometryFromText('POLYGON((12 34, 13 34, 13 35, 12 35, 12 34))'), 4326)) ;
area
------
1
(1 row)
But wasn't 4326 a SRID for longitude, latitude on spheroid? I don't understand.
Anyway, I think that these questions will remain without answer for a long. In fact in this period I am too busy and haven't free time to examine geospatial matters any more :-( .
