commit 12cca4176172fe57a74347b88829328dab25f8f5
parent 6b3bda0df41c824a9f5e1c2ac5eec7e65ae8c178
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Mon, 18 May 2020 19:11:10 +0200
Refactor htrdr_ran_lw_sample and accelerate it
Remove the now useless sampled_band_bounds argument.
Improve performances in "discrete" mode by sampling uniformly in the
sampled band instead of sampling the wavelength wrt the Planck function.
Diffstat:
5 files changed, 11 insertions(+), 52 deletions(-)
diff --git a/src/htrdr_compute_radiance_lw.c b/src/htrdr_compute_radiance_lw.c
@@ -139,7 +139,6 @@ htrdr_compute_radiance_lw
const double pos_in[3],
const double dir_in[3],
const double wlen, /* In nanometer */
- const double samp_band_bounds[2], /* In nanometer */
const size_t iband,
const size_t iquad)
{
@@ -154,22 +153,12 @@ htrdr_compute_radiance_lw
double range[2];
double pos_next[3];
double dir_next[3];
- double samp_band_bounds_m[2]; /* sub-interval bounds in meters */
- double delta_wlen_m; /* sub-interval size in meters */
double temperature;
+ double wlen_m = wlen * 1.e-9;
double g;
double w = 0; /* Weight */
ASSERT(htrdr && rng && pos_in && dir_in && ithread < htrdr->nthreads);
- ASSERT(samp_band_bounds);
-
- /* Convert to meters */
- samp_band_bounds_m[0] = samp_band_bounds[0] * 1e-9;
- samp_band_bounds_m[1] = samp_band_bounds[1] * 1e-9;
-
- /* If monochromatic, set delta_wlen to 1 m, otherwise lmax - lmin in meters */
- delta_wlen_m = samp_band_bounds_m[1] - samp_band_bounds_m[0];
- if(delta_wlen_m == 0) delta_wlen_m = 1.0;
/* Setup the phase function for this spectral band & quadrature point */
CHK(RES_OK == ssf_phase_create
@@ -224,8 +213,7 @@ htrdr_compute_radiance_lw
ASSERT(!SVX_HIT_NONE(&svx_hit));
temperature = htsky_fetch_temperature(htrdr->sky, pos_next);
/* weight is planck integrated over the spectral sub-interval */
- w = planck(samp_band_bounds_m[0], samp_band_bounds_m[1], temperature);
- w = w * delta_wlen_m;
+ w = planck_monochromatic(wlen_m, temperature);
break;
}
@@ -274,8 +262,7 @@ htrdr_compute_radiance_lw
w = 0;
} else {
/* weight is planck integrated over the spectral sub-interval */
- w = planck(samp_band_bounds_m[0], samp_band_bounds_m[1], temperature);
- w = w * delta_wlen_m;
+ w = planck_monochromatic(wlen_m, temperature);
}
break;
}
diff --git a/src/htrdr_draw_radiance.c b/src/htrdr_draw_radiance.c
@@ -605,7 +605,6 @@ draw_pixel_lw
size_t iband;
size_t iquad;
double usec;
- double samp_band_bounds[2];
double band_pdf;
/* Begin the registration of the time spent in the realisation */
@@ -624,7 +623,7 @@ draw_pixel_lw
r2 = ssp_rng_canonical(rng);
/* Sample a wavelength */
- wlen = htrdr_ran_lw_sample(htrdr->ran_lw, r0, r1, samp_band_bounds, &band_pdf);
+ wlen = htrdr_ran_lw_sample(htrdr->ran_lw, r0, r1, &band_pdf);
/* Select the associated band and sample a quadrature point */
iband = htsky_find_spectral_band(htrdr->sky, wlen);
@@ -632,7 +631,7 @@ draw_pixel_lw
/* Compute the integrated luminance in W/m^2/sr */
weight = htrdr_compute_radiance_lw(htrdr, ithread, rng, ray_org, ray_dir,
- wlen, samp_band_bounds, iband, iquad);
+ wlen, iband, iquad);
/* Importance sampling: correct weight with pdf */
weight /= band_pdf;
diff --git a/src/htrdr_ran_lw.c b/src/htrdr_ran_lw.c
@@ -158,7 +158,6 @@ ran_lw_sample_continue
const double r,
const double range[2], /* In nanometer */
const char* func_name,
- double sampled_band_bounds[2],
double* pdf)
{
/* Numerical parameters of the dichotomy search */
@@ -220,12 +219,6 @@ ran_lw_sample_continue
func_name, SPLIT2(range), ran_lw->ref_temperature);
}
- if(sampled_band_bounds) {
- const double lambda = lambda_m * 1.e+9;
- sampled_band_bounds[0] = lambda;
- sampled_band_bounds[1] = lambda;
- }
-
if(pdf) {
const double Tref = ran_lw->ref_temperature;
const double B_lambda = planck(lambda_m, lambda_m, Tref);
@@ -242,7 +235,6 @@ ran_lw_sample_discrete
const double r0,
const double r1,
const char* func_name,
- double sampled_band_bounds[2],
double* pdf)
{
const double* cdf = NULL;
@@ -278,23 +270,12 @@ ran_lw_sample_discrete
/* Fetch the pdf of the sampled band */
pdf_band = darray_double_cdata_get(&ran_lw->pdf)[i];
- /* Continously sample a wavelength in the sampled band */
- lambda = ran_lw_sample_continue
- (ran_lw, r1, band_range, func_name, sampled_band_bounds, &pdf_continue);
+ /* Uniformly sample a wavelength in the sampled band */
+ lambda = band_range[0] + (band_range[1] - band_range[0]) * r1;
+ pdf_continue = 1.0 / ((band_range[1] - band_range[0])*1.e-9);
if(pdf) {
- const double Tref = ran_lw->ref_temperature;
- const double lambda_m = lambda * 1.e-9;
- const double B_lambda = planck(lambda_m, lambda_m, Tref);
- double range_m[2];
- double B_mean;
-
- range_m[0] = ran_lw->range[0] * 1.e-9;
- range_m[1] = ran_lw->range[1] * 1.e-9;
-
- B_mean = planck(range_m[0], range_m[1], Tref);
- *pdf = B_lambda / (B_mean * (range_m[1]-range_m[0]));
- ASSERT(eq_eps(*pdf, pdf_continue * pdf_band, 1e-6));
+ *pdf = pdf_band * pdf_continue;
}
return lambda;
@@ -397,23 +378,17 @@ htrdr_ran_lw_sample
(const struct htrdr_ran_lw* ran_lw,
const double r0,
const double r1,
- double sampled_band_bounds[2],
double* pdf)
{
ASSERT(ran_lw);
if(ran_lw->nbands != HTRDR_RAN_LW_CONTINUE) { /* Discrete */
- return ran_lw_sample_discrete
- (ran_lw, r0, r1, FUNC_NAME, sampled_band_bounds, pdf);
+ return ran_lw_sample_discrete(ran_lw, r0, r1, FUNC_NAME, pdf);
} else if(eq_eps(ran_lw->range[0], ran_lw->range[1], 1.e-6)) {
- if(sampled_band_bounds) { /* Monochromatic */
- sampled_band_bounds[0] = ran_lw->range[0];
- sampled_band_bounds[1] = ran_lw->range[0];
- }
if(pdf) *pdf = 1;
return ran_lw->range[0];
} else { /* Continue */
return ran_lw_sample_continue
- (ran_lw, r0, ran_lw->range, FUNC_NAME, sampled_band_bounds, pdf);
+ (ran_lw, r0, ran_lw->range, FUNC_NAME, pdf);
}
}
diff --git a/src/htrdr_ran_lw.h b/src/htrdr_ran_lw.h
@@ -47,7 +47,6 @@ htrdr_ran_lw_sample
(const struct htrdr_ran_lw* ran_lw,
const double r0, /* Canonical number in [0, 1[ */
const double r1, /* Canonical number in [0, 1[ */
- double sampled_band_bounds[2], /* May be NULL */
double* pdf); /* May be NULL */
#endif /* HTRDR_RAN_LW_H */
diff --git a/src/htrdr_solve.h b/src/htrdr_solve.h
@@ -75,7 +75,6 @@ htrdr_compute_radiance_lw
const double pos[3],
const double dir[3],
const double wlen, /* In nanometer */
- const double samp_band_bounds[2], /* In nanometer */
const size_t iband, /* Index of the spectral band */
const size_t iquad); /* Index of the quadrature point into the band */
