[postgis-users] point_inside_circle(geometry,float,float,float) : units -> meters ?

[David Blasby]

Pedro Salazar wrote:

>Greetings,
>
>I would like if is it possible to use the
>point_inside_circle(geometry,float,float,float) function to obtain all
>the geometries (points) inside a circle but with a radius defined in
>meters and not in the native units?
>
There are lots of ways of doing this.  There are routines for doing 
distance calculations on an ellipsoid, but I think its more general to 
re-project your data.

Here's a fairly general method:

1. Find a cartesian projection thats accurate for where your data is 
(ie. UTM zone whatever)
2. make sure this projection is properly entered in your spatial_ref_sys 
table. Most likely its already there if you uploaded the standard epsg 
projections.
3. Optional (for indexing) : find a maximum conversion factor between 
your data's units and meters.  We;'ll be using this information to make 
a bounding box that's LARGER than the circle you want to enscribe

SELECT * FROM table WHERE
        the_geom && expand( <GEOMETRY>, <BOX SIZE> )
        AND distance(transform(the_geom,<SRID>), 
transform(<GEOMETRY>,<SRID>) <  <DIST>

This is a bit more general than what you're asking for, but:
   < GEOMETRY > = point thats the centre of your circle
   <BOX SIZE>      =  this is the step #3 conversion factor (degrees)
    <SRID>   =  step #2's SRID in the spatial_ref_sys table
    <DIST>  = distance (meters)

This is a bit complicated, but:

    the_geom && expand( <GEOMETRY>, <BOX SIZE> )
This uses the index to quickly find points that are nearby.  Box size 
will depend on step #3.  
NB: make sure this box size is ALWAYS bigger than the distance you're 
looking for - dont
worry if its too big.  Its in the units of your source data (ie. degrees)

   distance(transform(the_geom,<SRID>), transform(<GEOMETRY>,<SRID>) < 
 <DIST>
This calculates the distance between the centre of your circle 
("<GEOMETRY>") and the
data in the table.  
   transform( <geom>, <srid>)  - re-projects <geom> to the projection 
defined by <srid>

The expand() expression isnt really needed, but it will significantly 
speed this process up since you'll only be doing the expensive functions 
on a small portion of the data.

If you dont want to use transform() you can use the distance-on ellipse 
equations:
distance_spheroid(point, point, spheroid).  There's documentation on the 
postgis site (postgis.org) - but these function may not give the correct 
results (their accuracy has not been tested).


Hope this makes sense,
dave



>  
>