SPHERICALANGLE Compute angle between points on the sphere. ALPHA = sphericalAngle(P1, P2, P3) Computes angle (P1, P2, P3), i.e. the angle, measured at point P2, between the direction (P2, P1) and the direction (P2, P3). The result is given in radians, between 0 and 2*PI. Points are given either as [x y z] (there will be normalized to lie on the unit sphere), or as [phi theta], with phi being the longitude in [0 2*PI] and theta being the elevation on horizontal [-pi/2 pi/2]. NOTE: this is an 'oriented' version of the angle computation, that is, the result of sphericalAngle(P1, P2, P3) equals 2*pi-sphericalAngle(P3,P2,P1). To have the more classical relation (with results given betwen 0 and PI), it suffices to take the minimum of angle and 2*pi-angle. Examples % Use inputs as cartesian coordinates p1 = [0 1 0]; p2 = [1 0 0]; p3 = [0 0 1]; alpha = sphericalAngle(p1, p2, p3) alpha = 1.5708 % Use inputs as spherical coordinates sph1 = [.1 0]; sph2 = [0 0]; sph3 = [0 .1]; alphas = sphericalAngle(sph1, sph2, sph3) alphas = 1.5708 See also geom3d, angles3d, spheres, sph2cart
0001 function alpha = sphericalAngle(p1, p2, p3) 0002 %SPHERICALANGLE Compute angle between points on the sphere. 0003 % 0004 % ALPHA = sphericalAngle(P1, P2, P3) 0005 % Computes angle (P1, P2, P3), i.e. the angle, measured at point P2, 0006 % between the direction (P2, P1) and the direction (P2, P3). 0007 % The result is given in radians, between 0 and 2*PI. 0008 % 0009 % Points are given either as [x y z] (there will be normalized to lie on 0010 % the unit sphere), or as [phi theta], with phi being the longitude in [0 0011 % 2*PI] and theta being the elevation on horizontal [-pi/2 pi/2]. 0012 % 0013 % 0014 % NOTE: 0015 % this is an 'oriented' version of the angle computation, that is, the 0016 % result of sphericalAngle(P1, P2, P3) equals 0017 % 2*pi-sphericalAngle(P3,P2,P1). To have the more classical relation 0018 % (with results given betwen 0 and PI), it suffices to take the minimum 0019 % of angle and 2*pi-angle. 0020 % 0021 % Examples 0022 % % Use inputs as cartesian coordinates 0023 % p1 = [0 1 0]; 0024 % p2 = [1 0 0]; 0025 % p3 = [0 0 1]; 0026 % alpha = sphericalAngle(p1, p2, p3) 0027 % alpha = 0028 % 1.5708 0029 % 0030 % % Use inputs as spherical coordinates 0031 % sph1 = [.1 0]; 0032 % sph2 = [0 0]; 0033 % sph3 = [0 .1]; 0034 % alphas = sphericalAngle(sph1, sph2, sph3) 0035 % alphas = 0036 % 1.5708 0037 % 0038 % 0039 % See also 0040 % geom3d, angles3d, spheres, sph2cart 0041 0042 % ------ 0043 % Author: David Legland 0044 % E-mail: david.legland@inrae.fr 0045 % Created: 2005-02-21 0046 % Copyright 2005-2024 INRA - TPV URPOI - BIA IMASTE 0047 0048 % test if points are given as matlab spherical coordinates 0049 if size(p1, 2) == 2 0050 [x, y, z] = sph2cart(p1(:,1), p1(:,2), ones(size(p1,1), 1)); 0051 p1 = [x y z]; 0052 [x, y, z] = sph2cart(p2(:,1), p2(:,2), ones(size(p2,1), 1)); 0053 p2 = [x y z]; 0054 [x, y, z] = sph2cart(p3(:,1), p3(:,2), ones(size(p3,1), 1)); 0055 p3 = [x y z]; 0056 end 0057 0058 % normalize points 0059 p1 = normalizeVector3d(p1); 0060 p2 = normalizeVector3d(p2); 0061 p3 = normalizeVector3d(p3); 0062 0063 % create the plane tangent to the unit sphere and containing central point 0064 plane = createPlane(p2, p2); 0065 0066 % project the two other points on the plane 0067 pp1 = planePosition(projPointOnPlane(p1, plane), plane); 0068 pp3 = planePosition(projPointOnPlane(p3, plane), plane); 0069 0070 % compute angle on the tangent plane 0071 pp2 = zeros(max(size(pp1, 1), size(pp3,1)), 2); 0072 alpha = angle3Points(pp1, pp2, pp3); 0073