Implement simplified cone sampling algorithm.

* src/nd-random.sc (hollow-cone-random-vector): Implement simplified
algorithm that directly samples the surface of the sphere instead of
sampling a disk and projecting it onto the surface.

1 files changed, 8 insertions, 9 deletions

diff --git a/src/nd-random.sc b/src/nd-random.sc
index 02ad997..a9b08b3 100644
--- a/src/nd-random.sc
+++ b/src/nd-random.sc
@@ -98,16 +98,15 @@ dx. THETA should be in [0,pi]."
 (define (hollow-cone-random-vector r mean theta-min theta-max x)
   (void (const gsl-rng*) (const gsl-vector*) double double gsl-vector*)
   ;; Generate random vector around the nth canonical basis vector.
-  (let* ((n size-t (: x size)))
+  (let* ((n size-t (: x size))
+         (theta double (solid-angle->planar-angle
+                        (gsl-ran-flat r
+                                      (planar-angle->solid-angle theta-min n)
+                                      (planar-angle->solid-angle theta-max n))
+                        n)))
     (gsl-ran-dir-nd r (- n 1) (: x data))
-    (gsl-vector-scale x (* (cos theta-max)
-                           (tan (solid-angle->planar-angle
-                                 (gsl-ran-flat r
-                                               (planar-angle->solid-angle theta-min n)
-                                               (planar-angle->solid-angle theta-max n))
-                                 n))))
-    (gsl-vector-set x (- n 1) (cos theta-max)))
-  (gsl-vector-scale x (/ 1 (gsl-blas-dnrm2 x)))
+    (gsl-vector-scale x (sin theta))
+    (gsl-vector-set x (- n 1) (cos theta)))
 
   ;; Rotate to arbitrary basis.
   (rotate-from-nth-canonical x mean))