Surface area of a sphere sliced by two orthogonal planes

Question

I'm a programmer doing a graphics effect, and I've arrived at a math problem that my old dumb brain can't quite figure out the best way to solve.

Here's the problem: take a sphere and slice off a bit of it with a plane. You get a sphere cap and sliced sphere. Now, take another plane that is exactly orthogonal to the slicing plane and slice off another sphere cap. What is the surface area of the sliced sphere left over?

If the two sphere caps don't intersect, then the problem is simple: just subtract the surface area of the two sphere caps from the surface area of the sphere, and you're done. However, if the two planes are close enough to the center that the two sphere caps intersect, then there's a sliver of the sphere that gets sliced off twice, and you have to add it back. How do you do that?

In case anyone is curious about the actual thing I'm trying to do: I want to bake ambient occlusion into a texture that I'm going to use on a voxel mesh. The idea with AO is that you take a point on a surface, and a hemisphere centered around the surface normal, and you check how much of that surface area is occluded by objects. In a voxel world, that essentially reduces to the problem above, since the only thing that could occlude is neighboring voxels. I've googled this and the stuff I find are some variation of this method which uses vertex coloring (essentially) and looks quite bad. I want to do it properly.

Answer 1

Update: An exact formula for the area on a unit sphere with $x>a$, $y>b$ is $$A=\big(\frac 12-a-b\big)\pi + 2b\arctan(\frac ac)+2a\arctan(\frac bc)+\arctan(\frac{c^2-a^2 b^2}{2abc})$$ for positive $a$ and $b$ with $a^2+b^2<1$ and $c=\sqrt{1-a^2-b^2}$.

Finding this formula was difficult, requiring both some Mathematica and some guesswork for me. Fortunately we can verify it just by checking its values at endpoints and its changes in the middle.

In the limit as $a\rightarrow 0^+$, the last term goes to $\pi/2$, so this gives $A=(1-b)\pi$ and agrees with Archimedes's hatbox theorem. Similarly this is correct in the limit as $b\rightarrow 0^+$.

In the middle, $$\frac{d^2A}{da\ db}=\frac{2}{\sqrt{1-a^2-b^2}}=2\sqrt{1+f_x^2+f_y^2},$$ which is correctly twice the element of surface area for $f(x,y)=\sqrt{1-x^2-y^2}$, corresponding to the upper and lower hemispheres.

As an example, if $a=2/5$ and $b=3/4$, then $c=\sqrt{111}/20$ and $$A = \frac{-3\pi}{20} + \frac 45\arctan(\frac{15}{\sqrt{111}}) + \frac{3}{2}\arctan(\frac{8}{\sqrt{111}})-\arctan(\frac{4\sqrt{111}}{25}) \simeq .234333$$

Answer 2

Just to present another way with respect to the good answer above.

Let's refer to the attached 2-D orthogonal projection

We take a sphere of unitary radius, where the two cuts are individuated by the planes $x=a$ and $y=b$, or equivalently by the latitudinal angles $\alpha, \, \beta$.

The surface of the required "partial spherical wedge" will be described by the infinitesimal arc $d \phi$, rotated around the circle of radius $r$, at constant latitude $\phi$, by the logitudinal arc $2 r \theta$.

We have the following relations among the parameters. $$ \left\{ \matrix{ a = \cos \alpha \hfill \cr b = \sin \beta \hfill \cr r = \cos \phi \hfill \cr r\cos \theta = a \hfill \cr} \right. $$

So, clearly, we will have $$ S(\alpha ,\beta ) = 2\int_{\,\phi \, = \,\beta }^{\;\alpha } {r\,\theta \,d\phi } = 2\int_{\,\phi \, = \,\beta }^{\;\alpha } {\cos \phi \arccos \left( {\cos \alpha /\cos \phi } \right)d\phi } $$ where, if we want to preserve the positive sign for the integral we shall impose the inequality $\beta \le \alpha$.

I do not know (and doubt) whether the integral might be put into a closed form. Most probably you shall compute it numerically.

Answer 3

Area of piece $P$ in $\mathbb{S}^2$ :

i) Recall Gauss-Bonnet Theorem : $$ \int_{\Omega}\ K +\int_{\partial \Omega}\ k_g = 2\pi \chi(\Omega)$$ where $k_g$ is a geodesic curvature.

ii) Geodesic curvature on the boundary cap : $C_z$ is a cap wrt $z=\delta$ so that $$( {\rm area}\ C_z)_{=2\pi \delta} +\int_{\partial C_z}\ k_g =2\pi$$

Remark - Cone : And there is a cone $T_z$ which is tangent to $\mathbb{S}^2$ along $\partial C_z$. Note that a central angle of $T_z$ is $2\pi(1-\delta )$.

iii) Define two orthogonal planes $z=\delta,\ y=-1+\varepsilon$ in $\mathbb{E}^3$. Assume that two planes cut $P$ from $\mathbb{S}^2$

For $\partial C_z$, we have a parametrization $$ c_z(t)= \sqrt{2\delta-\delta^2 }(\cos\ t,\sin\ t,0) +(0,0,1-\delta) $$

Similarly, for $\partial C_y$ we have $$ c_y(s)=\sqrt{2\varepsilon -\varepsilon^2 }(\cos\ s,0,\sin\ s) +(0,-1+\varepsilon ,0) $$

Find intersection between $c_z$ and $c_y$ : $$ \sin\ t_i= \frac{-1+\varepsilon }{\sqrt{2\delta -\delta^2 }},\ \pi<t_1<3\pi/2 <t_2<2\pi,$$

$$ \sin\ s_i = \frac{1-\delta}{\sqrt{2\varepsilon -\varepsilon^2}},\ 0<s_2<\pi/2<s_1<\pi $$

Along smooth arc in $\partial P$, the integral of geodesic curvature is $$ |t_1-t_2|(1-\delta) + |s_1-s_2|(1-\varepsilon ) \ (:=T)$$

iv) Consider angle between two tangents at intersection point :

$$ c_z'(t_2)=\sqrt{2\delta-\delta^2 }(-\sin\ t_2,\cos\ t_2,0) = ( \varepsilon -1 ,\sqrt{ 2\delta-\delta^2 -(1-\varepsilon)^2 } ,0) $$

$$ c_y'(s_2)=\sqrt{2\varepsilon -\varepsilon^2 }(-\sin\ s_2,0,\cos\ s_2) = (\delta-1,0,\sqrt{2\varepsilon- \varepsilon^2 -(1-\delta)^2}) $$

Hence angle $\theta$ between two vectors is $$ \cos\ \theta=\frac{(1-\varepsilon )(1-\delta )}{\sqrt{2\varepsilon - \varepsilon^2}\sqrt{2 \delta -\delta^2}} $$

Hence Gauss-Bonnet theorem implies that $${\rm area}\ P= 2\pi- 2\theta-T $$