htrdr

Solving radiative transfer in heterogeneous media
git clone git://git.meso-star.fr/htrdr.git
Log | Files | Refs | README | LICENSE

commit 35c1cbf09445778096b7f0e54c16b6cec9949de9
parent cc419c395bf5eae2cef4b6805b8255e9bb52b6fd
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Fri,  6 Dec 2019 15:32:46 +0100

Merge branch 'release_0.2'

Diffstat:
MREADME.md | 7+++++++
Mcmake/CMakeLists.txt | 9++++-----
Mdoc/htrdr.1.txt.in | 14+++++++++++---
Msrc/htrdr.c | 2+-
Msrc/htrdr_args.c | 43+++++++++++++++++++++++++++++++++++++++++--
Msrc/htrdr_args.h.in | 4++++
Msrc/htrdr_compute_radiance_sw.c | 13+++++++------
Msrc/htrdr_ground.c | 102++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++-------
Msrc/htrdr_ground.h | 13++++++++++---
Msrc/htrdr_sky.c | 46++++++++++++++++++++++++++--------------------
10 files changed, 205 insertions(+), 48 deletions(-)

diff --git a/README.md b/README.md @@ -60,6 +60,13 @@ informations on CMake. ## Release notes +### Version 0.2 + +- Add the `-b` option that controls the BRDF of the ground geometry. +- Make optional the use of a ground geometry (option `-g`). +- Make optional the definition of the optical properties of water droplets + (option `-m`) when no cloud field is used. + ### Version 0.1 - Add the `-V` option that fixes the maximum definition of the octrees used to diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -52,6 +52,7 @@ include_directories( ${HTCP_INCLUDE_DIR} ${HTGOP_INCLUDE_DIR} ${HTMIE_INCLUDE_DIR} + ${HTRDR_SOURCE_DIR} ${MPI_INCLUDE_PATH} ${CMAKE_CURRENT_BINARY_DIR}) @@ -59,7 +60,7 @@ include_directories( # Generate files ################################################################################ set(VERSION_MAJOR 0) -set(VERSION_MINOR 1) +set(VERSION_MINOR 2) set(VERSION_PATCH 0) set(VERSION ${VERSION_MAJOR}.${VERSION_MINOR}.${VERSION_PATCH}) @@ -67,15 +68,13 @@ set(HTRDR_ARGS_DEFAULT_CAMERA_POS "0,0,0") set(HTRDR_ARGS_DEFAULT_CAMERA_TGT "0,1,0") set(HTRDR_ARGS_DEFAULT_CAMERA_UP "0,0,1") set(HTRDR_ARGS_DEFAULT_CAMERA_FOV "70") +set(HTRDR_ARGS_DEFAULT_GROUND_BSDF "HTRDR_BSDF_DIFFUSE") +set(HTRDR_ARGS_DEFAULT_GROUND_REFLECTIVITY "0.5") set(HTRDR_ARGS_DEFAULT_IMG_WIDTH "320") set(HTRDR_ARGS_DEFAULT_IMG_HEIGHT "240") set(HTRDR_ARGS_DEFAULT_IMG_SPP "1") -set(HTRDR_ARGS_DEFAULT_GROUND_REFLECTIVITY "0.5") set(HTRDR_ARGS_DEFAULT_OPTICAL_THICKNESS_THRESHOLD "1") -configure_file(${HTRDR_SOURCE_DIR}/htrdr_version.h.in - ${CMAKE_CURRENT_BINARY_DIR}/htrdr_version.h @ONLY) - configure_file(${HTRDR_SOURCE_DIR}/htrdr_args.h.in ${CMAKE_CURRENT_BINARY_DIR}/htrdr_args.h @ONLY) diff --git a/doc/htrdr.1.txt.in b/doc/htrdr.1.txt.in @@ -72,9 +72,16 @@ OPTIONS Path toward the file containing the gas optical properties of the atmosphere. Data must be formatted according to the fileformat described in [1]. +*-b* <diffuse|specular>:: + Define how light is reflected at the ground surface. When set to diffuse, the + reflections are lambertian while a specular surface simulates a perfect + mirror. By default, the ground surface is diffuse. Note that this option + controls only the directional part of the surface reflections; the amount of + reflected energy is defined by the *-e* option. + *-c* _clouds_:: Submit a *htcp*(5) file describing the properties of the clouds. If not - defined, only the atmopshere properties submitted through the *-a* option + defined, only the atmosphere properties submitted through the *-a* option are taken into account. *-C* <__camera-parameter__:...>:: @@ -107,7 +114,7 @@ OPTIONS {0,0,+1} regardless of _azimuth_ value. *-d*:: - Write in _output_ the space partitionning data structures used to speed up + Write in _output_ the space partitioning data structures used to speed up the radiative transfer computations in the clouds. The written data are octrees saved in the VTK file format [3]. Each octree node stores the minimum and the maximum of the extinction coefficients of the cloud cells overlapped @@ -117,7 +124,8 @@ OPTIONS *-e* _reflectivity_:: Reflectivity of the ground geometry in [0, 1]. By default it is set to @HTRDR_ARGS_DEFAULT_GROUND_REFLECTIVITY@. This parameter is fixed for the - while visible range. + while visible range and defines the amount of reflected energy. Refer to the + *-b* option to control how the light is reflected on the ground surface. *-f*:: Force overwrite of the _output_ file. diff --git a/src/htrdr.c b/src/htrdr.c @@ -514,7 +514,7 @@ htrdr_init goto error; } - res = htrdr_ground_create(htrdr, args->filename_obj, + res = htrdr_ground_create(htrdr, args->filename_obj, args->ground_bsdf_type, args->ground_reflectivity, args->repeat_ground, &htrdr->ground); if(res != RES_OK) goto error; diff --git a/src/htrdr_args.c b/src/htrdr_args.c @@ -27,6 +27,18 @@ /******************************************************************************* * Helper functions ******************************************************************************/ +static const char* +bsdf_type_to_string(const enum htrdr_bsdf_type type) +{ + const char* str = "<none>"; + switch(type) { + case HTRDR_BSDF_DIFFUSE: str = "diffuse"; break; + case HTRDR_BSDF_SPECULAR: str = "specular"; break; + default: FATAL("Unreachable code.\n"); break; + } + return str; +} + static void print_help(const char* cmd) { @@ -38,6 +50,10 @@ print_help(const char* cmd) printf( " -a ATMOSPHERE gas optical properties of the atmosphere.\n"); printf( +" -b <diffuse|specular>\n" +" BSDF of the ground. Default value is %s.\n", + bsdf_type_to_string(HTRDR_ARGS_DEFAULT.ground_bsdf_type)); + printf( " -c CLOUDS properties of the clouds.\n"); printf( " -C <camera> define the rendering point of view.\n"); @@ -354,6 +370,26 @@ error: goto exit; } +static res_T +parse_bsdf_type(struct htrdr_args* args, const char* str) +{ + res_T res = RES_OK; + if(!strcmp(str, "diffuse")) { + args->ground_bsdf_type = HTRDR_BSDF_DIFFUSE; + } else if(!strcmp(str, "specular")) { + args->ground_bsdf_type = HTRDR_BSDF_SPECULAR; + } else { + fprintf(stderr, "Invalid BRDF type `%s'.\n", str); + res = RES_BAD_ARG; + goto error; + } + +exit: + return res; +error: + goto exit; +} + /******************************************************************************* * Local functions ******************************************************************************/ @@ -378,9 +414,12 @@ htrdr_args_init(struct htrdr_args* args, int argc, char** argv) } } - while((opt = getopt(argc, argv, "a:C:c:D:de:fGg:hi:m:o:RrT:t:V:v")) != -1) { + while((opt = getopt(argc, argv, "a:b:C:c:D:de:fGg:hi:m:o:RrT:t:V:v")) != -1) { switch(opt) { case 'a': args->filename_gas = optarg; break; + case 'b': + res = parse_bsdf_type(args, optarg); + break; case 'C': res = parse_multiple_parameters (args, optarg, parse_camera_parameter); @@ -436,7 +475,7 @@ htrdr_args_init(struct htrdr_args* args, int argc, char** argv) res = RES_BAD_ARG; goto error; } - if(!args->filename_mie) { + if(args->filename_les && !args->filename_mie) { fprintf(stderr, "Missing the path toward the file of the Mie's data -- option '-m'\n"); res = RES_BAD_ARG; diff --git a/src/htrdr_args.h.in b/src/htrdr_args.h.in @@ -16,6 +16,8 @@ #ifndef HTRDR_ARGS_H #define HTRDR_ARGS_H +#include "htrdr_ground.h" + #include <limits.h> #include <rsys/rsys.h> @@ -48,6 +50,7 @@ struct htrdr_args { double sun_azimuth; /* In degrees */ double sun_elevation; /* In degrees */ double optical_thickness; /* Threshold used during octree building */ + enum htrdr_bsdf_type ground_bsdf_type; double ground_reflectivity; /* Reflectivity of the ground */ unsigned grid_max_definition[3]; /* Maximum definition of the grid */ @@ -84,6 +87,7 @@ struct htrdr_args { 0, /* Sun azimuth */ \ 90, /* Sun elevation */ \ @HTRDR_ARGS_DEFAULT_OPTICAL_THICKNESS_THRESHOLD@, /* Optical thickness */ \ + @HTRDR_ARGS_DEFAULT_GROUND_BSDF@, \ @HTRDR_ARGS_DEFAULT_GROUND_REFLECTIVITY@, /* Ground reflectivity */ \ {UINT_MAX, UINT_MAX, UINT_MAX}, /* Maximum definition of the grid */ \ (unsigned)~0, /* #threads */ \ diff --git a/src/htrdr_compute_radiance_sw.c b/src/htrdr_compute_radiance_sw.c @@ -278,21 +278,19 @@ htrdr_compute_radiance_sw ASSERT(htrdr && rng && pos_in && dir_in && ithread < htrdr->nthreads); - CHK(RES_OK == ssf_bsdf_create - (&htrdr->lifo_allocators[ithread], &ssf_lambertian_reflection, &bsdf)); CHK(RES_OK == ssf_phase_create (&htrdr->lifo_allocators[ithread], &ssf_phase_hg, &phase_hg)); CHK(RES_OK == ssf_phase_create (&htrdr->lifo_allocators[ithread], &ssf_phase_rayleigh, &phase_rayleigh)); - SSF(lambertian_reflection_setup - (bsdf, htrdr_ground_get_reflectivity(htrdr->ground))); - /* Setup the phase function for this spectral band & quadrature point */ g = htrdr_sky_fetch_particle_phase_function_asymmetry_parameter (htrdr->sky, iband, iquad); SSF(phase_hg_setup(phase_hg, g)); + /* Fetch the ground BSDF */ + bsdf = htrdr_ground_get_bsdf(htrdr->ground); + /* Fetch sun properties. Arbitrarily use the wavelength at the center of the * band to retrieve the sun radiance of the current band. Note that the sun * spectral data are defined by bands that, actually are the same of the SW @@ -397,6 +395,10 @@ htrdr_compute_radiance_sw bounce_reflectivity = ssf_bsdf_sample (bsdf, rng, wo, N, dir_next, &type, &pdf); + if(!(type & SSF_REFLECTION)) { /* Handle only reflections */ + bounce_reflectivity = 0; + } + if(d3_dot(N, sun_dir) < 0) { /* Below the ground */ R = 0; } else { @@ -439,7 +441,6 @@ htrdr_compute_radiance_sw d3_set(pos, pos_next); d3_set(dir, dir_next); } - SSF(bsdf_ref_put(bsdf)); SSF(phase_ref_put(phase_hg)); SSF(phase_ref_put(phase_rayleigh)); return w; diff --git a/src/htrdr_ground.c b/src/htrdr_ground.c @@ -26,6 +26,7 @@ #include <star/s3d.h> #include <star/s3daw.h> +#include <star/ssf.h> struct ray_context { float range[2]; @@ -48,7 +49,7 @@ struct htrdr_ground { struct s3d_scene_view* view; float lower[3]; /* Ground lower bound */ float upper[3]; /* Ground upper bound */ - double reflectivity; + struct ssf_bsdf* bsdf; int repeat; /* Make the ground infinite in X and Y */ struct htrdr* htrdr; @@ -136,7 +137,15 @@ setup_ground(struct htrdr_ground* ground, const char* obj_filename) size_t nshapes; size_t ishape; res_T res = RES_OK; - ASSERT(ground && obj_filename); + ASSERT(ground); + + if(!obj_filename) { + /* No geometry! Discard the creation of the scene view */ + htrdr_log_warn(ground->htrdr, + "%s: the scene does not have ground geometry.\n", + FUNC_NAME); + goto exit; + } res = s3daw_create(&ground->htrdr->logger, ground->htrdr->allocator, NULL, NULL, ground->htrdr->s3d, ground->htrdr->verbose, &s3daw); @@ -207,6 +216,63 @@ error: goto exit; } +static res_T +setup_bsdf_diffuse(struct htrdr_ground* ground, const double reflectivity) +{ + res_T res = RES_OK; + ASSERT(ground); + + res = ssf_bsdf_create + (ground->htrdr->allocator, &ssf_lambertian_reflection, &ground->bsdf); + if(res != RES_OK) goto error; + + res = ssf_lambertian_reflection_setup(ground->bsdf, reflectivity); + if(res != RES_OK) goto error; + +exit: + return res; +error: + htrdr_log_err(ground->htrdr, + "Could not setup the ground diffuse BSDF with a reflectivity of %g -- %s.\n", + reflectivity, res_to_cstr(res)); + if(ground->bsdf) { + SSF(bsdf_ref_put(ground->bsdf)); + ground->bsdf = NULL; + } + goto exit; +} + +static res_T +setup_bsdf_specular(struct htrdr_ground* ground, const double reflectivity) +{ + struct ssf_fresnel* fresnel = NULL; + res_T res = RES_OK; + ASSERT(ground); + + res = ssf_bsdf_create + (ground->htrdr->allocator, &ssf_specular_reflection, &ground->bsdf); + if(res != RES_OK) goto error; + res = ssf_fresnel_create + (ground->htrdr->allocator, &ssf_fresnel_constant, &fresnel); + if(res != RES_OK) goto error; + res = ssf_fresnel_constant_setup(fresnel, reflectivity); + if(res != RES_OK) goto error; + res = ssf_specular_reflection_setup(ground->bsdf, fresnel); + if(res != RES_OK) goto error; + +exit: + return res; +error: + htrdr_log_err(ground->htrdr, + "Could not setup the ground specular BSDF with a reflectivity of %g -- %s.\n", + reflectivity, res_to_cstr(res)); + if(ground->bsdf) { + SSF(bsdf_ref_put(ground->bsdf)); + ground->bsdf = NULL; + } + goto exit; +} + static void release_ground(ref_T* ref) { @@ -214,6 +280,7 @@ release_ground(ref_T* ref) ASSERT(ref); ground = CONTAINER_OF(ref, struct htrdr_ground, ref); if(ground->view) S3D(scene_view_ref_put(ground->view)); + if(ground->bsdf) SSF(bsdf_ref_put(ground->bsdf)); MEM_RM(ground->htrdr->allocator, ground); } @@ -223,7 +290,8 @@ release_ground(ref_T* ref) res_T htrdr_ground_create (struct htrdr* htrdr, - const char* obj_filename, + const char* obj_filename, /* May be NULL */ + const enum htrdr_bsdf_type bsdf_type, const double reflectivity, const int repeat_ground, /* Infinitely repeat the ground in X and Y */ struct htrdr_ground** out_ground) @@ -232,7 +300,7 @@ htrdr_ground_create struct htrdr_ground* ground = NULL; struct time t0, t1; res_T res = RES_OK; - ASSERT(htrdr && obj_filename && out_ground); + ASSERT(htrdr && out_ground); ASSERT(reflectivity >= 0 || reflectivity <= 1); ground = MEM_CALLOC(htrdr->allocator, 1, sizeof(*ground)); @@ -246,7 +314,20 @@ htrdr_ground_create ref_init(&ground->ref); ground->htrdr = htrdr; ground->repeat = repeat_ground; - ground->reflectivity = reflectivity; + f3_splat(ground->lower, (float)INF); + f3_splat(ground->upper,-(float)INF); + + switch(bsdf_type) { + case HTRDR_BSDF_DIFFUSE: + res = setup_bsdf_diffuse(ground, reflectivity); + break; + case HTRDR_BSDF_SPECULAR: + res = setup_bsdf_specular(ground, reflectivity); + break; + default: FATAL("Unreachable code\n"); + + } + if(res != RES_OK) goto error; time_current(&t0); res = setup_ground(ground, obj_filename); @@ -280,11 +361,11 @@ htrdr_ground_ref_put(struct htrdr_ground* ground) ref_put(&ground->ref, release_ground); } -double -htrdr_ground_get_reflectivity(const struct htrdr_ground* ground) +struct ssf_bsdf* +htrdr_ground_get_bsdf(const struct htrdr_ground* ground) { ASSERT(ground); - return ground->reflectivity; + return ground->bsdf; } res_T @@ -299,6 +380,11 @@ htrdr_ground_trace_ray res_T res = RES_OK; ASSERT(ground && org && dir && range && hit); + if(!ground->view) { /* No ground geometry */ + *hit = S3D_HIT_NULL; + goto exit; + } + if(!ground->repeat) { struct ray_context ray_ctx = RAY_CONTEXT_NULL; float ray_org[3]; diff --git a/src/htrdr_ground.h b/src/htrdr_ground.h @@ -18,15 +18,22 @@ #include <rsys/rsys.h> +enum htrdr_bsdf_type { + HTRDR_BSDF_DIFFUSE, + HTRDR_BSDF_SPECULAR +}; + /* Forward declarations */ struct htrdr; struct htrdr_ground; struct s3d_hit; +struct ssf_bsdf; extern LOCAL_SYM res_T htrdr_ground_create (struct htrdr* htrdr, - const char* obj_filename, + const char* obj_filename, /* May be NULL <=> No ground geometry */ + const enum htrdr_bsdf_type bsdf_type, const double reflectivity, /* In [0, 1] */ const int repeat_ground, /* Infinitely repeat the ground in X and Y */ struct htrdr_ground** ground); @@ -39,8 +46,8 @@ extern LOCAL_SYM void htrdr_ground_ref_put (struct htrdr_ground* ground); -extern LOCAL_SYM double -htrdr_ground_get_reflectivity +extern LOCAL_SYM struct ssf_bsdf* +htrdr_ground_get_bsdf (const struct htrdr_ground* ground); extern LOCAL_SYM res_T diff --git a/src/htrdr_sky.c b/src/htrdr_sky.c @@ -209,7 +209,7 @@ struct htrdr_sky { ******************************************************************************/ static INLINE int aabb_intersect - (const double aabb0_low[3], + (const double aabb0_low[3], const double aabb0_upp[3], const double aabb1_low[3], const double aabb1_upp[3]) @@ -1667,10 +1667,6 @@ setup_sw_bands_properties(struct htrdr_sky* sky) size_t i; ASSERT(sky); - /* Fetch short wave bands range */ - res = htgop_get_sw_spectral_intervals_CIE_XYZ(sky->htgop, sky->sw_bands_range); - if(res != RES_OK) goto error; - nbands = htrdr_sky_get_sw_spectral_bands_count(sky); ASSERT(nbands); sky->sw_bands = MEM_CALLOC @@ -1774,7 +1770,7 @@ htrdr_sky_create int htgop_upd = 1; int force_upd = 1; res_T res = RES_OK; - ASSERT(htrdr && sun && htgop_filename && htmie_filename && out_sky); + ASSERT(htrdr && sun && htgop_filename && out_sky); ASSERT(optical_thickness_threshold >= 0); sky = MEM_CALLOC(htrdr->allocator, 1, sizeof(*sky)); @@ -1790,6 +1786,8 @@ htrdr_sky_create sky->repeat_clouds = repeat_clouds; sky->is_cloudy = htcp_filename != NULL; darray_split_init(htrdr->allocator, &sky->svx2htcp_z); + sky->sw_bands_range[0] = 1; + sky->sw_bands_range[1] = 0; /* Load the gas optical properties */ res = htgop_create @@ -1806,6 +1804,21 @@ htrdr_sky_create goto error; } + /* Fetch short wave bands range */ + res = htgop_get_sw_spectral_intervals_CIE_XYZ(sky->htgop, sky->sw_bands_range); + if(res != RES_OK) goto error; + + /* Setup the atmopshere */ + time_current(&t0); + res = setup_atmosphere(sky, optical_thickness_threshold); + if(res != RES_OK) goto error; + time_sub(&t0, time_current(&t1), &t0); + time_dump(&t0, TIME_ALL, NULL, buf, sizeof(buf)); + htrdr_log(htrdr, "Setup atmosphere in %s\n", buf); + + /* Nothing more to do */ + if(!sky->is_cloudy) goto exit; + /* Load MIE data */ res = htmie_create(&htrdr->logger, htrdr->allocator, htrdr->verbose, &sky->htmie); if(res != RES_OK) { @@ -1822,17 +1835,6 @@ htrdr_sky_create res = setup_sw_bands_properties(sky); if(res != RES_OK) goto error; - /* Setup the atmopshere */ - time_current(&t0); - res = setup_atmosphere(sky, optical_thickness_threshold); - if(res != RES_OK) goto error; - time_sub(&t0, time_current(&t1), &t0); - time_dump(&t0, TIME_ALL, NULL, buf, sizeof(buf)); - htrdr_log(htrdr, "Setup atmosphere in %s\n", buf); - - /* Nothing more to do */ - if(!sky->is_cloudy) goto exit; - /* Load clouds properties */ res = htcp_create(&htrdr->logger, htrdr->allocator, htrdr->verbose, &sky->htcp); if(res != RES_OK) { @@ -1915,8 +1917,12 @@ htrdr_sky_fetch_particle_phase_function_asymmetry_parameter ASSERT(ispectral_band >= sky->sw_bands_range[0]); ASSERT(ispectral_band <= sky->sw_bands_range[1]); (void)iquad; - i = ispectral_band - sky->sw_bands_range[0]; - return sky->sw_bands[i].g_avg; + if(!sky->is_cloudy) { + return 0; + } else { + i = ispectral_band - sky->sw_bands_range[0]; + return sky->sw_bands[i].g_avg; + } } double @@ -2276,7 +2282,7 @@ htrdr_sky_fetch_svx_voxel_property size_t htrdr_sky_get_sw_spectral_bands_count(const struct htrdr_sky* sky) { - ASSERT(sky); + ASSERT(sky && sky->sw_bands_range[0] <= sky->sw_bands_range[1]); return sky->sw_bands_range[1] - sky->sw_bands_range[0] + 1; }