star-sf

commit 8a5c5d917397aa7471c11b08d8db6bf5f7742d12
parent a568b9b7b1603905a49ee0dfdbadd2fafe0fd010
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Tue,  9 Feb 2021 10:22:26 +0100

Test the ssf_phase_rdgfa_sample function

Diffstat:
M
src/ssf_phase_rdgfa.c
 | 
6
+++++-
M
src/test_ssf_phase_rdgfa.c
 | 
20
++++++++++++++++++++

2 files changed, 25 insertions(+), 1 deletion(-)
diff --git a/src/ssf_phase_rdgfa.c b/src/ssf_phase_rdgfa.c
@@ -264,6 +264,7 @@ rdgfa_sample
   double thetas[2];
   double theta;
   double sin_theta;
+  double cos_theta;
   double phi;
   double r;
   double u;
@@ -297,12 +298,13 @@ rdgfa_sample
   /* Uniformly sample a phi angle in [0, 2PI[ */
   phi = ssp_rng_uniform_double(rng, 0, 2*PI);
   sin_theta = sin(theta);
+  cos_theta = cos(theta);
 
   /* Compute the cartesian coordinates of the sampled direction in the _local_
    * phase function space */
   wi[0] = cos(phi) * sin_theta;
   wi[1] = sin(phi) * sin_theta;
-  wi[2] = cos(theta);
+  wi[2] = cos_theta;
 
   /* Compute the transformation matrix from local phase function to world
    * space. Note that by convention, in Star-SF the directions point outward
@@ -311,6 +313,8 @@ rdgfa_sample
   d33_basis(frame, d3_minus(w, wo));
   d33_muld3(wi, frame, wi);
 
+  ASSERT(eq_eps(d3_dot(wi, w), cos_theta, fabs(cos_theta*1.e-6)));
+
   if(pdf) *pdf = rdgfa_eval(rdgfa, wo, wi);
 }
 
diff --git a/src/test_ssf_phase_rdgfa.c b/src/test_ssf_phase_rdgfa.c
@@ -20,6 +20,8 @@
 int
 main(int argc, char** argv)
 {
+  static const size_t NSAMPS = 10000;
+
   struct mem_allocator allocator;
   struct ssf_phase_rdgfa_setup_args args = SSF_PHASE_RDGFA_SETUP_ARGS_DEFAULT;
   struct ssf_phase_rdgfa_desc desc = SSF_PHASE_RDGFA_DESC_NULL;
@@ -27,6 +29,7 @@ main(int argc, char** argv)
   struct ssf_phase* phase;
   struct ssf_phase* dummy;
   struct ssp_rng* rng;
+  double wo[3];
   size_t i;
   (void)argc, (void)argv;
 
@@ -63,12 +66,28 @@ main(int argc, char** argv)
 
   FOR_EACH(i, 0, desc.nintervals) {
     const struct ssf_phase_rdgfa_interval* ref = phase_rdgfa_ref_cumul + i;
+
     CHK(ssf_phase_rdgfa_get_interval(phase, i, &interval) == RES_OK);
+
     CHK(eq_eps(interval.range[0], ref->range[0], fabs(ref->range[0]*1.e-6)));
     CHK(eq_eps(interval.range[1], ref->range[1], fabs(ref->range[1]*1.e-6)));
     CHK(eq_eps(interval.cumulative, ref->cumulative, fabs(ref->cumulative*1.e-6)));
   }
 
+  ssp_ran_sphere_uniform(rng, wo, NULL);
+  FOR_EACH(i, 0, NSAMPS) {
+    double wi[3];
+    double pdf;
+    ssf_phase_sample(phase, rng, wo, wi, &pdf);
+    CHK(eq_eps(pdf, ssf_phase_eval(phase, wo, wi), fabs(pdf*1.e-6)));
+    CHK(d3_is_normalized(wi));
+
+#if 0
+    fprintf(stderr, "v %g %g %g\n", wi[0]*pdf, wi[1]*pdf, wi[2]*pdf);
+    fprintf(stderr, "p %lu\n", (unsigned long)(i+1));
+#endif
+  }
+
   CHK(ssf_phase_ref_put(phase) == RES_OK);
   CHK(ssf_phase_ref_put(dummy) == RES_OK);
   CHK(ssp_rng_ref_put(rng) == RES_OK);
@@ -78,3 +97,4 @@ main(int argc, char** argv)
   CHK(mem_allocated_size() == 0);
   return 0;
 }
+