Calculate Earth global distances
I wrote this sp to be able to calculate the distance between two locations on the earth if the latitude and longitude is known for each. Could be useful for GPS work or maps.
Usage:
If I know the following city locations:
Sydney 151.2 E 33.87 S
Canberra 149.15 E 35.24 S
DECLARE @RC int
DECLARE @Dist decimal(28,20)
EXEC @RC = [dbo].[sp_EarthDistance] 149.15, -35.24, 151.2, -33.87, @Dist OUTPUT
Select @RC, 'Distance from Sydney and Canberra = ', @Dist
Result (kilometers):
0 Distance between Sydney and Canberra = 242.04260703665335000000
You can change the precision or distance units yourself.
SET QUOTED_IDENTIFIER ON
GO
SET ANSI_NULLS ON
GO
if exists (select * from dbo.sysobjects where id = object_id(N'[dbo].[sp_EarthDistance]') and OBJECTPROPERTY(id, N'IsProcedure') = 1)
drop procedure [dbo].[sp_EarthDistance]
GO
create proc sp_EarthDistance
@lng1 Decimal(18,6) = NULL
, @lat1 Decimal(18,6) = NULL
, @lng2 Decimal(18,6) = NULL
, @lat2 Decimal(18,6) = NULL
, @Dist Decimal(28,20) OUTPUT
/********1*********2*********3*********4*********5*********6*********8*********9*********0*********1*********2*****
**
** $Header: $
**
*******************************************************************************************************************
** $Revision: $
** $Author: $
** $History: $
**
*******************************************************************************************************************
**
** Name: sp_EarthDistance
** Desc: Stored procedure sp_EarthDistance calculates the Great Circle distance between two points on
** the earth, given the longitude and latitude for each in degrees.
** Usage:
** EXEC sp_EarthDistance @lng1, @lat2, @lng2, @lat2, @Dist OUTPUT
**
** Note that negative values represent West and South
**
** Return values:
** 4 - Error Calculating Distance
** 3 - Root Sum Square of delta coordinates Calculation Error
** 2 - Cartesian Coordinate Calculation Error
** 1 - Invalid Coordinates
** 0 - Completed OK
*******************************************************************************************************************
** Change History - All Author comments below this point.
*******************************************************************************************************************
** Author Date Description
** ------- -------- -------------------------------------------
** Neil Jacobson 10-Jan-2003 Original Issue
******************************************************************************************************************/
as
Declare @err Int
, @radlng1 Decimal(28,20)
, @radlng2 Decimal(28,20)
, @radcolat1 Decimal(28,20)
, @radcolat2 Decimal(28,20)
, @EarthRadius Decimal(28,20)
, @x1 Decimal(28,20)
, @y1 Decimal(28,20)
, @z1 Decimal(28,20)
, @x2 Decimal(28,20)
, @y2 Decimal(28,20)
, @z2 Decimal(28,20)
, @alpha Decimal(28,20)
/*
** Step 1 - Convert Latitude and Longitude
** to radians and get Co-Latitude
** Latitude is the angle from the Equator to the position
** Co-Latitude is the angle from the pole to the position
*/
If (@lng1 is null)
Begin
Return 1
End
Else
Begin
Select @radlng1 = @lng1 * PI()/180.0
End
If (@lat1 is null)
Begin
Return 1
End
Else
Begin
Select @radcolat1 = (PI()/2) - (@lat1 * PI()/180.0)
End
If (@lng2 is null)
Begin
Return 1
End
Else
Begin
Select @radlng2 = @lng2 * PI()/180.0
End
If (@lat2 is null)
Begin
Return 1
End
Else
Begin
Select @radcolat2 = (PI()/2) - (@lat2 * PI()/180.0)
End
Select @EarthRadius = 6378.0000 -- Kilometres, Change this to convert units.
Select @Err = 0
/*
** Step 2 Get the Cartesian Coordinates
** x = radius Cos(T) Sin(F)
** y = radius Sin(T) Sin(F)
** z = radius Cos(F)
** T (Theta) is the Longitude angle in radians
** F (Phi) is the Co-Latitude angle in radians (90 degrees - Latitude)
*/
Select @x1 = @EarthRadius * COS(@radlng1) * SIN(@radcolat1)
, @y1 = @EarthRadius * SIN(@radlng1) * SIN(@radcolat1)
, @z1 = @EarthRadius * COS(@radcolat1)
, @x2 = @EarthRadius * COS(@radlng2) * SIN(@radcolat2)
, @y2 = @EarthRadius * SIN(@radlng2) * SIN(@radcolat2)
, @z2 = @EarthRadius * COS(@radcolat2)
Select @Err = @@Error
If @Err <> 0
Begin
return 2
End
/*
** Step 3 Calculate the square root of the sum
** of the deltas of x, y and z
** ______________________
** a = \/ dx.dx + dy.dy + dz.dz
**
** where dx = x2 - x1, dy = y2 - y1, dz = z2 - z1
*/
Select @alpha = SQRT(((@x2-@x1)*(@x2-@x1))+((@y2-@y1)*(@y2-@y1))+((@z2-@z1)*(@z2-@z1)))
Select @Err = @@Error
If @Err <> 0
Begin
return 3
End
/*
** Step 4 Calculate the distance with the law of cosines
**
** distance = r * ArcCos((r*r + r*r - a*a)/(2*r*r))
**
** r = earth radius
**
*/
Select @Dist = @EarthRadius * (ACOS(( (@EarthRadius*@EarthRadius) + (@EarthRadius*@EarthRadius) - (@alpha*@alpha) )/(2*@EarthRadius*@EarthRadius)) )
Select @Err = @@Error
If @Err <> 0
Begin
return 4
End
Return 0
GO
SET QUOTED_IDENTIFIER OFF
GO
SET ANSI_NULLS ON
GO
GRANT EXECUTE ON [dbo].[sp_EarthDistance] TO [public]
GO
