[postgis-users]

David Blasby dblasby at refractions.net
Tue Sep 16 09:29:10 PDT 2003


Mark Cave-Ayland wrote:
> Hi everyone,
> 
> I have a database containing a number of polygons in lat/lon coordinates
> for which I would like to determine the area in metres. Unforunately the
> area() function only calculates area in the given coordinate system. How
> can I reproject the data to get the answer in metres? Is it possible to
> do this using the length_spheroid() / transform() functions?

Your best bet is to transform the polygon to a meter-based coordinate 
system, then compute the area of the polygon.

If you know an approriate one for your data, ensure that both your
lat/lon coordinate system and new (meters-based) system are in
the spatial_ref_sys table.

(There should be a .sql file on the website (or the download) that 
should define several 1000 coordinate systems).

Then, to do the transform, make sure your original polygon has the 
correct SRID (so it points to the correct lat/lon coordinate system in 
spatial_ref_sys).

SELECT area(tranform(the_geom, ###)) FROM <table>;

Where ### is the SRID (in spatial_ref_sys) of your meter's coordinate 
system.

------

If you really want to get fancy, you can make a custom projection for 
each polygon.  There are some transformations that you can 
build-on-the-fly that are most accurate at a specific location on the 
globe.  I cannot remember the one that is best for doing area 
calculations, but alber equal area project comes to mind.  You can read 
up at the PROJ4 homepage ("http://www.remotesensing.org/proj").

To do this, you'll probably need to have the lat/lon location of the 
middle of your polygon.

(xmin(the_geom::box3d) + xmax(the_geom::box3d))/2,
(ymin(the_geom::box3d) + ymax(the_geom::box3d))/2

Then build a PROJ4 string that defines a projection centered at that 
location.

Since you're using a custom proj4 string (instead of one in 
spatial_ref_sys), you'll have to use the tranform_geom() function::

tranform_geom( GEOMETRY, TEXT (input proj4), TEXT (output proj4), INT 
(output srid) )


dave






More information about the postgis-users mailing list