Technical Article

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        DateDescription
**  -------        ---------------------------------------------------
**  Neil Jacobson10-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

Rate

5 (1)

You rated this post out of 5. Change rating

Share

Share

Rate

5 (1)

You rated this post out of 5. Change rating