[OSM-dev] Help needed: postgis database transform()
Richard Duivenvoorde
rdmailings at duif.net
Mon Dec 1 20:42:47 GMT 2008
Hi Joerg,
apparently the schema you use is not properly 'spatial enabled'.
If you have a postgresql database, and installed postgis properly, you
either have to 'spatial enable' a database by using the right
'gis'-template, OR you run two big sql file which give your database:
- a lot of spatial functions
- a table called 'geometry_columns' (to 'registre' a geometry column:
give it the right srid/proj etc)
- a table called 'spatial_ref_sys' (which gives postgis the proper
parameters for the proj-functions to do reprojection etc).
So if your spatial_ref_sys is empty, you should run a script which is
called spatial_ref_sys.sql (here with me it's in
/usr/share/postgresql-8.2-postgis/spatial_ref_sys.sql)
The other script by the way is in the same directory and called
'lwpostgis.sql'.
select * from spatial_ref_sys where srid = 4326
gives me something like:
4326;"EPSG";4326;"GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS
84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]]";"+proj=longlat
+ellps=WGS84 +datum=WGS84 +no_defs "
So IF you want to transform to or from 4326, this record should be
available....
But besides this, if the import of your osm-data went ok, there is also
an srid in the geometry itself defined:
select srid(way) from planet_osm_point limit 5
you can also view this with:
select AsEWKT(way) from planet_osm_point limit 5
Now the clue is that to use spatial functions, both need a or the same srid.
By the way:
900913 is normally NOT in the spatial_ref_sys because it's not an
'official epgs' one. But you can put it in there, see:
http://spatialreference.org/ref/user/6/
INSERT into spatial_ref_sys (srid, auth_name, auth_srid, proj4text,
srtext) values ( 900913, 'spatialreference.org', 6, '+proj=merc
+a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0
+units=m +nadgrids=@null +wktext +no_defs',
'PROJCS["unnamed",GEOGCS["unnamed
ellipse",DATUM["unknown",SPHEROID["unnamed",6378137,0]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]],PROJECTION["Mercator_2SP"],PARAMETER["standard_parallel_1",0],PARAMETER["central_meridian",0],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1],EXTENSION["PROJ4","+proj=merc
+a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0
+units=m +nadgrids=@null +wktext +no_defs"]]');
now you can do:
select osm_id,
asEWKT(way),
asEWKT(transform(way,900913)),
asEWKT(transform(way,4326))
from planet_osm_point limit 1
Hope this helps,
Regards,
Richard Duivenvoorde
Joerg Ostertag (OSM Tettnang/Germany) wrote:
> Probably my problem boils down to :
> Which SRID can I use to get data into the coordinate system used inside the
> planet.osm?
>
> But here in Detail:
> I need a little help in reading from the mapnik posgis database. I plan to
> expand the osmtrackfilter.pl to compare against the local data already
> existing in my mapnik-postgis-database. So I tried to get all street-segments
> inside a bounding box. For this i tried to select from the database. But
> there I have the problem that the coordintesystems are different and I have
> no clue how to convert between them.
> I want to operate with coordinates compatible to the planet.osm File.
> Something like
> lat="48.8046469" lon="9.0401476"
>
> So i thought i'll start with
> echo "select osm_id,asText(way) from planet_osm_point limit 5;" | psql gis
> Result:
> osm_id | astext
> -----------+-------------------------------------------
> 293106148 | POINT(-2876147.15643226 4560658.25114169)
> 293105894 | POINT(-2874584.97662211 4563768.28055896)
> 293107003 | POINT(-2871694.69962705 4559831.77530905)
> 292977203 | POINT(-2871581.82166339 4556767.11400322)
> 292976398 | POINT(-2870606.29556972 4556208.15857941)
> (5 rows)
>
> Well at least I get some data from the database.
> Then the next step is to get transformed Data in the coordinate system I
> desire ... I did remember WGS-84 was SRID 4326, but either I'm wrong or
> whatever. Doing the following
> echo "select osm_id,asText(transform(way,4326)) from planet_osm_point limit
> 5;" | psql gis
> Results is an:
> ERROR: AddToPROJ4SRSCache: Cannot find SRID (4326) in spatial_ref_sys
>
>
> Well I thought I'll have a look into this in "spatial_ref_sys".... but ....
> echo "select * from spatial_ref_sys;" | psql gis
> Result:
> srid | auth_name | auth_srid | srtext | proj4text
> ------+-----------+-----------+--------+-----------
> (0 rows)
>
> Then I found SRID 900913 inside osm2pgsql and tried:
> echo "select osm_id,asText(transform(way,900913)) from planet_osm_point limit
> 5;" | psql gis
> Result:
> osm_id | astext
> -----------+-------------------------------------------
> 293106148 | POINT(-2876147.15643226 4560658.25114169)
> 293105894 | POINT(-2874584.97662211 4563768.28055896)
> 293107003 | POINT(-2871694.69962705 4559831.77530905)
> 292977203 | POINT(-2871581.82166339 4556767.11400322)
> 292976398 | POINT(-2870606.29556972 4556208.15857941)
> (5 rows)
>
> Well, this seems to be the same SRID which is used inside the Database.
>
> So my Question is:
> Which SRID can I use to get into the coordinate system I want?
>
>
> Jörg
>
> http://www.ostertag.name/
More information about the dev
mailing list