distance script

http://www.zipwise.com/online-zip-code-search.php

<%

':::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
'::: :::
'::: This routine calculates the distance between two points (given the :::
'::: latitude/longitude of those points). It is being used to calculate :::
'::: the distance between two ZIP Codes. :::
'::: :::
'::: Definitions: :::
'::: South latitudes are negative, east longitudes are positive :::
'::: :::
'::: Passed to function: :::
'::: lat1, lon1 = Latitude and Longitude of point 1 (in decimal degrees) :::
'::: lat2, lon2 = Latitude and Longitude of point 2 (in decimal degrees) :::
'::: unit = the unit you desire for results :::
'::: where: 'M' is statute miles :::
'::: 'K' is kilometers (default) :::
'::: 'N' is nautical miles :::
'::: :::
'::: United States ZIP Code databases with latitude :::
'::: & longitude are available at http://www.zipwise.com :::
'::: :::
':::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::

const pi = 3.14159265358979323846

Function distance(lat1, lon1, lat2, lon2, unit)
Dim theta, dist
theta = lon1 - lon2
dist = sin(deg2rad(lat1)) * sin(deg2rad(lat2)) + cos(deg2rad(lat1)) * cos(deg2rad(lat2)) * cos(deg2rad(theta))
dist = acos(dist)
dist = rad2deg(dist)
distance = dist * 60 * 1.1515
Select Case ucase(unit)
Case "K"
distance = distance * 1.609344
Case "N"
distance = distance * 0.8684
End Select
End Function


'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
'::: This function get the arccos function from arctan function :::
'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
Function acos(rad)
If Abs(rad) <> 1 Then
acos = pi/2 - Atn(rad / Sqr(1 - rad * rad))
ElseIf rad = -1 Then
acos = pi
End If
End function


'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
'::: This function converts decimal degrees to radians :::
'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
Function deg2rad(Deg)
deg2rad = cdbl(Deg * pi / 180)
End Function

'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
'::: This function converts radians to decimal degrees :::
'::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
Function rad2deg(Rad)
rad2deg = cdbl(Rad * 180 / pi)
End Function

response.write distance(32.9697, -96.80322, 29.46786, -98.53506, "M") & " Miles<br>"
response.write distance(32.9697, -96.80322, 29.46786, -98.53506, "K") & " Kilometers<br>"
response.write distance(32.9697, -96.80322, 29.46786, -98.53506, "N") & " Nautical Miles<br>"

%>