Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: compute surface area within closed loop on sphere
Date: Tue, 8 Jan 2013 20:22:09 +0000 (UTC)
Organization: The MathWorks, Inc.
Lines: 19
Message-ID: <kchv5h$rqv$1@newscl01ah.mathworks.com>
References: <kcdovp$lo7$1@newscl01ah.mathworks.com>
Reply-To: <HIDDEN>
NNTP-Posting-Host: www-00-blr.mathworks.com
Content-Type: text/plain; charset=UTF-8; format=flowed
Content-Transfer-Encoding: 8bit
X-Trace: newscl01ah.mathworks.com 1357676529 28511 172.30.248.45 (8 Jan 2013 20:22:09 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Tue, 8 Jan 2013 20:22:09 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1187260
Xref: news.mathworks.com comp.soft-sys.matlab:786208

"Peter Schreiber" <schreiber.peter15@gmail.com> wrote in message <kcdovp$lo7$1@newscl01ah.mathworks.com>...
> I have a closed loop on a sphere and would like to compute the surface area enclosed within the closed loop. I don't have an analytical expression for the curve (closed loop) on the sphere surface.
> How would I go about calculating the surface area?
- - - - - - - - - -
  Let p0 = [x0,y0,z0] be a point on the sphere's surface preferably interior to the closed loop, and let P1 be an n x 3 array of coordinates of your discrete closed loop points, where we suppose that the sphere has its center at the origin and radius R.  (If the center is not at the origin, subtract it from p0 and P1.)  Then you can break up the region inside the loop into n spherical triangles where p0 is always one vertex and two successive points on P1 are the other two vertices.  For each triangle you can use Lhuillier's formula to compute its surface area.  Their sum will then give you the total area.  The following code provides a signed area for each triangle in order to allow triangles to cancel out that may have overlapped due to the nature of the P1 curve.  Consequently P1 should make a net total of one loop around p0 in a counterclockwise direction in order to obtain a 
positive result. 

 P0 = repmat(p0,size(P1,1),1);
 P2 = P1([2:end,1],:);
 P01 = atan2(sqrt(sum(cross(P0,P1,2).^2,2)),dot(P0,P1,2));
 P12 = atan2(sqrt(sum(cross(P1,P2,2).^2,2)),dot(P1,P2,2));
 P20 = atan2(sqrt(sum(cross(P2,P0,2).^2,2)),dot(P2,P0,2));
 S = (P01+P12+P20)/2;
 A = R^2*sum( sign(dot(P0,cross(P1,P2,2),2)) .* ...
  (4*atan(sqrt(tan(S/2).*tan((S-P01)/2).*tan((S-P12)/2).*tan((S-P20)/2)))) );

A is the area.

Roger Stafford