SQL Server Forums
Profile | Register | Active Topics | Members | Search | Forum FAQ
 
Register Now and get your question answered!
Username:
Password:
Save Password
Forgot your Password?

 All Forums
 General SQL Server Forums
 Script Library
 Great Circle Distance Function - Haversine Formula
 New Topic  Reply to Topic
 Printer Friendly
Next Page
Author Previous Topic Topic Next Topic
Page: of 3

Michael Valentine Jones
Yak DBA Kernel (pronounced Colonel)

USA
7020 Posts

Posted - 03/28/2007 :  22:43:31  Show Profile  Reply with Quote
This function computes the great circle distance in Kilometers using the Haversine formula distance calculation.

If you want it in miles, change the average radius of Earth to miles in the function.

create function dbo.F_GREAT_CIRCLE_DISTANCE
	(
	@Latitude1  float,
	@Longitude1 float,
	@Latitude2  float,
	@Longitude2 float
	)
returns float
as
/*
fUNCTION: F_GREAT_CIRCLE_DISTANCE

	Computes the Great Circle distance in kilometers
	between two points on the Earth using the
	Haversine formula distance calculation.

Input Parameters:
	@Longitude1 - Longitude in degrees of point 1
	@Latitude1  - Latitude  in degrees of point 1
	@Longitude2 - Longitude in degrees of point 2
	@Latitude2  - Latitude  in degrees of point 2

*/
begin
declare @radius float

declare @lon1  float
declare @lon2  float
declare @lat1  float
declare @lat2  float

declare @a float
declare @distance float

-- Sets average radius of Earth in Kilometers
set @radius = 6371.0E

-- Convert degrees to radians
set @lon1 = radians( @Longitude1 )
set @lon2 = radians( @Longitude2 )
set @lat1 = radians( @Latitude1 )
set @lat2 = radians( @Latitude2 )

set @a = sqrt(square(sin((@lat2-@lat1)/2.0E)) + 
	(cos(@lat1) * cos(@lat2) * square(sin((@lon2-@lon1)/2.0E))) )

set @distance =
	@radius * ( 2.0E *asin(case when 1.0E < @a then 1.0E else @a end ))

return @distance

end


Edit: corrected spelling


CODO ERGO SUM

Edited by - Michael Valentine Jones on 04/01/2007 01:33:57

SwePeso
Patron Saint of Lost Yaks

Sweden
30265 Posts

Posted - 03/29/2007 :  08:44:40  Show Profile  Visit SwePeso's Homepage  Reply with Quote
Well done Michael!

Here is some more information about Haversine formula and also the average radius in miles (3956)
http://www.movable-type.co.uk/scripts/GIS-FAQ-5.1.html

It's interesting to see that the old Pythagorean Theorem is only a few percent wrong
quote:
less than 30 meters (100 ft) for latitudes less than 70 degrees
less than 20 meters ( 66 ft) for latitudes less than 50 degrees
less than 9 meters ( 30 ft) for latitudes less than 30 degrees


Haversine formula can error as much as 2 km (1 mi). But only under the context "half around the world", 20000 km / 12000 miles.


Peter Larsson
Helsingborg, Sweden

Edited by - SwePeso on 03/29/2007 08:52:27
Go to Top of Page

SwePeso
Patron Saint of Lost Yaks

Sweden
30265 Posts

Posted - 03/29/2007 :  08:49:15  Show Profile  Visit SwePeso's Homepage  Reply with Quote
Just for the curiosity, here is Vincenty's formula (which is accurate to the millimeter)
http://www.movable-type.co.uk/scripts/LatLongVincenty.html

Haversine is much faster and accurate enough
quote:
Vincentyfs formula is accurate to within 0.5mm, or 0.000015 (!), on the ellipsoid being used. Calculations based on a spherical model, such as the (much simpler) Haversine, are accurate to around 0.3% (which is still good enough for most purposes, of course).



Peter Larsson
Helsingborg, Sweden

Edited by - SwePeso on 03/29/2007 08:57:30
Go to Top of Page

Michael Valentine Jones
Yak DBA Kernel (pronounced Colonel)

USA
7020 Posts

Posted - 03/29/2007 :  18:15:55  Show Profile  Reply with Quote
quote:
Originally posted by Peso

Just for the curiosity, here is Vincenty's formula (which is accurate to the millimeter)
http://www.movable-type.co.uk/scripts/LatLongVincenty.html

Haversine is much faster and accurate enough
quote:
Vincentyfs formula is accurate to within 0.5mm, or 0.000015 (!), on the ellipsoid being used. Calculations based on a spherical model, such as the (much simpler) Haversine, are accurate to around 0.3% (which is still good enough for most purposes, of course).



Peter Larsson
Helsingborg, Sweden



I looked into that before, and my head started hurting when I saw the formula.

