AndyH79
asked on
Conversion from Excel VBA to T-SQL not generating the same results.
I'm needing to convert the following code from Excel VBA to T-SQL. The problem that I am running into is that I can't get the system to generate the same result when run in T-SQL as when run in VBA. Both processes should return a result of 193585.85. Currently the T-SQL process is returning a result of -99 which is due to the debugging process identifying that all 20 iterations have occurred prior to a result being generated which makes Abs(@lambda - @lambdaP) > @EPSILON = FALSE
Since the VBA code uses a constant value for Pi I have used the same constant value as a declared variable rather than use the PI() SQL function since I'm unsure of the number of significant digits the SQL function uses.
The first issue that I run into is the calculation of:
@U1
Excel VBA code Line 62 U1 = Atn((1 - f) * Tan(toRad(lat1)))
returns 0.791295464197368
SQL Server code Line 68 SET @U1 = Atan((1 - @f) * Tan(RADIANS(@lat1)))
returns 0.791295464197410
@U2
Excel VBA code Line 62 U2 = Atn((1 - f) * Tan(toRad(lat2)))
returns 0.821498274972993
SQL Server code Line 69 -- SET @U2 = Atan((1 - @f) * Tan(RADIANS(@lat2)))
returns 0.821498274973035
Neither difference seems to be significant, I'm just not sure what is causing them to happen. Originally I was also manually converting the degrees to radians in SQL Server using the following function. Since the result was the same using the Radian() sql function I removed the references to this udf in the code.
CREATE FUNCTION [dbo].[udfMath_DegreesToRa dian]
(
-- Add the parameters for the function here
@degrees Decimal(38,30)
)
RETURNS Decimal(38,30)
AS
BEGIN
DECLARE @res Decimal(38,30)
DECLARE @PI Decimal(38,30)
SET @PI = 3.14159265358979
SELECT @res = (@degrees * (@PI /180))
RETURN @res
END
Complete VBA Code
Complete T-SQL Code:
Since the VBA code uses a constant value for Pi I have used the same constant value as a declared variable rather than use the PI() SQL function since I'm unsure of the number of significant digits the SQL function uses.
The first issue that I run into is the calculation of:
@U1
Excel VBA code Line 62 U1 = Atn((1 - f) * Tan(toRad(lat1)))
returns 0.791295464197368
SQL Server code Line 68 SET @U1 = Atan((1 - @f) * Tan(RADIANS(@lat1)))
returns 0.791295464197410
@U2
Excel VBA code Line 62 U2 = Atn((1 - f) * Tan(toRad(lat2)))
returns 0.821498274972993
SQL Server code Line 69 -- SET @U2 = Atan((1 - @f) * Tan(RADIANS(@lat2)))
returns 0.821498274973035
Neither difference seems to be significant, I'm just not sure what is causing them to happen. Originally I was also manually converting the degrees to radians in SQL Server using the following function. Since the result was the same using the Radian() sql function I removed the references to this udf in the code.
CREATE FUNCTION [dbo].[udfMath_DegreesToRa
(
-- Add the parameters for the function here
@degrees Decimal(38,30)
)
RETURNS Decimal(38,30)
AS
BEGIN
DECLARE @res Decimal(38,30)
DECLARE @PI Decimal(38,30)
SET @PI = 3.14159265358979
SELECT @res = (@degrees * (@PI /180))
RETURN @res
END
Complete VBA Code
Option Explicit
Private Const PI = 3.14159265358979
Private Const EPSILON As Double = 0.000000000001
Public Sub test()
Dim metersDistance As Double
metersDistance = Geo.GetDistance(45.434094, -122.5791498, 47.1643339, -122.293198)
End Sub
Public Function GetDistance(ByVal lat1 As Double, ByVal lon1 As Double, ByVal lat2 As Double, ByVal lon2 As Double) As Double
'=================================================================================
' Calculate geodesic distance (in m) between two points specified by latitude/longitude (in numeric degrees)
' using Vincenty inverse formula for ellipsoids
Dim low_a As Double
Dim low_b As Double
Dim f As Double
Dim L As Double
Dim U1 As Double
Dim U2 As Double
Dim sinU1 As Double
Dim sinU2 As Double
Dim cosU1 As Double
Dim cosU2 As Double
Dim lambda As Double
Dim lambdaSub1 As Double
Dim lambdaSub2 As Double
Dim lambdaSub3 As Double
Dim lambdaP As Double
Dim iterLimit As Integer
Dim sinLambda As Double
Dim cosLambda As Double
Dim sinSigma As Double
Dim cosSigma As Double
Dim sigma As Double
Dim sinAlpha As Double
Dim cosSqAlpha As Double
Dim cos2SigmaM As Double
Dim C As Double
Dim uSq As Double
Dim upper_A As Double
Dim upper_ASub1 As Double
Dim upper_ASub2 As Double
Dim upper_ASub3 As Double
Dim upper_B As Double
Dim upper_BSub1 As Double
Dim upper_BSub2 As Double
Dim upper_BSub3 As Double
Dim deltaSigma As Double
Dim deltaSigmaSub1 As Double
Dim deltaSigmaSub2 As Double
Dim deltaSigmaSub3 As Double
Dim s As Double
low_a = 6378137
low_b = 6356752.31414
f = 1 / 298.257222101 'GRS-80 ellipsiod
L = toRad(lon2 - lon1)
U1 = Atn((1 - f) * Tan(toRad(lat1)))
U2 = Atn((1 - f) * Tan(toRad(lat2)))
sinU1 = Sin(U1)
cosU1 = Cos(U1)
sinU2 = Sin(U2)
cosU2 = Cos(U2)
lambda = L
lambdaP = 2 * PI
iterLimit = 20
While (Abs(lambda - lambdaP) > EPSILON) And (iterLimit > 0)
iterLimit = iterLimit - 1
sinLambda = Sin(lambda)
cosLambda = Cos(lambda)
sinSigma = Sqr(((cosU2 * sinLambda) ^ 2) + ((cosU1 * sinU2 - sinU1 * cosU2 * cosLambda) ^ 2))
If sinSigma = 0 Then
GetDistance = 0 'co-incident points
Exit Function
End If
cosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLambda
sigma = Atan2(cosSigma, sinSigma)
sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma
cosSqAlpha = 1 - sinAlpha * sinAlpha
If cosSqAlpha = 0 Then 'check for a divide by zero
cos2SigmaM = 0 '2 points on the equator
Else
cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlpha
End If
C = f / 16 * cosSqAlpha * (4 + f * (4 - 3 * cosSqAlpha))
lambdaP = lambda
'Had to split the calc of lambda since excel says is too complex
'lambda = L + (1 - C) * f * sinAlpha * (sigma + C * sinSigma * (cos2SigmaM + C * cosSigma * (-1 + 2 * (cos2SigmaM ^ 2))))
lambdaSub1 = (1 - C) * f * sinAlpha
lambdaSub2 = (-1 + 2 * (cos2SigmaM ^ 2))
lambdaSub3 = (cos2SigmaM + C * cosSigma * lambdaSub2)
lambda = L + lambdaSub1 * (sigma + C * sinSigma * lambdaSub3)
Wend
If iterLimit < 1 Then
MsgBox "iteration limit has been reached, something didn't work."
Exit Function
End If
uSq = cosSqAlpha * (low_a ^ 2 - low_b ^ 2) / (low_b ^ 2)
'Had to split the calc of upper_A since excel says is too complex
'upper_A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)))
upper_ASub1 = 320 - 175 * uSq
upper_ASub2 = (-768 + uSq * (upper_ASub1))
upper_ASub3 = (4096 + uSq * upper_ASub2)
upper_A = 1 + uSq / 16384 * upper_ASub3
'Had to split the calc of upper_B since excel says is too complex
'upper_B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)))
upper_BSub1 = (74 - 47 * uSq)
upper_BSub2 = (-128 + uSq * upper_BSub1)
upper_BSub3 = (256 + uSq * upper_BSub2)
upper_B = uSq / 1024 * upper_BSub3
'Had to split the calc of deltaSigma since excel says is too complex
'deltaSigma = upper_B * sinSigma * (cos2SigmaM + upper_B / 4 * (cosSigma * (-1 + 2 * cos2SigmaM ^ 2) - upper_B / 6 * cos2SigmaM * (-3 + 4 * sinSigma ^ 2) * (-3 + 4 * cos2SigmaM ^ 2))
deltaSigmaSub1 = (-3 + 4 * sinSigma ^ 2) * (-3 + 4 * cos2SigmaM ^ 2)
deltaSigmaSub2 = (cosSigma * (-1 + 2 * cos2SigmaM ^ 2) - upper_B / 6 * cos2SigmaM * deltaSigmaSub1)
deltaSigmaSub3 = (cos2SigmaM + upper_B / 4 * deltaSigmaSub2)
deltaSigma = upper_B * sinSigma * deltaSigmaSub3
s = low_b * upper_A * (sigma - deltaSigma)
GetDistance = Round(s, 3) 'Rounds distance to nearest millimeter.
End Function
Private Function toRad(ByVal degrees As Double) As Double
toRad = degrees * (PI / 180)
End Function
Private Function Atan2(ByVal X As Double, ByVal y As Double) As Double
If y > 0 Then
If X >= y Then
Atan2 = Atn(y / X)
ElseIf X <= -y Then
Atan2 = Atn(y / X) + PI
Else
Atan2 = PI / 2 - Atn(X / y)
End If
Else
If X >= -y Then
Atan2 = Atn(y / X)
ElseIf X <= y Then
Atan2 = Atn(y / X) - PI
Else
Atan2 = -Atn(X / y) - PI / 2
End If
End If
End Function
Complete T-SQL Code:
DECLARE @lat1 FLOAT
DECLARE @lon1 FLOAT
DECLARE @lat2 FLOAT
DECLARE @lon2 FLOAT
-- 0027
SET @lat1 = 45.4340940
SET @lon1 = -122.5791498
-- 0601
SET @lat2 = 47.1643339
SET @lon2 =-122.2931980
DECLARE @EPSILON FLOAT
DECLARE @PI FLOAT
SET @EPSILON = 0.000000000001
SET @PI = 3.14159265358979
DECLARE @low_a FLOAT,
@low_b FLOAT,
@f FLOAT,
@L FLOAT,
@U1 FLOAT,
@U2 FLOAT,
@sinU1 FLOAT,
@sinU2 FLOAT,
@cosU1 FLOAT,
@cosU2 FLOAT,
@lambda FLOAT,
@lambdaSub1 FLOAT,
@lambdaSub2 FLOAT,
@lambdaSub3 FLOAT,
@lambdaP FLOAT,
@iterLimit TinyInt,
@sinLambda FLOAT,
@cosLambda FLOAT,
@sinSigma FLOAT,
@cosSigma FLOAT,
@sigma FLOAT,
@sinAlpha FLOAT,
@cosSqAlpha FLOAT,
@cos2SigmaM FLOAT,
@C FLOAT,
@uSq FLOAT,
@upper_A FLOAT,
@upper_ASub1 FLOAT,
@upper_ASub2 FLOAT,
@upper_ASub3 FLOAT,
@upper_B FLOAT,
@upper_BSub1 FLOAT,
@upper_BSub2 FLOAT,
@upper_BSub3 FLOAT,
@deltaSigma FLOAT,
@deltaSigmaSub1 FLOAT,
@deltaSigmaSub2 FLOAT,
@deltaSigmaSub3 FLOAT,
@s FLOAT
DECLARE @res FLOAT
SET @low_a = 6378137
SET @low_b = 6356752.31414
SET @f = 1 / 298.257222101 --GRS-80 ellipsiod
SET @L = RADIANS(@lon2 - @lon1)
SET @U1 = Atan((1 - @f) * Tan(RADIANS(@lat1)))
SET @U2 = Atan((1 - @f) * Tan(RADIANS(@lat2)))
--Included for debugging purposes
SELECT @L as L, @U1 as U1, @U2 as U2
SET @sinU1 = Sin(@U1)
SET @cosU1 = Cos(@U1)
SET @sinU2 = Sin(@U2)
SET @cosU2 = Cos(@U2)
SET @lambda = @L
SET @lambdaP = 2 * @PI
SET @iterLimit = 20
--Included for debugging purposes
--SELECT @sinU1 as sinU1, @cosU1 as cosU1, @sinU2 as sinU2,@cosU2 as cosU2, @lambda as lambda, @lambdaP as lambdaP, @iterLimit as iterLimit
While (Abs(@lambda - @lambdaP) > @EPSILON) AND (@iterLimit > 0)
BEGIN
SET @iterLimit = @iterLimit - 1
SET @sinLambda = Sin(@lambda)
SET @cosLambda = Cos(@lambda)
SET @sinSigma = Square((Square((@cosU2 * @sinLambda))) + (Square(@cosU1 * @sinU2 - @sinU1 * @cosU2 * @cosLambda)))
IF @sinSigma = 0
BEGIN
SELECT @iterLimit as Iterations, -98 AS ErrorCode--co-incident points
--RETURN @res
END
SET @cosSigma = @sinU1 * @sinU2 + @cosU1 * @cosU2 * @cosLambda
SET @sigma = atn2(@cosSigma, @sinSigma)
SET @sinAlpha = @cosU1 * @cosU2 * @sinLambda / @sinSigma
SET @cosSqAlpha = 1 - @sinAlpha * @sinAlpha
If @cosSqAlpha = 0 --check for a divide by zero
Begin
SET @cos2SigmaM = 0 --2 points on the equator
END
Else
Begin
SET @cos2SigmaM = @cosSigma - 2 * @sinU1 * @sinU2 / @cosSqAlpha
End
SET @C = @f / 16 * @cosSqAlpha * (4 + @f * (4 - 3 * @cosSqAlpha))
SET @lambdaP = @lambda
--Had to split the calc of lambda since excel says is too complex
--lambda = L + (1 - C) * f * sinAlpha * (sigma + C * sinSigma * (cos2SigmaM + C * cosSigma * (-1 + 2 * (Power(cos2SigmaM,2)))))
SET @lambdaSub1 = (1 - @C) * @f * @sinAlpha
SET @lambdaSub2 = (-1 + 2 * (Square(@cos2SigmaM)))
SET @lambdaSub3 = (@cos2SigmaM + @C * @cosSigma * @lambdaSub2)
SET @lambda = @L + @lambdaSub1 * (@sigma + @C * @sinSigma * @lambdaSub3)
END
If @iterLimit < 1
BEGIN
SET @res = -99 AS ErrorCode
SELECT @res --iteration limit has been reached, something didn't work."
--RETURN @res
END
SET @uSq = @cosSqAlpha * (Square(@low_a) - Square(@low_b) / Square(@low_b))
--Had to split the calc of upper_A since excel says is too complex
--upper_A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)))
SET @upper_ASub1 = 320 - 175 * @uSq
SET @upper_ASub2 = (-768 + @uSq * (@upper_ASub1))
SET @upper_ASub3 = (4096 + @uSq * @upper_ASub2)
SET @upper_A = 1 + @uSq / 16384 * @upper_ASub3
--Had to split the calc of upper_B since excel says is too complex
--upper_B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)))
SET @upper_BSub1 = (74 - 47 * @uSq)
SET @upper_BSub2 = (-128 + @uSq * @upper_BSub1)
SET @upper_BSub3 = (256 + @uSq * @upper_BSub2)
SET @upper_B = @uSq / 1024 * @upper_BSub3
--Had to split the calc of deltaSigma since excel says is too complex
--deltaSigma = upper_B * sinSigma * (cos2SigmaM + upper_B / 4 * (cosSigma * (-1 + 2 * Power(cos2SigmaM,2)) - upper_B / 6 * cos2SigmaM * (-3 + 4 * Power(sinSigma,2)) * (-3 + 4 * Power(cos2SigmaM,2)))
SET @deltaSigmaSub1 = (-3 + 4 * Square(@sinSigma)) * (-3 + 4 * Square(@cos2SigmaM))
SET @deltaSigmaSub2 = (@cosSigma * (-1 + 2 * Square(@cos2SigmaM)) - @upper_B / 6 * @cos2SigmaM * @deltaSigmaSub1)
SET @deltaSigmaSub3 = (@cos2SigmaM + @upper_B / 4 * @deltaSigmaSub2)
SET @deltaSigma = @upper_B * @sinSigma * @deltaSigmaSub3
SET @s = @low_b * @upper_A * (@sigma - @deltaSigma)
If @res <> -98 and @res <> -99
BEGIN
SET @res = Round(@s, 3)--Rounds distance to nearest millimeter.
SELECT @res
END
ASKER CERTIFIED SOLUTION
membership
This solution is only available to members.
To access this solution, you must be a member of Experts Exchange.
ASKER
Both the SQL and VBA code make use of a locally declared value for PI = 3.14159265358979. For sql it is @PI. As a test when I update the VBA PI constant to 3.1415926535897931 it only results in a change in value of 0.0000000000000001where the values for @U1 and @U2 are 0.000000000000042 different from their VBA counterparts. Could this be due to some rounding error because I'm using float values instead of decimal values? I may look at doing this as a CLR in the future.
Hi,
I've been working on an iterative process-calculation. What I found is that is division using decimal only returns 6 decimal places, so I'm using float, which equates to double in VBA.
Potentially there maybe a little difference in the returns from tan and arctan functions. Have you checked them individually?
Something to try is in every place you want float values returned specify float constants ie 180.0 instead of 180 and see if that makes a difference. It made a significant difference at one point for my project.
Regards
David
I've been working on an iterative process-calculation. What I found is that is division using decimal only returns 6 decimal places, so I'm using float, which equates to double in VBA.
Potentially there maybe a little difference in the returns from tan and arctan functions. Have you checked them individually?
Something to try is in every place you want float values returned specify float constants ie 180.0 instead of 180 and see if that makes a difference. It made a significant difference at one point for my project.
Regards
David
ASKER
I took your advice and converted it into a CLR function.
[Microsoft.SqlServer.Serve r.SqlFunct ion]
public static SqlString udfMath_InverseSolutionOfG eodesics
(
Double lat1,
Double lon1,
Double lat2,
Double lon2
)
{
//Calculate geodesic distance (in m) between two points specified by latitude/longitude (in numeric degrees)
//using Vincenty inverse formula for ellipsoids
Double result;
Double a;
Double b;
Double f;
Double L;
Double U1;
Double U2;
Double sinU1;
Double sinU2;
Double cosU1;
Double cosU2;
Double lambda;
Double lambdaP;
int iterLimit;
Double sinLambda;
Double cosLambda;
Double sinSigma = 0;
Double cosSigma = 0;
Double sigma = 0;
Double sinAlpha;
Double cosSqAlpha = 0;
Double cos2SigmaM = 0;
Double C;
Double uSq;
Double A;
Double B;
Double deltaSigma;
Double s;
Double PI = 3.14159265358979;
Double EPSILON = 0.000000000001;
a = 6378137;
b = 6356752.31414;
f = 1 / 298.257222101; //GRS-80 ellipsiod
L = PI * (lon2-lon1) / 180.0;
U1 = Math.Atan((1 - f) * Math.Tan(PI * (lat1) / 180.0));
U2 = Math.Atan((1 - f) * Math.Tan(PI * (lat2) / 180.0));
sinU1 = Math.Sin(U1);
cosU1 = Math.Cos(U1);
sinU2 = Math.Sin(U2);
cosU2 = Math.Cos(U2);
lambda = L;
lambdaP = 2 * PI;
iterLimit = 20;
while (Math.Abs(lambda - lambdaP) > EPSILON && iterLimit > 0)
{
iterLimit = iterLimit - 1;
sinLambda = Math.Sin(lambda);
cosLambda = Math.Cos(lambda);
sinSigma = Math.Sqrt(Math.Pow((cosU2 * sinLambda),2) + Math.Pow((cosU1 * sinU2 - sinU1 * cosU2 * cosLambda),2));
if (sinSigma == 0)
{
result = -98; //co-incident points
}
cosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLambda;
sigma = Math.Atan2(sinSigma,cosSig ma);
sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma;
cosSqAlpha = 1 - sinAlpha * sinAlpha;
if (cosSqAlpha == 0) //check for a divide by zero
{
cos2SigmaM = 0; //2 points on the equator
}
else
{
cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlpha;
}
C = f / 16 * cosSqAlpha * (4 + f * (4 - 3 * cosSqAlpha));
lambdaP = lambda;
lambda = L + (1 - C) * f * sinAlpha * (sigma + C * sinSigma * (cos2SigmaM + C * cosSigma * (-1 + 2 * Math.Pow(cos2SigmaM,2))));
}
if (iterLimit < 1)
{
result = -99; //iteration limit has been reached, something didn't work.
}
uSq = cosSqAlpha * (Math.Pow(a, 2) - Math.Pow(b, 2)) / Math.Pow(b, 2);
A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)));
B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)));
deltaSigma = B * sinSigma * (cos2SigmaM + B / 4 * (cosSigma * (-1 + 2 * (Math.Pow(cos2SigmaM,2) - B / 6 * cos2SigmaM * (-3 + 4 * Math.Pow(sinSigma,2) * (-3 + 4 * Math.Pow(cos2SigmaM,2))))) ));
s = b * A * (sigma - deltaSigma);
result = Math.Round(s, 3); //Rounds distance to nearest millimeter
return new SqlString(result.ToString( ));
}
[Microsoft.SqlServer.Serve
public static SqlString udfMath_InverseSolutionOfG
(
Double lat1,
Double lon1,
Double lat2,
Double lon2
)
{
//Calculate geodesic distance (in m) between two points specified by latitude/longitude (in numeric degrees)
//using Vincenty inverse formula for ellipsoids
Double result;
Double a;
Double b;
Double f;
Double L;
Double U1;
Double U2;
Double sinU1;
Double sinU2;
Double cosU1;
Double cosU2;
Double lambda;
Double lambdaP;
int iterLimit;
Double sinLambda;
Double cosLambda;
Double sinSigma = 0;
Double cosSigma = 0;
Double sigma = 0;
Double sinAlpha;
Double cosSqAlpha = 0;
Double cos2SigmaM = 0;
Double C;
Double uSq;
Double A;
Double B;
Double deltaSigma;
Double s;
Double PI = 3.14159265358979;
Double EPSILON = 0.000000000001;
a = 6378137;
b = 6356752.31414;
f = 1 / 298.257222101; //GRS-80 ellipsiod
L = PI * (lon2-lon1) / 180.0;
U1 = Math.Atan((1 - f) * Math.Tan(PI * (lat1) / 180.0));
U2 = Math.Atan((1 - f) * Math.Tan(PI * (lat2) / 180.0));
sinU1 = Math.Sin(U1);
cosU1 = Math.Cos(U1);
sinU2 = Math.Sin(U2);
cosU2 = Math.Cos(U2);
lambda = L;
lambdaP = 2 * PI;
iterLimit = 20;
while (Math.Abs(lambda - lambdaP) > EPSILON && iterLimit > 0)
{
iterLimit = iterLimit - 1;
sinLambda = Math.Sin(lambda);
cosLambda = Math.Cos(lambda);
sinSigma = Math.Sqrt(Math.Pow((cosU2 * sinLambda),2) + Math.Pow((cosU1 * sinU2 - sinU1 * cosU2 * cosLambda),2));
if (sinSigma == 0)
{
result = -98; //co-incident points
}
cosSigma = sinU1 * sinU2 + cosU1 * cosU2 * cosLambda;
sigma = Math.Atan2(sinSigma,cosSig
sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma;
cosSqAlpha = 1 - sinAlpha * sinAlpha;
if (cosSqAlpha == 0) //check for a divide by zero
{
cos2SigmaM = 0; //2 points on the equator
}
else
{
cos2SigmaM = cosSigma - 2 * sinU1 * sinU2 / cosSqAlpha;
}
C = f / 16 * cosSqAlpha * (4 + f * (4 - 3 * cosSqAlpha));
lambdaP = lambda;
lambda = L + (1 - C) * f * sinAlpha * (sigma + C * sinSigma * (cos2SigmaM + C * cosSigma * (-1 + 2 * Math.Pow(cos2SigmaM,2))));
}
if (iterLimit < 1)
{
result = -99; //iteration limit has been reached, something didn't work.
}
uSq = cosSqAlpha * (Math.Pow(a, 2) - Math.Pow(b, 2)) / Math.Pow(b, 2);
A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)));
B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)));
deltaSigma = B * sinSigma * (cos2SigmaM + B / 4 * (cosSigma * (-1 + 2 * (Math.Pow(cos2SigmaM,2) - B / 6 * cos2SigmaM * (-3 + 4 * Math.Pow(sinSigma,2) * (-3 + 4 * Math.Pow(cos2SigmaM,2)))))
s = b * A * (sigma - deltaSigma);
result = Math.Round(s, 3); //Rounds distance to nearest millimeter
return new SqlString(result.ToString(
}
ASKER
This is used for determining the distance between to coordinates in meters. Credit for the original code goes to http://www.movable-type.co.uk/scripts/latlong-vincenty.html
In SQL PI is
3.1415926535897931
so try that in your VBA code
select convert( decimal( 38, 30 ), pi()) as pie_as_decimal, convert( float( 53 ), pi()) as pie_as_float;
Your differences start in the 13th digit and to me this is the difference in the two values of pi.
HTH
David