I am reading a matlab function for calculating great circle distance written by my senior collegue. The distance between two points on the earth surface should be calculated using this formula:

d = r * arccos[(sin(lat1) * sin(lat2))   cos(lat1) * cos(lat2) * cos(long2 – long1)]

However, the script has the code like this:

dist = (acos(cos(pi/180*(90-lat2)).*cos(pi/180*(90-lat1)) sin(pi/180*(90-lat2)).*sin(pi/180*(90-lat1)).*cos(pi/180*(diff_long)))) .* r_local;
(-180 < long1,long2 <= 180, -90 < lat1,lat2 <= 90)

Why are sin(pi/2-A) and cos(pi/2-A) used to replace cos(A) and sin(A)？ Doesn't it introduced more error source by using the constant pi? Since lat1, lat2 might be very close to zero in my work, is this a trick on the numerical accuracy of MATLAB's sin() and cos() function?

Look forward to answers that explain how trigonometric functions in MATLAB work and analyze the error of these functions when the argument is close or equal to 0 and pi/2.

CodePudding user response：

If the purpose is to increase accuracy, this seems a very poor idea. When the angle is small, 90-A spoils any accuracy. That even makes tiny angles vanish (90-ε=90).

On the opposite, the sine of tiny angles is very close to the angle itself (radians) and for this reason quite accurately computed, while the cosine is virtually 1 or 1-A²/2. For top accuracy on tiny angles, you may resort to the versine, using versin(A):= 1-cos(A) = 2 sin²(A/2) and rework the equations in terms of 1-versin(A) instead of cos(A).

If the angle is close to 90°, accuracy is lost anyway, 90°-A will not restore it.

CodePudding user response：

I very much doubt this has to do with accuracy. Or at least, I don't think this helps any when it comes to accuracy.

The maximum difference between both sin(pi/2-A) - cos(A) and cos(pi/2-A) - sin(A) is 1.1102e-16, which is very small. This is just basic floating point accuracy, and there's really no way of telling which of the numbers is more correct. Note that cos(pi/2) = 6.1232e-17. So, if theta = 0, your colleague's code cos(pi/2-0) will give an error of 6.1232e-17, while simply doing the obvious sin(0) will be correct.

If you need numbers that are more accurate than this then you can try vpa.

I guess this is either because your colleague found another formula and implemented that, or he/she's confused and has attempted to increase the accuracy.

The latter might be the case if he/she tried to avoid the approximations sin(theta) ≈ theta and cos(theta) ≈ 1 for small values of theta. However, this doesn't make sense, since cos(pi/2-theta) ≈ theta and sin(pi/2-theta) ≈ 1 for small values of theta.