In any case, there are so many other factors that introduce distance inaccuracy, such as difference in altitude between locations, that it seemed unlikely that it would even produce a more accurate result. Basically, the input data is likely to have so much inaccuracy built into it that anything more complex than the Haversine distance calculation is unlikely to produce better results. Garbage in, garbage out.


When I get to it, I want to post a follow-up on the method of calculating a square around the circle in order to limit the number of locations to be searched with an index lookup on longitude and latitude. For database applications, this is likely to have the most benefit. Basically you calculate a maximum and minimum longitude and latitude, and use that to limit the locations to be searched. Using this method, you eliminate most locations because they are outside the square, and use the Haversine function on the remaining locations. By chance, about 75% on the locations inside the square will be inside the search circle.

I did some work on this before, and found that calculating the east/west limits was much more complex than it may seem at first glance. The kilometers per degree along the latitude circle vary according to latitude. There is also an error that increases the further you are from the equator, because the distance east/west is shorter along the great circle distance, than just calculating the distance along the latitude circle. The magnitude of this error also increases as the size of the search circle increases. This means that you can leave out locations by making the square too small. More complications happen if your search circle is so large that one of the earth’s poles is inside the circle, or if the circle crosses the International Date Line.





CODO ERGO SUM

Edited by - Michael Valentine Jones on 03/29/2007 20:14:35
Go to Top of Page

SwePeso
Patron Saint of Lost Yaks

Sweden
30265 Posts

Posted - 03/29/2007 :  18:59:20  Show Profile  Visit SwePeso's Homepage  Reply with Quote
I did a test measuring the distance between the two cities Stockholm and Kiruna in Sweden.

The Haversine formula said the distance is 954.2 kilometers (593 miles).
Pythagorean Theorem said the distance is 957.7 kilometers (595 miles).

The error is only 0.35% between these two methods!
A mere 3.5 kilometer (2 mile) stray when travelling more than 950 kilometers (590 miles)!


Peter Larsson
Helsingborg, Sweden
Go to Top of Page

Michael Valentine Jones
Yak DBA Kernel (pronounced Colonel)

USA
7020 Posts

Posted - 03/29/2007 :  20:38:36  Show Profile  Reply with Quote
Well, what I really want is a simple, accurate method to compute the min/max longitude to get the limits for the query. For performance, I think it is far more important to limit the number of points for which you do a computation, than to use a less CPU intensive calculation.

Instead of "How far away is this point?", I would like to find, "What are the longitudes that are exactly X number of distance units east or west of a given latitude/longitude using the great circle distance?"

I found the following formula in Wikipedia that may be what I need. Unfortunately, they failed to post a TSQL implementation, so I guess I have to do that myself.

http://en.wikipedia.org/wiki/Longitude
"As opposed to a degree of latitude, which always corresponds almost exactly to sixty nautical miles or about 111 km (69 statute miles, each of 5280 feet), a degree of longitude corresponds to a distance that varies from 0 to 111 km: it is 111 km times the cosine of the latitude, when the distance is laid out on a circle of constant latitude; if the shortest distance, on a great circle were used, the distance would be even a little less. More precisely, one degree of longitude = (111.320 + 0.373sin²ö)cosö km, where ö is latitude)."








CODO ERGO SUM
Go to Top of Page

jezemine
Flowing Fount of Yak Knowledge

USA
2886 Posts

Posted - 03/31/2007 :  16:13:37  Show Profile  Visit jezemine's Homepage  Reply with Quote
minor typo in your function: @Logitude ought to be @Longitude.

that is unless you did it that way on purpose so it would line up nicely with @Latitude

also I agree Haversine will be plenty accurate enough for most all applications. You might notice a small discrepancy going from the equator north 90 degrees, as opposed to east or west 90 degrees. but even then the difference would small.

if the earth's period were shorter, like spinning at 100Hz or so, there would be more of a difference because it would be much more ellipsoidal. of course then we'd be too dizzy to care...



www.elsasoft.org

Edited by - jezemine on 03/31/2007 16:38:56
Go to Top of Page

Kristen
Test

United Kingdom
22415 Posts

Posted - 04/02/2007 :  03:32:04  Show Profile  Reply with Quote
I've got these two links in my file of "Ready to add to FAQ", included here for completeness. The discussion is mostly about Accuracy v. Speed

http://www.sqlteam.com/forums/topic.asp?TOPIC_ID=40295
http://www.sqlteam.com/forums/topic.asp?TOPIC_ID=57242

Kristen
Go to Top of Page

7695af9b
Starting Member

1 Posts

Posted - 10/24/2007 :  04:03:28  Show Profile  Reply with Quote
What do you mean?
Go to Top of Page

AndrewMurphy
Flowing Fount of Yak Knowledge

Ireland
2916 Posts

