Discussion:
SRID for analyzing a USA national data set in Meters
mfrumin
2007-04-02 13:11:55 UTC
Permalink
I have a data set that includes points all over the United States in lat/lng
(say, SRID = 4326) but I need to do some analysis on these data in real
units, say meters. So, my task is to select an appropriate SRID that will
work reasonably well for the whole USA to reproject lat/lng into meters.
There are a bazillion SRID's available. Any suggestions?

This analysis doesn't have to be super-precise, so I'm not so concerned with
the limitations of a national projection (as opposed to state/region).

thanks,
Michael
--
View this message in context: http://www.nabble.com/SRID-for-analyzing-a-USA-national-data-set-in-Meters-tf3505664.html#a9790471
Sent from the PostGIS - User mailing list archive at Nabble.com.
Pedro Doria Meunier
2007-04-02 13:29:42 UTC
Permalink
Hello Michael,

You determine the correct SRID with this bit of code (PHP):
(the $loc was obtained from a SELECT with the astext(geometry) function)

The idea behind this is summarized in the UtmZoneFromLong() function.
Pick a Long value and calculate the utm zone for it.
Then use the SridFromUtmZone() function to get the srid from the given mask
(see below) --> this queries spatial_ref_sys with the calculated $utm param.
SELECT srid FROM spatial_ref_sys WHERE srtext LIKE 'PROJCS[\"WGS 84 / UTM
zone $utm%' LIMIT 1

// extract the POINT geometry
$geometry = str_replace("POINT(","", $loc);
$geometry = str_replace(")","",$geometry);
$coordinate = explode(" ",$geometry);
$lon=$coordinate[0];
$lat=$coordinate[1];

// get utm zone from coordinates
$utm=UtmZoneFromLong($lon);
// get srid from the spatial_ref_sys table for the given utm zone
// add N or S
if($lat>=0) $utm.="N";
else $utm.="S";
$srid=SridFromUtmZone($utm, $connection);


// determine UTM zone from LonLat coordinates
// params $long: longitude
// returns (int) utm
function UtmZoneFromLong($lon) {
$utm = floor(($lon + 180) / 6) + 1;
return $utm;
}

