commit 90e03015f2c2405d5af63b4cc6f68d7a61660ff3
parent c12efa8c9aff11e3c5e0cbf79ec821632e53ca36
Author: Vincent Forest <vincent.forest@meso-star.com>
Date: Mon, 31 Oct 2022 14:53:36 +0100
htrdr-planeto: add a dummy compute radiance function
Start debugging the entire building
Diffstat:
5 files changed, 181 insertions(+), 12 deletions(-)
diff --git a/cmake/planeto/CMakeLists.txt b/cmake/planeto/CMakeLists.txt
@@ -53,6 +53,7 @@ include_directories(
set(HTRDR_PLANETO_FILES_SRC
htrdr_planeto.c
htrdr_planeto_args.c
+ htrdr_planeto_compute_radiance.c
htrdr_planeto_draw_map.c
htrdr_planeto_main.c
htrdr_planeto_source.c)
diff --git a/src/planeto/htrdr_planeto.c b/src/planeto/htrdr_planeto.c
@@ -67,11 +67,13 @@ compute_min_band_len(const struct htrdr_planeto* cmd)
/* At least one band must be overlaped by the spectral domain */
ASSERT(ibands[0]<=ibands[1]);
FOR_EACH(i, ibands[0], ibands[1]+1) {
+ struct rnatm_band_desc band;
double band_range[2];
- RNATM(band_get_range(cmd->atmosphere, i, band_range));
+ RNATM(band_get_desc(cmd->atmosphere, i, &band));
/* Make the upper bound inclusive */
- band_range[1] = nextafter(band_range[1], 0);
+ band_range[0] = band.lower;
+ band_range[1] = nextafter(band.upper, 0);
/* Clamp the band range to the spectral domain */
band_range[0] = MMAX(band_range[0], range[0]);
@@ -535,7 +537,7 @@ htrdr_planeto_run(struct htrdr_planeto* cmd)
switch(cmd->output_type) {
case HTRDR_PLANETO_ARGS_OUTPUT_IMAGE:
- htrdr_log_warn(cmd->htrdr, "image rendering is not yet implemented\n");
+ res = planeto_draw_map(cmd);
break;
case HTRDR_PLANETO_ARGS_OUTPUT_OCTREES:
res = write_vtk_octrees(cmd);
diff --git a/src/planeto/htrdr_planeto_c.h b/src/planeto/htrdr_planeto_c.h
@@ -62,7 +62,20 @@ struct planeto_pixel_image {
HTRDR_ACCUM_NULL__ \
}
-struct planeto_pixel_xwave;
+struct planeto_compute_radiance_args {
+ struct ssp_rng* rng;
+ size_t ithread; /* Index of the thread executing the function */
+
+ double path_org[3]; /* Origin of the path to trace */
+ double path_dir[3]; /* Initial direction of the path to trace */
+
+ double wlen; /* In nm */
+ size_t iband; /* Spectral band index */
+ size_t iquad; /* Quadrature point */
+};
+#define PLANETO_COMPUTE_RADIANCE_ARGS_NULL__ {NULL, 0, {0,0,0}, {0,0,0}, 0, 0, 0}
+static const struct planeto_compute_radiance_args
+PLANETO_COMPUTE_RADIANCE_ARGS_NULL = PLANETO_COMPUTE_RADIANCE_ARGS_NULL__;
struct htrdr_planeto {
struct rnatm* atmosphere;
@@ -98,4 +111,10 @@ planeto_get_pixel_format
(const struct htrdr_planeto* cmd,
struct htrdr_pixel_format* fmt);
+/* Return the radiance in W/m²/sr/m */
+extern LOCAL_SYM double
+planeto_compute_radiance
+ (struct htrdr_planeto* cmd,
+ const struct planeto_compute_radiance_args* args);
+
#endif /* HTRDR_PLANETO_C_H */
diff --git a/src/planeto/htrdr_planeto_compute_radiance.c b/src/planeto/htrdr_planeto_compute_radiance.c
@@ -0,0 +1,142 @@
+/* Copyright (C) 2018, 2019, 2020, 2021 |Meso|Star> (contact@meso-star.com)
+ * Copyright (C) 2018, 2019, 2021 CNRS
+ * Copyright (C) 2018, 2019, Université Paul Sabatier
+ *
+ * This program is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. */
+
+#include "planeto/htrdr_planeto_c.h"
+#include "planeto/htrdr_planeto_source.h"
+
+#include <rad-net/rnatm.h>
+#include <rad-net/rngrd.h>
+
+#include <star/s3d.h>
+#include <star/ssf.h>
+#include <star/ssp.h>
+
+#include <rsys/double2.h>
+#include <rsys/double3.h>
+
+/*******************************************************************************
+ * Helper functions
+ ******************************************************************************/
+static INLINE res_T
+check_planeto_compute_radiance_args
+ (const struct htrdr_planeto* cmd,
+ const struct planeto_compute_radiance_args* args)
+{
+ struct rnatm_band_desc band = RNATM_BAND_DESC_NULL;
+ res_T res = RES_OK;
+
+ if(!args || !args->rng)
+ return RES_BAD_ARG;
+
+ /* Invalid thread index */
+ if(args->ithread >= htrdr_get_threads_count(cmd->htrdr))
+ return RES_BAD_ARG;
+
+ /* Invalid input direction */
+ if(!d3_is_normalized(args->path_dir))
+ return RES_BAD_ARG;
+
+ /* Invalid wavelength */
+ if(args->wlen < cmd->spectral_domain.wlen_range[0]
+ || args->wlen > cmd->spectral_domain.wlen_range[1])
+ return RES_BAD_ARG;
+
+ res = rnatm_band_get_desc(cmd->atmosphere, args->iband, &band);
+ if(res != RES_OK) return res;
+
+ /* Inconsistent spectral dimension */
+ if(args->wlen < band.lower
+ || args->wlen >= band.upper /* Exclusive */
+ || args->iquad >= band.quad_pts_count)
+ return RES_BAD_ARG;
+
+ return RES_OK;
+}
+
+/*******************************************************************************
+ * Local functions
+ ******************************************************************************/
+double
+planeto_compute_radiance
+ (struct htrdr_planeto* cmd,
+ const struct planeto_compute_radiance_args* args)
+{
+ struct rngrd_trace_ray_args rt_args = RNGRD_TRACE_RAY_ARGS_DEFAULT;
+ struct rngrd_create_bsdf_args bsdf_args = RNGRD_CREATE_BSDF_ARGS_NULL;
+ struct s3d_hit s3d_hit0 = S3D_HIT_NULL;
+ struct s3d_hit s3d_hit1 = S3D_HIT_NULL;
+ struct ssf_bsdf* bsdf = NULL;
+ double pos[3];
+ double N[3];
+ double wo[3];
+ double wi[3];
+ double L = 0;
+ double Li = 0;
+ double cos_wi_N = 0;
+ double rho = 0;
+ ASSERT(cmd && check_planeto_compute_radiance_args(cmd, args) == RES_OK);
+
+ /* Look for a surface intersection */
+ d3_set(rt_args.ray_org, args->path_org);
+ d3_set(rt_args.ray_dir, args->path_dir);
+ d2(rt_args.ray_range, 0, INF);
+ RNGRD(trace_ray(cmd->ground, &rt_args, &s3d_hit0));
+
+ /* No surface intersection */
+ if(S3D_HIT_NONE(&s3d_hit0)) goto exit; /* Nothing more to do */
+
+ /* Calculate the position of intersection */
+ pos[0] = s3d_hit0.distance * rt_args.ray_dir[0] + rt_args.ray_org[0];
+ pos[1] = s3d_hit0.distance * rt_args.ray_dir[1] + rt_args.ray_org[1];
+ pos[2] = s3d_hit0.distance * rt_args.ray_dir[2] + rt_args.ray_org[2];
+
+ /* Sample a direction toward the source */
+ htrdr_planeto_source_sample_direction(cmd->source, args->rng, pos, wi);
+ Li = htrdr_planeto_source_get_radiance(cmd->source, args->wlen);
+
+ /* The source is not facing the current position */
+ d3_normalize(N, d3_set_f3(N, s3d_hit0.normal));
+ if(d3_dot(N, rt_args.ray_dir) > 0) d3_minus(N, N);
+ cos_wi_N = d3_dot(wi, N);
+ if(cos_wi_N < 0) goto exit; /* Nothing more to do */
+
+ /* Check the source visibility */
+ d3_set(rt_args.ray_org, pos);
+ d3_set(rt_args.ray_dir, wi);
+ d2(rt_args.ray_range, 0, INF);
+ rt_args.hit_from = s3d_hit0;
+ RNGRD(trace_ray(cmd->ground, &rt_args, &s3d_hit1));
+ if(!S3D_HIT_NONE(&s3d_hit1)) goto exit; /* Nothing more to do */
+
+ /* Fetch the surface BSDF */
+ bsdf_args.prim = s3d_hit0.prim;
+ bsdf_args.barycentric_coords[0] = s3d_hit0.uv[0];
+ bsdf_args.barycentric_coords[1] = s3d_hit0.uv[1];
+ bsdf_args.barycentric_coords[2] = 1 - s3d_hit0.uv[0] - s3d_hit0.uv[1];
+ bsdf_args.wavelength = args->wlen;
+ bsdf_args.r = ssp_rng_canonical(args->rng);
+ RNGRD(create_bsdf(cmd->ground, &bsdf_args, &bsdf));
+
+ d3_minus(wo, rt_args.ray_org);
+ rho = ssf_bsdf_eval(bsdf, wo, N, wi);
+ SSF(bsdf_ref_put(bsdf));
+
+ L = rho * Li * cos_wi_N;
+
+exit:
+ return L;
+}
diff --git a/src/planeto/htrdr_planeto_draw_map.c b/src/planeto/htrdr_planeto_draw_map.c
@@ -49,6 +49,9 @@ draw_pixel_image
const struct htrdr_draw_pixel_args* args,
void* data)
{
+ struct planeto_compute_radiance_args rad_args =
+ PLANETO_COMPUTE_RADIANCE_ARGS_NULL;
+
struct htrdr_accum XYZ[3]; /* X, Y, and Z */
struct htrdr_accum time;
struct htrdr_planeto* cmd;
@@ -95,7 +98,7 @@ draw_pixel_image
sample.lens[1] = ssp_rng_canonical(args->rng);
/* Generate a camera ray */
- scam_generate_ray(cmd->camera, &sample, &ray);
+ SCAM(generate_ray(cmd->camera, &sample, &ray));
r0 = ssp_rng_canonical(args->rng);
r1 = ssp_rng_canonical(args->rng);
@@ -122,13 +125,15 @@ draw_pixel_image
RNATM(band_sample_quad_pt(cmd->atmosphere, r2, iband[0], &iquad));
ASSERT(iband[0] == iband[1]);
- /* TODO Compute the radiance in W/m²/sr/m */
-#if 0
- weight = planeto_compute_radiance(cmd, args->ithread, args->rng,
- ray.org, ray.dir, wlen, iband, iquad);
-#else
- weight = 0;
-#endif
+ /* Compute the radiance in W/m²/sr/m */
+ d3_set(rad_args.path_org, ray.org);
+ d3_set(rad_args.path_dir, ray.dir);
+ rad_args.rng = args->rng;
+ rad_args.ithread = args->ithread;
+ rad_args.wlen = wlen[0];
+ rad_args.iband = iband[0];
+ rad_args.iquad = iquad;
+ weight = planeto_compute_radiance(cmd, &rad_args);
ASSERT(weight >= 0);
weight /= pdf; /* In W/m²/sr */