Posted - 10/24/2007 :  07:26:38  Show Profile  Reply with Quote
"What do you mean?"....This must qualify as one of the more obscure questions ever posted here. To whom is this latest question posed? And to what are you referring?
Go to Top of Page

Kristen
Test

United Kingdom
22415 Posts

Posted - 10/24/2007 :  07:50:53  Show Profile  Reply with Quote
not to mention the 6.5 elapsed months
Go to Top of Page

Kristen
Test

United Kingdom
22415 Posts

Posted - 11/02/2007 :  09:15:08  Show Profile  Reply with Quote
"( 2.0E *asin(case when 1.0E < @a then 1.0E else @a end ))"

I'm wondering if this is a robust alternative to

2 * atan2(sqrt(a), sqrt(1-a))

because in the tests I've been doing this seems to be where I'm getting some error creeping in (i.e. comparing SQL's calculation above to a hand-calculation using atan2.

Obviously SQL doesn't have atan2(x, y) ...

Kristen
Go to Top of Page

Zoroaster
Aged Yak Warrior

USA
702 Posts

Posted - 11/02/2007 :  10:08:15  Show Profile  Reply with Quote
quote:
Originally posted by 7695af9b

What do you mean?



42



Future guru in the making.
Go to Top of Page

SwePeso
Patron Saint of Lost Yaks

Sweden
30265 Posts

Posted - 11/06/2007 :  04:42:01  Show Profile  Visit SwePeso's Homepage  Reply with Quote
quote:
Originally posted by Kristen

Obviously SQL doesn't have atan2(x, y) ...

http://msdn2.microsoft.com/en-us/library/ms173854.aspx

ATN2



E 12°55'05.25"
N 56°04'39.16"
Go to Top of Page

Kristen
Test

United Kingdom
22415 Posts

Posted - 11/06/2007 :  10:32:05  Show Profile  Reply with Quote
Would you Adam and Eve it? Its even in SQL2000 ... I wonder how I missed that :(

I'll retro fit that and see if it helps my accuracy problem. Thanks.

Kristen
Go to Top of Page

Van
Constraint Violating Yak Guru

458 Posts

Posted - 11/06/2007 :  10:50:18  Show Profile  Reply with Quote
In other news .9999999 repeating = 1
Go to Top of Page

SwePeso
Patron Saint of Lost Yaks

Sweden
30265 Posts

Posted - 11/06/2007 :  11:45:27  Show Profile  Visit SwePeso's Homepage  Reply with Quote
Try these constants for earth radius

km - 6366.70701949371
mi - 3956.0883313286096695299



E 12°55'05.25"
N 56°04'39.16"
Go to Top of Page

jezemine
Flowing Fount of Yak Knowledge

USA
2886 Posts

Posted - 11/06/2007 :  15:17:39  Show Profile  Visit jezemine's Homepage  Reply with Quote
quote:
Originally posted by Van

In other news .9999999 repeating = 1



http://qntm.org/pointnine


elsasoft.org
Go to Top of Page

jezemine
Flowing Fount of Yak Knowledge

USA
2886 Posts

Posted - 11/06/2007 :  15:22:16  Show Profile  Visit jezemine's Homepage  Reply with Quote
quote:
Originally posted by Peso

Try these constants for earth radius

km - 6366.70701949371
mi - 3956.0883313286096695299



E 12°55'05.25"
N 56°04'39.16"




those are some very accurate measurements! the last one has the radius of the earth to within a millionth of the size of a hydrogen atom!




elsasoft.org
Go to Top of Page

Kristen
Test

United Kingdom
22415 Posts

Posted - 11/06/2007 :  17:17:16  Show Profile  Reply with Quote
Yeah, but by the time you've sum'd the squares of the other two sides to get the hypotenuse you're up to nearly a 100th of an electron ...
Go to Top of Page

Michael Valentine Jones
Yak DBA Kernel (pronounced Colonel)

USA
7020 Posts

Posted - 11/06/2007 :  17:36:34  Show Profile  Reply with Quote
quote:
Originally posted by jezemine

quote:
Originally posted by Peso

Try these constants for earth radius

km - 6366.70701949371
mi - 3956.0883313286096695299



E 12°55'05.25"
N 56°04'39.16"




those are some very accurate measurements! the last one has the radius of the earth to within a millionth of the size of a hydrogen atom!




elsasoft.org




Yes, but that's gigantic compared to the size of the strings a hydrogen atom is (allegedly) made from.

CODO ERGO SUM
Go to Top of Page
Page: of 3 Previous Topic Topic Next Topic  
Next Page
 New Topic  Reply to Topic
 Printer Friendly
Jump To:
SQL Server Forums © 2000-2009 SQLTeam Publishing, LLC Go To Top Of Page
This page was generated in 0.21 seconds. Powered By: Snitz Forums 2000