commit bddd7fdd86dfbbeb73ef43b974ef6384538ce0ff
parent 9fb45f8fb33244c05c12f37a110a09249b790b48
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Wed, 17 Apr 2019 09:54:31 +0200
Estimate the time per realisation for each pixel
Diffstat:
3 files changed, 78 insertions(+), 31 deletions(-)
diff --git a/src/htrdr.c b/src/htrdr.c
@@ -102,6 +102,7 @@ static res_T
dump_accum_buffer
(struct htrdr* htrdr,
struct htrdr_buffer* buf,
+ struct htrdr_accum* time_acc, /* May be NULL */
const char* stream_name,
FILE* stream)
{
@@ -112,34 +113,34 @@ dump_accum_buffer
(void)stream_name;
htrdr_buffer_get_layout(buf, &layout);
- if(layout.elmt_size != sizeof(struct htrdr_accum[3])/*#channels*/
- || layout.alignment < ALIGNOF(struct htrdr_accum[3])) {
+ if(layout.elmt_size != sizeof(struct htrdr_accum[4])/*#channels*/
+ || layout.alignment < ALIGNOF(struct htrdr_accum[4])) {
htrdr_log_err(htrdr,
"%s: invalid buffer layout. "
- "The pixel size must be the size of an accumulator.\n",
+ "The pixel size must be the size of 4 accumulators.\n",
FUNC_NAME);
res = RES_BAD_ARG;
goto error;
}
fprintf(stream, "%lu %lu\n", layout.width, layout.height);
+
+ if(time_acc) *time_acc = HTRDR_ACCUM_NULL;
FOR_EACH(y, 0, layout.height) {
FOR_EACH(x, 0, layout.width) {
const struct htrdr_accum* accums = htrdr_buffer_at(buf, x, y);
int i;
- FOR_EACH(i, 0, 3) {
- const double N = (double)accums[i].nweights;
- double E = 0;
- double V = 0;
- double SE = 0;
-
- if(accums[i].nweights) {
- E = accums[i].sum_weights / N;
- V = MMAX(accums[i].sum_weights_sqr / N - E*E, 0);
- SE = sqrt(V/N);
- }
+ FOR_EACH(i, 0, 4) {
+ double E, SE;
+
+ htrdr_accum_get_estimation(&accums[i], &E, &SE);
fprintf(stream, "%g %g ", E, SE);
}
+ if(time_acc) {
+ time_acc->sum_weights += accums[3].sum_weights;
+ time_acc->sum_weights_sqr += accums[3].sum_weights_sqr;
+ time_acc->nweights += accums[3].nweights;
+ }
fprintf(stream, "\n");
}
fprintf(stream, "\n");
@@ -527,12 +528,15 @@ htrdr_init
/* Create the image buffer only on the master process; the image parts
* rendered by the processes are gathered onto the master process. */
if(!htrdr->dump_vtk && htrdr->mpi_rank == 0) {
+ /* 4 accums: X, Y, Z components and one more for the per realisation time */
+ const size_t pixsz = sizeof(struct htrdr_accum[4]);
+ const size_t pixal = ALIGNOF(struct htrdr_accum[4]);
res = htrdr_buffer_create(htrdr,
args->image.definition[0], /* Width */
args->image.definition[1], /* Height */
- args->image.definition[0]*sizeof(struct htrdr_accum[3]), /* Pitch */
- sizeof(struct htrdr_accum[3]),
- ALIGNOF(struct htrdr_accum[3]), /* Alignment */
+ args->image.definition[0] * pixsz, /* Pitch */
+ pixsz, /* Size of a pixel */
+ pixal, /* Alignment of a pixel */
&htrdr->buf);
if(res != RES_OK) goto error;
}
@@ -629,9 +633,16 @@ htrdr_run(struct htrdr* htrdr)
htrdr->height, htrdr->spp, htrdr->buf);
if(res != RES_OK) goto error;
if(htrdr->mpi_rank == 0) {
- res = dump_accum_buffer
- (htrdr, htrdr->buf, str_cget(&htrdr->output_name), htrdr->output);
+ struct htrdr_accum path_time_acc = HTRDR_ACCUM_NULL;
+ double E, SE;
+
+ res = dump_accum_buffer(htrdr, htrdr->buf, &path_time_acc,
+ str_cget(&htrdr->output_name), htrdr->output);
if(res != RES_OK) goto error;
+
+ htrdr_accum_get_estimation(&path_time_acc, &E, &SE);
+ htrdr_log(htrdr,
+ "Time per radiative path (in micro seconds): %g +/- %g\n", E, SE);
}
}
exit:
diff --git a/src/htrdr_draw_radiance_sw.c b/src/htrdr_draw_radiance_sw.c
@@ -92,7 +92,7 @@ tile_create(struct mem_allocator* allocator)
const size_t tile_sz =
sizeof(struct tile) - sizeof(struct htrdr_accum)/*rm dummy accum*/;
const size_t buf_sz = /* Flexiblbe array element */
- TILE_SIZE*TILE_SIZE*sizeof(struct htrdr_accum)*3/*#channels*/;
+ TILE_SIZE*TILE_SIZE*sizeof(struct htrdr_accum[4/*#channels*/]);
ASSERT(allocator);
tile = MEM_ALLOC(allocator, tile_sz+buf_sz);
@@ -135,7 +135,7 @@ tile_at
const size_t y) /* In tile space */
{
ASSERT(tile && x < TILE_SIZE && y < TILE_SIZE);
- return tile->data.accums + (y*TILE_SIZE + x)*3/*#channels*/;
+ return tile->data.accums + (y*TILE_SIZE + x)*4/*#channels*/;
}
static void
@@ -150,6 +150,8 @@ write_tile_data(struct htrdr_buffer* buf, const struct tile_data* tile_data)
htrdr_buffer_get_layout(buf, &layout);
buf_mem = htrdr_buffer_get_data(buf);
+ ASSERT(layout.elmt_size == sizeof(struct htrdr_accum[4]));
+
/* Compute the row/column of the tile origin into the buffer */
icol = tile_data->x * (size_t)TILE_SIZE;
@@ -163,9 +165,9 @@ write_tile_data(struct htrdr_buffer* buf, const struct tile_data* tile_data)
FOR_EACH(irow_tile, 0, nrows_tile) {
char* buf_row = buf_mem + (irow + irow_tile) * layout.pitch;
const struct htrdr_accum* tile_row =
- tile_data->accums + irow_tile*TILE_SIZE*3/*#channels*/;
- memcpy(buf_row + icol*sizeof(struct htrdr_accum)*3, tile_row,
- ncols_tile*sizeof(struct htrdr_accum)*3/*#channels*/);
+ tile_data->accums + irow_tile*TILE_SIZE*4/*#channels*/;
+ memcpy(buf_row + icol*sizeof(struct htrdr_accum[4/*#channels*/]), tile_row,
+ ncols_tile*sizeof(struct htrdr_accum[4/*#channels*/]));
}
}
@@ -391,7 +393,7 @@ mpi_gather_tiles
/* Compute the size of the tile_data */
const size_t msg_sz =
sizeof(struct tile_data) - sizeof(struct htrdr_accum)/*dummy*/
- + TILE_SIZE*TILE_SIZE*sizeof(struct htrdr_accum)*3/*#channels*/;
+ + TILE_SIZE*TILE_SIZE*sizeof(struct htrdr_accum[4/*#channels*/]);
struct list_node* node = NULL;
struct tile* tile = NULL;
@@ -480,7 +482,8 @@ draw_tile
/* Fetch and reset the pixel accumulator */
pix_accums = tile_at(tile, ipix_tile[0], ipix_tile[1]);
-
+ pix_accums[3] = HTRDR_ACCUM_NULL; /* Reset time per radiative path */
+
/* Compute the pixel coordinate */
ipix[0] = tile_org[0] + ipix_tile[0];
ipix[1] = tile_org[1] + ipix_tile[1];
@@ -490,12 +493,14 @@ draw_tile
pix_accums[ichannel] = HTRDR_ACCUM_NULL;
FOR_EACH(isamp, 0, spp) {
+ struct time t0, t1;
double pix_samp[2];
double ray_org[3];
double ray_dir[3];
double weight;
size_t iband;
size_t iquad;
+ double usec;
/* Sample a position into the pixel, in the normalized image plane */
pix_samp[0] = ((double)ipix[0] + ssp_rng_canonical(rng)) * pix_sz[0];
@@ -523,14 +528,22 @@ draw_tile
}
/* Compute the radiance that reach the pixel through the ray */
+ time_current(&t0);
weight = htrdr_compute_radiance_sw
(htrdr, ithread, rng, ray_org, ray_dir, iband, iquad);
ASSERT(weight >= 0);
+ time_sub(&t0, time_current(&t1), &t0);
+ usec = (double)time_val(&t0, TIME_NSEC) * 0.001;
- /* Update the pixel accumulator */
+ /* Update the pixel accumulator of the current channel */
pix_accums[ichannel].sum_weights += weight;
pix_accums[ichannel].sum_weights_sqr += weight*weight;
pix_accums[ichannel].nweights += 1;
+
+ /* Update the pixel accumulator of per realisation time */
+ pix_accums[3].sum_weights += usec;
+ pix_accums[3].sum_weights_sqr += usec*usec;
+ pix_accums[3].nweights += 1;
}
}
}
@@ -724,11 +737,11 @@ htrdr_draw_radiance_sw
ASSERT(layout.width || layout.height || layout.elmt_size);
ASSERT(layout.width == width && layout.height == height);
- if(layout.elmt_size != sizeof(struct htrdr_accum[3])/*#channels*/
- || layout.alignment < ALIGNOF(struct htrdr_accum[3])) {
+ if(layout.elmt_size != sizeof(struct htrdr_accum[4])/*#channels*/
+ || layout.alignment < ALIGNOF(struct htrdr_accum[4])) {
htrdr_log_err(htrdr,
"%s: invalid buffer layout. "
- "The pixel size must be the size of 3 * accumulators.\n",
+ "The pixel size must be the size of 4 accumulators.\n",
FUNC_NAME);
res = RES_BAD_ARG;
goto error;
diff --git a/src/htrdr_solve.h b/src/htrdr_solve.h
@@ -51,7 +51,7 @@ htrdr_draw_radiance_sw
const size_t width, /* Image width */
const size_t height, /* Image height */
const size_t spp, /* #samples per pixel, i.e. #realisations */
- /* Buffer of struct htrdr_accum[3]. May be NULL on on non master processes */
+ /* Buffer of struct htrdr_accum[4]. May be NULL on non master processes */
struct htrdr_buffer* buf);
extern LOCAL_SYM int
@@ -62,4 +62,27 @@ htrdr_ground_filter
void* ray_data,
void* filter_data);
+static FINLINE void
+htrdr_accum_get_estimation
+ (const struct htrdr_accum* acc,
+ double* expected_value,
+ double* std_err)
+{
+ ASSERT(acc && expected_value && std_err);
+
+ if(!acc->nweights) {
+ *expected_value = 0;
+ *std_err = 0;
+ } else {
+ const double N = (double)acc->nweights;
+ double E, V, SE;
+ E = acc->sum_weights / N;
+ V = MMAX(acc->sum_weights_sqr / N - E*E, 0);
+ SE = sqrt(V/N);
+
+ *expected_value = E;
+ *std_err = SE;
+ }
+}
+
#endif /* HTRDR_SOLVE_H */