// query the spatial_ref_sys table for the srid of the given utm zone
// params $utm: utm zone string (MUST BE TRAILED WITH 'N' OR 'S'!), $conn:
connection to the db
// returns (int) srid
function SridFromUtmZone($utm, $conn){
$myresult = pg_exec($conn, "SELECT srid FROM spatial_ref_sys WHERE srtext
LIKE 'PROJCS[\"WGS 84 / UTM zone $utm%' LIMIT 1");
$srid=pg_result($myresult, 0, 0);
return $srid;
}

HTH,
Pedro Doria Meunier

-----Original Message-----
From: postgis-users-***@postgis.refractions.net
[mailto:postgis-users-***@postgis.refractions.net] On Behalf Of mfrumin
Sent: segunda-feira, 2 de Abril de 2007 14:12
To: postgis-***@postgis.refractions.net
Subject: [postgis-users] SRID for analyzing a USA national data set in
Meters


I have a data set that includes points all over the United States in lat/lng
(say, SRID = 4326) but I need to do some analysis on these data in real
units, say meters. So, my task is to select an appropriate SRID that will
work reasonably well for the whole USA to reproject lat/lng into meters.
There are a bazillion SRID's available. Any suggestions?

This analysis doesn't have to be super-precise, so I'm not so concerned with
the limitations of a national projection (as opposed to state/region).

thanks,
Michael
--
View this message in context:
http://www.nabble.com/SRID-for-analyzing-a-USA-national-data-set-in-Meters-t
f3505664.html#a9790471
Sent from the PostGIS - User mailing list archive at Nabble.com.
Michael Frumin
2007-04-02 15:13:02 UTC
Permalink
pedro,

This appears to suggest that i should pick an SRID based on one point in
my data set. However, since the data set is national, it includes
points from all around the USA (i.e. at least 4 UTM zones). So, if I
pick the SRID for, say, Mountain Time, will that horribly distort points
on the east and/or west coast? is there no SRID that is perhaps less
accurate in any one time zone, but more accurate across all time zones
in the [continental] US?

thanks,
mike
Post by Pedro Doria Meunier
Hello Michael,
(the $loc was obtained from a SELECT with the astext(geometry) function)
The idea behind this is summarized in the UtmZoneFromLong() function.
Pick a Long value and calculate the utm zone for it.
Then use the SridFromUtmZone() function to get the srid from the given mask
(see below) --> this queries spatial_ref_sys with the calculated $utm param.
SELECT srid FROM spatial_ref_sys WHERE srtext LIKE 'PROJCS[\"WGS 84 / UTM
zone $utm%' LIMIT 1
// extract the POINT geometry
$geometry = str_replace("POINT(","", $loc);
$geometry = str_replace(")","",$geometry);
$coordinate = explode(" ",$geometry);
$lon=$coordinate[0];
$lat=$coordinate[1];
// get utm zone from coordinates
$utm=UtmZoneFromLong($lon);
// get srid from the spatial_ref_sys table for the given utm zone
// add N or S
if($lat>=0) $utm.="N";
else $utm.="S";
$srid=SridFromUtmZone($utm, $connection);
// determine UTM zone from LonLat coordinates
// params $long: longitude
// returns (int) utm
function UtmZoneFromLong($lon) {
$utm = floor(($lon + 180) / 6) + 1;
return $utm;
}
// query the spatial_ref_sys table for the srid of the given utm zone
connection to the db
// returns (int) srid
function SridFromUtmZone($utm, $conn){
$myresult = pg_exec($conn, "SELECT srid FROM spatial_ref_sys WHERE srtext
LIKE 'PROJCS[\"WGS 84 / UTM zone $utm%' LIMIT 1");
$srid=pg_result($myresult, 0, 0);
return $srid;
}
HTH,
Pedro Doria Meunier
-----Original Message-----
Sent: segunda-feira, 2 de Abril de 2007 14:12
Subject: [postgis-users] SRID for analyzing a USA national data set in
Meters
I have a data set that includes points all over the United States in lat/lng
(say, SRID = 4326) but I need to do some analysis on these data in real
units, say meters. So, my task is to select an appropriate SRID that will
work reasonably well for the whole USA to reproject lat/lng into meters. There are a bazillion SRID's available. Any suggestions?
This analysis doesn't have to be super-precise, so I'm not so concerned with
the limitations of a national projection (as opposed to state/region).
thanks,
Michael
http://www.nabble.com/SRID-for-analyzing-a-USA-national-data-set-in-Meters-t
f3505664.html#a9790471
Sent from the PostGIS - User mailing list archive at Nabble.com.
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
Stephen Woodbridge
2007-04-02 15:25:44 UTC
Permalink
Michael,

The Simple answer is that there is no simple answer. You have not given
us enough information about what kind of calculations you need to do and
what are acceptable errors. Are you worried about distortions in
distance, area, or direction/heading for example.

Why don't you do all you analysis in lat/lon and convert to meters based
on the fact that 1 deg is roughly 69..... nautical miles and you can
convert miles to meters? You this would be problematic over large
changes in latitude. Or use distance_sphere() which uses another
approximation based on a spherical earth.

Any suggestion we can offer will have problems. You need to look at the
problems in light of you analysis constraints and pick a compomise that
you can live with. If you give us additional information me might be
able to make additional suggestions.

Also variants on this request can be found in the archives as this is a
fairly common question.

Hope this helps,
-Steve W
Post by Michael Frumin
pedro,
This appears to suggest that i should pick an SRID based on one point in
my data set. However, since the data set is national, it includes
points from all around the USA (i.e. at least 4 UTM zones). So, if I
pick the SRID for, say, Mountain Time, will that horribly distort points
on the east and/or west coast? is there no SRID that is perhaps less
accurate in any one time zone, but more accurate across all time zones
in the [continental] US?
thanks,
mike
Post by Pedro Doria Meunier
Hello Michael,
(the $loc was obtained from a SELECT with the astext(geometry) function)
The idea behind this is summarized in the UtmZoneFromLong() function.
Pick a Long value and calculate the utm zone for it.
Then use the SridFromUtmZone() function to get the srid from the given mask
(see below) --> this queries spatial_ref_sys with the calculated $utm param.
SELECT srid FROM spatial_ref_sys WHERE srtext LIKE 'PROJCS[\"WGS 84 / UTM
zone $utm%' LIMIT 1
// extract the POINT geometry
$geometry = str_replace("POINT(","", $loc);
$geometry = str_replace(")","",$geometry);
$coordinate = explode(" ",$geometry);
$lon=$coordinate[0];
$lat=$coordinate[1];
// get utm zone from coordinates
$utm=UtmZoneFromLong($lon);
// get srid from the spatial_ref_sys table for the given utm zone
// add N or S
if($lat>=0) $utm.="N";
else $utm.="S";
$srid=SridFromUtmZone($utm, $connection);
// determine UTM zone from LonLat coordinates
// params $long: longitude
// returns (int) utm
function UtmZoneFromLong($lon) {
$utm = floor(($lon + 180) / 6) + 1;
return $utm;
}
// query the spatial_ref_sys table for the srid of the given utm zone
connection to the db
// returns (int) srid
function SridFromUtmZone($utm, $conn){
$myresult = pg_exec($conn, "SELECT srid FROM spatial_ref_sys WHERE srtext
LIKE 'PROJCS[\"WGS 84 / UTM zone $utm%' LIMIT 1");
$srid=pg_result($myresult, 0, 0);
return $srid;
}
HTH,
Pedro Doria Meunier
-----Original Message-----
Sent: segunda-feira, 2 de Abril de 2007 14:12
Subject: [postgis-users] SRID for analyzing a USA national data set in
Meters
I have a data set that includes points all over the United States in lat/lng
(say, SRID = 4326) but I need to do some analysis on these data in real
units, say meters. So, my task is to select an appropriate SRID that will
work reasonably well for the whole USA to reproject lat/lng into
meters. There are a bazillion SRID's available. Any suggestions?
This analysis doesn't have to be super-precise, so I'm not so
concerned with
the limitations of a national projection (as opposed to state/region).
thanks,
Michael
http://www.nabble.com/SRID-for-analyzing-a-USA-national-data-set-in-Meters-t
f3505664.html#a9790471
Sent from the PostGIS - User mailing list archive at Nabble.com.
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
Pedro Doria Meunier
2007-04-02 15:30:18 UTC
Permalink
Hey Mike,

The only global system, adopted by 1984, is WGS84 (this has been revised in
2004 and it's valid until 2010).
This, as you well know, is a system whose units are expressed in degrees.
If you wish to work in metres you must use a 'localized' system (projected).
The correct geoid must be used for accurate calculations.

So my advice to you is defining 'work areas' (as in regions) if you really
need to work with metres..

This comes from my notes:
"an inaccurate way for calculating the distance would be:
distance [m] = 6378137.0 [m] * Pi * distance [degree] / 180.0
6378137.0 = earth radius

accurate measurement involves the correct geoid to be used (UTM projection)"

INACURRACY ALARM! The above method is horribly inaccurate...

Depending on your task you can programmatically (as in using PHP or whatever
lang) process your points.
Could you elaborate on the task at hand?

Pedro.


-----Original Message-----
From: postgis-users-***@postgis.refractions.net
[mailto:postgis-users-***@postgis.refractions.net] On Behalf Of Michael
Frumin
Sent: segunda-feira, 2 de Abril de 2007 16:13
To: postgis-***@postgis.refractions.net
Subject: Re: [postgis-users] SRID for analyzing a USA national data set in
Meters

pedro,

This appears to suggest that i should pick an SRID based on one point in
my data set. However, since the data set is national, it includes
points from all around the USA (i.e. at least 4 UTM zones). So, if I
pick the SRID for, say, Mountain Time, will that horribly distort points
on the east and/or west coast? is there no SRID that is perhaps less
accurate in any one time zone, but more accurate across all time zones
in the [continental] US?

thanks,
mike
Paul Ramsey
2007-04-02 15:35:19 UTC
Permalink
Post by Pedro Doria Meunier
Hey Mike,
The only global system, adopted by 1984, is WGS84 (this has been revised in
2004 and it's valid until 2010).
This, as you well know, is a system whose units are expressed in degrees.
If you wish to work in metres you must use a 'localized' system (projected).
The correct geoid must be used for accurate calculations.
So my advice to you is defining 'work areas' (as in regions) if you really
need to work with metres..
distance [m] = 6378137.0 [m] * Pi * distance [degree] / 180.0
6378137.0 = earth radius
Use distance_sphere() to get this as a built-in postgis function.
Distance_spheroid() returns an accurate distance, but at 10 times the
CPU cost of the sphere calc.

P
Post by Pedro Doria Meunier
accurate measurement involves the correct geoid to be used (UTM projection)"
INACURRACY ALARM! The above method is horribly inaccurate...
Depending on your task you can programmatically (as in using PHP or whatever
lang) process your points.
Could you elaborate on the task at hand?
Pedro.
-----Original Message-----
Frumin
Sent: segunda-feira, 2 de Abril de 2007 16:13
Subject: Re: [postgis-users] SRID for analyzing a USA national data set in
Meters
pedro,
This appears to suggest that i should pick an SRID based on one point in
my data set. However, since the data set is national, it includes
points from all around the USA (i.e. at least 4 UTM zones). So, if I
pick the SRID for, say, Mountain Time, will that horribly distort points
on the east and/or west coast? is there no SRID that is perhaps less
accurate in any one time zone, but more accurate across all time zones
in the [continental] US?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
Paul Ramsey
Refractions Research
http://www.refractions.net
***@refractions.net
Phone: 250-383-3022
Cell: 250-885-0632
Michael Frumin
2007-04-02 15:39:09 UTC
Permalink
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?

As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the
map these sets, with a 25 mile buffer around the first set, using
Geoserver. So, distance and area are both somewhat important, heading
not at all. distance_sphere(oid), sounds good for the calculation, but
won't help with the buffering because it doesn't tell me the 'distance'
in lat/lng space that would equate to 25 miles (because of course this
varies over the globe). To achieve this I would need to reproject into
something that is in meters, and buffer around that.

How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?

thanks,
mike
Post by Pedro Doria Meunier
Michael,
The Simple answer is that there is no simple answer. You have not given us enough information about what kind of calculations you need to do and what are acceptable errors. Are you worried about distortions in distance, area, or direction/heading for example.
Why don't you do all you analysis in lat/lon and convert to meters based on the fact that 1 deg is roughly 69..... nautical miles and you can convert miles to meters? You this would be problematic over large changes in latitude. Or use distance_sphere() which uses another approximation based on a spherical earth.
Any suggestion we can offer will have problems. You need to look at the problems in light of you analysis constraints and pick a compomise that you can live with. If you give us additional information me might be able to make additional suggestions.
Also variants on this request can be found in the archives as this is a fairly common question.
Hope this helps,
-Steve W
Post by Michael Frumin
pedro,
This appears to suggest that i should pick an SRID based on one point in > my data set. However, since the data set is national, it includes > points from all around the USA (i.e. at least 4 UTM zones). So, if I > pick the SRID for, say, Mountain Time, will that horribly distort points > on the east and/or west coast? is there no SRID that is perhaps less > accurate in any one time zone, but more accurate across all time zones > in the [continental] US?
thanks,
mike
Hello Michael,
(the $loc was obtained from a SELECT with the astext(geometry) function)
The idea behind this is summarized in the UtmZoneFromLong() function.
Pick a Long value and calculate the utm zone for it.
Then use the SridFromUtmZone() function to get the srid from the given >> mask
(see below) --> this queries spatial_ref_sys with the calculated $utm >> param.
SELECT srid FROM spatial_ref_sys WHERE srtext LIKE 'PROJCS[\"WGS 84 / UTM
zone $utm%' LIMIT 1
// extract the POINT geometry
$geometry = str_replace("POINT(","", $loc);
$geometry = str_replace(")","",$geometry);
$coordinate = explode(" ",$geometry);
$lon=$coordinate[0];
$lat=$coordinate[1];
// get utm zone from coordinates
$utm=UtmZoneFromLong($lon);
// get srid from the spatial_ref_sys table for the given utm zone
// add N or S
if($lat>=0) $utm.="N";
else $utm.="S";
$srid=SridFromUtmZone($utm, $connection);
// determine UTM zone from LonLat coordinates
// params $long: longitude
// returns (int) utm
function UtmZoneFromLong($lon) {
$utm = floor(($lon + 180) / 6) + 1;
return $utm;
}
// query the spatial_ref_sys table for the srid of the given utm zone
connection to the db
// returns (int) srid
function SridFromUtmZone($utm, $conn){
$myresult = pg_exec($conn, "SELECT srid FROM spatial_ref_sys WHERE srtext
LIKE 'PROJCS[\"WGS 84 / UTM zone $utm%' LIMIT 1");
$srid=pg_result($myresult, 0, 0);
return $srid;
}
HTH,
Pedro Doria Meunier
-----Original Message-----
Sent: segunda-feira, 2 de Abril de 2007 14:12
Subject: [postgis-users] SRID for analyzing a USA national data set in
Meters
I have a data set that includes points all over the United States in >> lat/lng
(say, SRID = 4326) but I need to do some analysis on these data in real
units, say meters. So, my task is to select an appropriate SRID that >> will
work reasonably well for the whole USA to reproject lat/lng into >> meters. There are a bazillion SRID's available. Any suggestions?
This analysis doesn't have to be super-precise, so I'm not so >> concerned with
the limitations of a national projection (as opposed to state/region).
thanks,
Michael
http://www.nabble.com/SRID-for-analyzing-a-USA-national-data-set-in-Meters-t >>
f3505664.html#a9790471
Sent from the PostGIS - User mailing list archive at Nabble.com.
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
Post by Pedro Doria Meunier
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
David William Bitner
2007-04-02 19:03:08 UTC
Permalink
It is a common mis-assumption that you need to do a buffer calculation. To
do this type of analysis, all you need is the distance_sphere() calculation:
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;

Just add a limit 1 to the end if all you want is the closest.

No buffer necessary at all.


David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
Jeremy Dunck
2007-04-02 20:03:26 UTC
Permalink
On 4/2/07, David William Bitner <***@gmail.com> wrote:
...
Post by David William Bitner
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order
What's the 1609 for?
David William Bitner
2007-04-02 20:14:00 UTC
Permalink
He wanted 25 miles, distance_sphere returns meters, it's a rough conversion
factor.
Post by Jeremy Dunck
...
Post by David William Bitner
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order
What's the 1609 for?
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
************************************
David William Bitner
Michael Frumin
2007-04-02 19:17:14 UTC
Permalink
correct. I wanted to calculate the buffer so that I can render a map
using GeoServer, which sits right on top of PostGIS, which would show
the first set of hospitals, the 25 mile radius, and the second set.

thanks,
mike
Post by David William Bitner
It is a common mis-assumption that you need to do a buffer calculation. To
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
------------------------------------------------------------------------
It is a common mis-assumption that you need to do a buffer
calculation. To do this type of analysis, all you need is the
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2
where distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I
was looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals,
and I just want to be able to identify for each hospital in the
first set, the hospitals in the second set within approximately 25
miles. I will the map these sets, with a 25 mile buffer around
the first set, using Geoserver. So, distance and area are both
somewhat important, heading not at all. distance_sphere(oid),
sounds good for the calculation, but won't help with the buffering
because it doesn't tell me the 'distance' in lat/lng space that
would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something
that is in meters, and buffer around that.
How egregious would you expect the errors to be if I simply use
the projection for the UTM zone that represents, say, Central time?
thanks,
mike
David William Bitner
2007-04-02 20:17:23 UTC
Permalink
If you're only concerned with 25 miles, you could use the method to
determine which UTM zone your data is in for your first hospital and then
use that for the calculation and buffer creation. You could convert Pedro's
PHP script into PL/PGSQL to make it easy to work with in the db.
Post by Michael Frumin
correct. I wanted to calculate the buffer so that I can render a map
using GeoServer, which sits right on top of PostGIS, which would show the
first set of hospitals, the 25 mile radius, and the second set.
thanks,
mike
It is a common mis-assumption that you need to do a buffer calculation. To
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
------------------------------
It is a common mis-assumption that you need to do a buffer calculation.
To do this type of analysis, all you need is the distance_sphere()
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2
where distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
************************************
David William Bitner
David William Bitner
2007-04-02 20:42:44 UTC
Permalink
This function should work to return the correct UTM zone in WGS84 (unless I
messed something up quickly pulling this together -- verify before you use
it).

CREATE OR REPLACE FUNCTION utmzone(geometry)
RETURNS integer AS
$BODY$
declare
geomgeog geometry;
zone int;
pref int;
begin
geomgeog:=transform($1,4326);
if (y(geomgeog))>0 then
pref:=32600;
else
pref:=32700;
end if;
zone:=floor((x(geomgeog)+180)/6)+1;
return zone+pref;
end;
$BODY$
LANGUAGE 'plpgsql' VOLATILE;
Post by David William Bitner
If you're only concerned with 25 miles, you could use the method to
determine which UTM zone your data is in for your first hospital and then
use that for the calculation and buffer creation. You could convert Pedro's
PHP script into PL/PGSQL to make it easy to work with in the db.
Post by Michael Frumin
correct. I wanted to calculate the buffer so that I can render a map
using GeoServer, which sits right on top of PostGIS, which would show the
first set of hospitals, the 25 mile radius, and the second set.
thanks,
mike
It is a common mis-assumption that you need to do a buffer calculation. To
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
------------------------------
It is a common mis-assumption that you need to do a buffer calculation.
To do this type of analysis, all you need is the distance_sphere()
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2
where distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
************************************
David William Bitner
--
************************************
David William Bitner
Jeremy Dunck
2007-04-02 20:55:18 UTC
Permalink
Post by David William Bitner
CREATE OR REPLACE FUNCTION utmzone(geometry)
...
Post by David William Bitner
LANGUAGE 'plpgsql' VOLATILE;
This function is guaranteed to return the same result as long as
spatial_ref_sys doesn't change, and so could be classified as:
LANGUAGE 'plpgsql' IMMUTABLE;
for better performance.

Cheers,
Jeremy
Michael Frumin
2007-04-02 20:50:11 UTC
Permalink
this actualy returns the SRID, not the UTMZone; even better!

totally awesome, thanks.

-mike
Post by David William Bitner
This function should work to return the correct UTM zone in WGS84 (unless I
messed something up quickly pulling this together -- verify before you use
it).
CREATE OR REPLACE FUNCTION utmzone(geometry)
RETURNS integer AS
$BODY$
declare
geomgeog geometry;
zone int;
pref int;
begin
geomgeog:=transform($1,4326);
if (y(geomgeog))>0 then
pref:=32600;
else
pref:=32700;
end if;
zone:=floor((x(geomgeog)+180)/6)+1;
return zone+pref;
end;
$BODY$
LANGUAGE 'plpgsql' VOLATILE;
Post by David William Bitner
If you're only concerned with 25 miles, you could use the method to
determine which UTM zone your data is in for your first hospital and then
use that for the calculation and buffer creation. You could convert Pedro's
PHP script into PL/PGSQL to make it easy to work with in the db.
Post by Michael Frumin
correct. I wanted to calculate the buffer so that I can render a map
using GeoServer, which sits right on top of PostGIS, which would show the
first set of hospitals, the 25 mile radius, and the second set.
thanks,
mike
It is a common mis-assumption that you need to do a buffer calculation. To
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
------------------------------
It is a common mis-assumption that you need to do a buffer calculation.
To do this type of analysis, all you need is the distance_sphere()
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2
where distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
************************************
David William Bitner
-- ************************************
David William Bitner
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
------------------------------------------------------------------------
This function should work to return the correct UTM zone in WGS84
(unless I messed something up quickly pulling this together -- verify
before you use it).
CREATE OR REPLACE FUNCTION utmzone(geometry)
RETURNS integer AS
$BODY$
declare
geomgeog geometry;
zone int;
pref int;
begin
geomgeog:=transform($1,4326);
if (y(geomgeog))>0 then
pref:=32600;
else
pref:=32700;
end if;
zone:=floor((x(geomgeog)+180)/6)+1;
return zone+pref;
end;
$BODY$
LANGUAGE 'plpgsql' VOLATILE;
If you're only concerned with 25 miles, you could use the method
to determine which UTM zone your data is in for your first
hospital and then use that for the calculation and buffer
creation. You could convert Pedro's PHP script into PL/PGSQL to
make it easy to work with in the db.
correct. I wanted to calculate the buffer so that I can
render a map using GeoServer, which sits right on top of
PostGIS, which would show the first set of hospitals, the 25
mile radius, and the second set.
thanks,
mike
Post by David William Bitner
It is a common mis-assumption that you need to do a buffer calculation. To
select *,distance_sphere(hospital1geom,hospitals2geom) from hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order by
distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Post by Michael Frumin
Right, I should be more specific from the outset. I did do some
searching thru the PostGIS archives, and didn't find the answer I was
looking for; is there a PostGIS FAQ somewhere?
As for my problem, my inputs are two sets of geocoded hospitals, and I
just want to be able to identify for each hospital in the first set, the
hospitals in the second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set, using Geoserver.
So, distance and area are both somewhat important, heading not at all.
distance_sphere(oid), sounds good for the calculation, but won't help with
the buffering because it doesn't tell me the 'distance' in lat/lng space
that would equate to 25 miles (because of course this varies over the
globe). To achieve this I would need to reproject into something that is in
meters, and buffer around that.
How egregious would you expect the errors to be if I simply use the
projection for the UTM zone that represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
------------------------------------------------------------------------
It is a common mis-assumption that you need to do a buffer
calculation. To do this type of analysis, all you need is
select *,distance_sphere(hospital1geom,hospitals2geom) from
hospitals2 where
distance_sphere(hospital1geom,hospitals2geom)<25*1609 order
by distance_sphere(hospital1geom,hospitals2geom) asc;
Just add a limit 1 to the end if all you want is the closest.
No buffer necessary at all.
David
Right, I should be more specific from the outset. I did
do some searching thru the PostGIS archives, and didn't
find the answer I was looking for; is there a PostGIS FAQ
somewhere?
As for my problem, my inputs are two sets of geocoded
hospitals, and I just want to be able to identify for
each hospital in the first set, the hospitals in the
second set within approximately 25 miles. I will the map
these sets, with a 25 mile buffer around the first set,
using Geoserver. So, distance and area are both somewhat
important, heading not at all. distance_sphere(oid),
sounds good for the calculation, but won't help with the
buffering because it doesn't tell me the 'distance' in
lat/lng space that would equate to 25 miles (because of
course this varies over the globe). To achieve this I
would need to reproject into something that is in meters,
and buffer around that.
How egregious would you expect the errors to be if I
simply use the projection for the UTM zone that
represents, say, Central time?
thanks,
mike
_______________________________________________
postgis-users mailing list
http://postgis.refractions.net/mailman/listinfo/postgis-users
--
************************************
David William Bitner
--
************************************
David William Bitner
Michael Frumin
2007-04-03 02:51:07 UTC
Permalink
_______________________________________________
postgis-users mailing list
postgis-***@postgis.refractions.net
http://postgis.refractions.net/mailman/listinfo/postgis-users

Loading...