commit f65b00ea6f4468e2ce4d02b1dd7aa9cce7c4b325
parent ec361da113dfe30dda367335e435cd7905113d6e
Author: Eduardo Fontana Lazzari <edufonlaz@gmail.com>
Date: Thu, 5 Mar 2026 17:27:36 +0100
Add tests for internally defined BTDF models
Implement one test for each internally defined BTDF: KEEP_CURRENT_DIR,
LAMBERT, and SNELL_DIELECTRIC.
Each test builds a three-surface configuration (emitter -> transmitter
-> absorber) where the photon density reaching the absorber can be
computed analytically. The Monte Carlo results, computed with the LOSSES
reported by the MVREA algorithm, are compared against the expected
solution within statistical tolerance.
Diffstat:
5 files changed, 882 insertions(+), 2 deletions(-)
diff --git a/.gitignore b/.gitignore
@@ -9,6 +9,7 @@ spec.txt
*.so
*.swp
*.stl
+*.dat
tags
input
output
diff --git a/Makefile b/Makefile
@@ -174,7 +174,10 @@ uninstall:
################################################################################
TEST_SRC =\
src/test_sphor_lib.c\
- src/test_sphor_MVREA_analytical1.c
+ src/test_sphor_MVREA_analytical1.c\
+ src/test_sphor_MVREA_btdf_keep_current_dir.c\
+ src/test_sphor_MVREA_btdf_snell_dielectric.c\
+ src/test_sphor_MVREA_btdf_lambertian.c
TEST_OBJ = $(TEST_SRC:.c=.o)
TEST_DEP = $(TEST_SRC:.c=.d)
TEST_TGT = $(TEST_SRC:.c=.t)
@@ -222,11 +225,14 @@ $(TEST_OBJ): config.mk sphor-local.pc
test_sphor_lib\
test_sphor_MVREA_analytical1\
+test_sphor_MVREA_btdf_keep_current_dir\
+test_sphor_MVREA_btdf_snell_dielectric\
+test_sphor_MVREA_btdf_lambertian\
: config.mk sphor-local.pc $(LIBNAME)
$(CC) $(CFLAGS_TEST) -o $@ src/$@.o $(LDFLAGS_TEST)
clean_test:
rm -f $(TEST_DEP) $(TEST_OBJ) $(TEST_TGT) \
- input output cube.stl cube_back.stl cube_front.stl \
+ input output *.stl *.dat \
spec.txt rho.txt
for i in $(TEST_SRC); do rm -f "$$(basename "$${i}" ".c")"; done
diff --git a/src/test_sphor_MVREA_btdf_keep_current_dir.c b/src/test_sphor_MVREA_btdf_keep_current_dir.c
@@ -0,0 +1,272 @@
+/* Copyright (C) 2024-2026 Centre National de la Recherche Scientifique
+ * Copyright (C) 2024-2026 Clermont Auvergne INP
+ * Copyright (C) 2024-2026 INSA Lyon
+ * Copyright (C) 2024-2026 Institut Mines Télécom Albi-Carmaux
+ * Copyright (C) 2024-2026 Institut National Polytechnique de Toulouse
+ * Copyright (C) 2024-2026 |Méso|Star> (contact@meso-star.com)
+ * Copyright (C) 2024-2026 PhotonLyX (info@photonlyx.com)
+ * Copyright (C) 2024-2026 Université de Lorraine
+ * Copyright (C) 2024-2026 Université Paul Sabatier
+ * Copyright (C) 2024-2026 Université Toulouse - Jean Jaurès
+ *
+ * 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/>. */
+#define _POSIX_C_SOURCE 200112L
+
+#include "sphor.h"
+
+#include <rsys/cstr.h>
+#include <rsys/logger.h>
+#include <rsys/str.h>
+#include <rsys/text_reader.h>
+
+#include <math.h>
+#include <stdio.h>
+
+static void
+write_emitting_surface()
+{
+ FILE* fp = fopen("emitting_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid emitting_surface\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -0.5773502691896257 -2 -1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -1.4433756729740643 -1.5 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid emitting_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf_surface()
+{
+ FILE* fp = fopen("btdf_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid btdf_surface\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 1\n"
+ " vertex 0 -1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 -1\n"
+ " vertex 0 1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid btdf_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_absorbing_surface()
+{
+ FILE* fp = fopen("absorbing_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid absorbing_surface\n"
+ " facet normal -0.5 -0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex 1.4433756729740643 1.5 -1\n"
+ " vertex 0.5773502691896257 2 1\n"
+ " vertex 0.5773502691896257 2 -1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal -0.5 -0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex 0.5773502691896257 2 1\n"
+ " vertex 1.4433756729740643 1.5 -1\n"
+ " vertex 1.4433756729740643 1.5 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid absorbing_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_spectrum()
+{
+ FILE* fp = fopen("emission_spectrum.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf()
+{
+ FILE* fp = fopen("btdf.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_bsdf_null()
+{
+ FILE* fp = fopen("bsdf_null.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 0\n");
+ fprintf(fp, "401 0\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_input_file(FILE* fp)
+{
+ write_emitting_surface();
+ write_btdf_surface();
+ write_absorbing_surface();
+ write_spectrum();
+ write_btdf();
+ write_bsdf_null();
+
+ fprintf(fp, "surface: \"source\"\n");
+ fprintf(fp, "\tgeometry: FRONT emitting_surface.stl\n");
+ fprintf(fp, "\tsource:\n");
+ fprintf(fp, "\tflux_density: 1 umol/m^2/s emission_spectrum.dat nm nm^-1\n");
+ fprintf(fp, "\t\tdirection: COLLIM NORMAL\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"transmitter\"\n");
+ fprintf(fp, "\tgeometry: BACK btdf_surface.stl\n");
+ fprintf(fp, "\tbtdf: KEEP_CURRENT_DIR btdf.dat\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"absorber\"\n");
+ fprintf(fp, "\tgeometry: FRONT absorbing_surface.stl\n");
+ fprintf(fp, "\tbrdf: SPECULAR bsdf_null.dat\n");
+ fprintf(fp, "\tsensor:\n");
+ fprintf(fp, "\t\tresponse_function: 1\n");
+
+ CHK(fflush(fp) == 0);
+}
+
+int main(int argc, char** argv)
+{
+ struct sphor* sphor = NULL;
+ struct sphor_create_args args = SPHOR_CREATE_ARGS_DEFAULT;
+ struct txtrdr* txtrdr = NULL;
+ struct str line;
+ char* token = NULL;
+ char* token_ptr = NULL;
+ char input_filename[] = "input";
+ char output_filename[] = "output";
+ double avg = 0;
+ double std = 0;
+ double nsamples = 0;
+
+ FILE* fp = NULL;
+
+ fp = fopen(input_filename, "w+");
+ CHK(NULL != fp);
+
+ write_input_file(fp);
+ CHK(fclose(fp) == 0);
+
+ args.input_filename = input_filename;
+ args.output_filename = output_filename;
+ args.force = 1;
+
+ (void)argc;
+ (void)argv;
+
+ /* Test sphor_create with NULL logger and allocator */
+ CHK(sphor_create(&args, &sphor) == RES_OK);
+
+ /* Run sphor with the config file */
+ CHK(sphor_run(sphor) == RES_OK);
+
+ CHK(sphor_ref_put(sphor) == RES_OK);
+
+ fp = fopen(output_filename, "r");
+ CHK(NULL != fp);
+
+ str_init(NULL, &line);
+ CHK(txtrdr_stream(NULL, fp, output_filename, '#', &txtrdr) == RES_OK);
+
+ nsamples = (double)args.samples;
+
+ /* Verify first level: total scene MVREA */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "MVREA") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 0, 2 * std / sqrt(nsamples)));
+
+ /* Verify first level: total scene losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 1, 2 * std / sqrt(nsamples)));
+
+ /* Verify second level: per surface losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(strcmp(token, "absorber") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 1, 2 * std / sqrt(nsamples)));
+
+ CHK(fclose(fp) == 0);
+
+ str_release(&line);
+ txtrdr_ref_put(txtrdr);
+
+ CHK(mem_allocated_size() == 0);
+
+ return 0;
+}
diff --git a/src/test_sphor_MVREA_btdf_lambertian.c b/src/test_sphor_MVREA_btdf_lambertian.c
@@ -0,0 +1,292 @@
+/* Copyright (C) 2024-2026 Centre National de la Recherche Scientifique
+ * Copyright (C) 2024-2026 Clermont Auvergne INP
+ * Copyright (C) 2024-2026 INSA Lyon
+ * Copyright (C) 2024-2026 Institut Mines Télécom Albi-Carmaux
+ * Copyright (C) 2024-2026 Institut National Polytechnique de Toulouse
+ * Copyright (C) 2024-2026 |Méso|Star> (contact@meso-star.com)
+ * Copyright (C) 2024-2026 PhotonLyX (info@photonlyx.com)
+ * Copyright (C) 2024-2026 Université de Lorraine
+ * Copyright (C) 2024-2026 Université Paul Sabatier
+ * Copyright (C) 2024-2026 Université Toulouse - Jean Jaurès
+ *
+ * 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/>. */
+#define _POSIX_C_SOURCE 200112L
+
+#include "sphor.h"
+
+#include <rsys/cstr.h>
+#include <rsys/logger.h>
+#include <rsys/str.h>
+#include <rsys/text_reader.h>
+
+#include <math.h>
+#include <stdio.h>
+
+static void
+write_emitting_surface()
+{
+ FILE* fp = fopen("emitting_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid emitting_surface\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -0.5773502691896257 -2 -1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -1.4433756729740643 -1.5 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid emitting_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf_surface()
+{
+ FILE* fp = fopen("btdf_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid btdf_surface\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 1\n"
+ " vertex 0 -1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 -1\n"
+ " vertex 0 1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid btdf_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_absorbing_surface()
+{
+ FILE* fp = fopen("absorbing_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid absorbing_surface\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 1 -1 -1\n"
+ " vertex 1 1 1\n"
+ " vertex 1 -1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 1 -1 -1\n"
+ " vertex 1 1 -1\n"
+ " vertex 1 1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid absorbing_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_spectrum()
+{
+ FILE* fp = fopen("emission_spectrum.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf()
+{
+ FILE* fp = fopen("btdf.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_bsdf_null()
+{
+ FILE* fp = fopen("bsdf_null.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 0\n");
+ fprintf(fp, "401 0\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_input_file(FILE* fp)
+{
+ write_emitting_surface();
+ write_btdf_surface();
+ write_absorbing_surface();
+ write_spectrum();
+ write_btdf();
+ write_bsdf_null();
+
+ fprintf(fp, "surface: \"source\"\n");
+ fprintf(fp, "\tgeometry: FRONT emitting_surface.stl\n");
+ fprintf(fp, "\tsource:\n");
+ fprintf(fp, "\tflux_density: 2 umol/m^2/s emission_spectrum.dat nm nm^-1\n");
+ fprintf(fp, "\t\tdirection: COLLIM NORMAL\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"transmitter\"\n");
+ fprintf(fp, "\tgeometry: BACK btdf_surface.stl\n");
+ fprintf(fp, "\tbtdf: LAMBERT btdf.dat\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"absorber\"\n");
+ fprintf(fp, "\tgeometry: BACK absorbing_surface.stl\n");
+ fprintf(fp, "\tbrdf: SPECULAR bsdf_null.dat\n");
+ fprintf(fp, "\tsensor:\n");
+ fprintf(fp, "\t\tresponse_function: 1\n");
+
+ CHK(fflush(fp) == 0);
+}
+
+int main(int argc, char** argv)
+{
+ struct sphor* sphor = NULL;
+ struct sphor_create_args args = SPHOR_CREATE_ARGS_DEFAULT;
+ struct txtrdr* txtrdr = NULL;
+ struct str line;
+ char* token = NULL;
+ char* token_ptr = NULL;
+ char input_filename[] = "input";
+ char output_filename[] = "output";
+ double avg = 0;
+ double std = 0;
+ double nsamples = 0;
+
+ FILE* fp = NULL;
+
+ fp = fopen(input_filename, "w+");
+ CHK(NULL != fp);
+
+ write_input_file(fp);
+ CHK(fclose(fp) == 0);
+
+ args.input_filename = input_filename;
+ args.output_filename = output_filename;
+ args.force = 1;
+
+ (void)argc;
+ (void)argv;
+
+ /* Test sphor_create with NULL logger and allocator */
+ CHK(sphor_create(&args, &sphor) == RES_OK);
+
+ /* Run sphor with the config file */
+ CHK(sphor_run(sphor) == RES_OK);
+
+ CHK(sphor_ref_put(sphor) == RES_OK);
+
+ fp = fopen(output_filename, "r");
+ CHK(NULL != fp);
+
+ str_init(NULL, &line);
+ CHK(txtrdr_stream(NULL, fp, output_filename, '#', &txtrdr) == RES_OK);
+
+ nsamples = (double)args.samples;
+
+ /* Verify first level: total scene MVREA */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "MVREA") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 0, 2 * std / sqrt(nsamples)));
+
+ /* Verify first level: total scene losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare against the analytical view factor.
+ *
+ * This computes the exact configuration (view) factor F12 between two
+ * coaxial, parallel square surfaces of side length 2 separated by a
+ * distance l = 1.
+ *
+ * The expression evaluated is:
+ *
+ * F12 = (1 / A1) int_{A1} int_{A2} ( l^2 / (pi r^4) ) dA2 dA1
+ *
+ * where:
+ * - A1 and A2 are the two square surfaces,
+ * - l is the axial separation,
+ * - r = ||p2 - p1|| is the distance between differential elements,
+ * - the integrand results from the general view factor formula
+ * (cos theta_1 cos theta_2) / (pi r^2),
+ * using cos theta_1 = cos theta_2 = l / r for parallel planes.
+ *
+ * The result is the density of diffuse power leaving the btdf surface
+ * that reaches the absorbing surface, weighted by the surface area ratio of
+ * the emitting surface and the absorbing surface.
+ */
+ CHK(eq_eps(avg, 0.4152532835771472, 2 * std / sqrt(nsamples)));
+
+ /* Verify second level: per surface losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(strcmp(token, "absorber") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ CHK(eq_eps(avg, 0.4152532835771472, 2 * std / sqrt(nsamples)));
+
+ CHK(fclose(fp) == 0);
+
+ str_release(&line);
+ txtrdr_ref_put(txtrdr);
+
+ CHK(mem_allocated_size() == 0);
+
+ return 0;
+}
diff --git a/src/test_sphor_MVREA_btdf_snell_dielectric.c b/src/test_sphor_MVREA_btdf_snell_dielectric.c
@@ -0,0 +1,309 @@
+/* Copyright (C) 2024-2026 Centre National de la Recherche Scientifique
+ * Copyright (C) 2024-2026 Clermont Auvergne INP
+ * Copyright (C) 2024-2026 INSA Lyon
+ * Copyright (C) 2024-2026 Institut Mines Télécom Albi-Carmaux
+ * Copyright (C) 2024-2026 Institut National Polytechnique de Toulouse
+ * Copyright (C) 2024-2026 |Méso|Star> (contact@meso-star.com)
+ * Copyright (C) 2024-2026 PhotonLyX (info@photonlyx.com)
+ * Copyright (C) 2024-2026 Université de Lorraine
+ * Copyright (C) 2024-2026 Université Paul Sabatier
+ * Copyright (C) 2024-2026 Université Toulouse - Jean Jaurès
+ *
+ * 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/>. */
+#define _POSIX_C_SOURCE 200112L
+
+#include "sphor.h"
+
+#include <rsys/cstr.h>
+#include <rsys/logger.h>
+#include <rsys/str.h>
+#include <rsys/text_reader.h>
+
+#include <math.h>
+#include <stdio.h>
+
+static void
+write_emitting_surface()
+{
+ FILE* fp = fopen("emitting_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid emitting_surface\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -0.5773502691896257 -2 -1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 0.5 0.8660254037844386 0\n"
+ " outer loop\n"
+ " vertex -0.5773502691896257 -2 1\n"
+ " vertex -1.4433756729740643 -1.5 -1\n"
+ " vertex -1.4433756729740643 -1.5 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid emitting_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf_surface()
+{
+ FILE* fp = fopen("btdf_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid btdf_surface\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 1\n"
+ " vertex 0 -1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ " facet normal 1 0 0\n"
+ " outer loop\n"
+ " vertex 0 -1 -1\n"
+ " vertex 0 1 -1\n"
+ " vertex 0 1 1\n"
+ " endloop\n"
+ " endfacet\n"
+ "endsolid btdf_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_absorbing_surface()
+{
+ FILE* fp = fopen("absorbing_surface.stl", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp,
+ "solid absorbing_surface\n"
+ "facet normal 0 0 0\n"
+ "outer loop\n"
+ "vertex 1.7320508075688772 0 -1\n"
+ "vertex 1.7320508075688772 0 1\n"
+ "vertex 0.8660254037844386 1.5 1\n"
+ "endloop\n"
+ "endfacet\n"
+ "facet normal 0 0 0\n"
+ "outer loop\n"
+ "vertex 0.8660254037844386 1.5 1\n"
+ "vertex 0.8660254037844386 1.5 -1\n"
+ "vertex 1.7320508075688772 0 -1\n"
+ "endloop\n"
+ "endfacet\n"
+ "endsolid absorbing_surface\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_spectrum()
+{
+ FILE* fp = fopen("emission_spectrum.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_btdf()
+{
+ FILE* fp = fopen("btdf.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_bsdf_null()
+{
+ FILE* fp = fopen("bsdf_null.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 0\n");
+ fprintf(fp, "401 0\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_n_back()
+{
+ FILE* fp = fopen("n_back.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1\n");
+ fprintf(fp, "401 1\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_n_front()
+{
+ FILE* fp = fopen("n_front.dat", "w");
+ CHK(fp != NULL);
+
+ fprintf(fp, "400 1.7320508075688772\n");
+ fprintf(fp, "401 1.7320508075688772\n");
+
+ CHK(fclose(fp) == 0);
+}
+
+static void
+write_input_file(FILE* fp)
+{
+ write_emitting_surface();
+ write_btdf_surface();
+ write_absorbing_surface();
+ write_spectrum();
+ write_btdf();
+ write_bsdf_null();
+ write_n_front();
+ write_n_back();
+
+ fprintf(fp, "surface: \"source\"\n");
+ fprintf(fp, "\tgeometry: FRONT emitting_surface.stl\n");
+ fprintf(fp, "\tsource:\n");
+ fprintf(fp, "\tflux_density: 1.7320508075688772 umol/m^2/s");
+ fprintf(fp, " emission_spectrum.dat nm nm^-1\n");
+ fprintf(fp, "\t\tdirection: COLLIM NORMAL\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"transmitter\"\n");
+ fprintf(fp, "\tgeometry: BACK btdf_surface.stl\n");
+ fprintf(fp, "\tbtdf: SNELL_DIELECTRIC btdf.dat\n");
+ fprintf(fp, "\n");
+
+ fprintf(fp, "surface: \"absorber\"\n");
+ fprintf(fp, "\tgeometry: FRONT absorbing_surface.stl\n");
+ fprintf(fp, "\tbrdf: SPECULAR bsdf_null.dat\n");
+ fprintf(fp, "\tsensor:\n");
+ fprintf(fp, "\t\tresponse_function: 1\n");
+
+ fprintf(fp, "volume: \"n_front\"\n");
+ fprintf(fp, "\tgeometry: FRONT btdf_surface.stl\n");
+ fprintf(fp, "\trefractive_index:\n");
+ fprintf(fp, "\t\tn_real: n_front.dat\n");
+
+ fprintf(fp, "volume: \"n_back\"\n");
+ fprintf(fp, "\tgeometry: BACK btdf_surface.stl\n");
+ fprintf(fp, "\trefractive_index:\n");
+ fprintf(fp, "\t\tn_real: n_back.dat\n");
+
+ CHK(fflush(fp) == 0);
+}
+
+int main(int argc, char** argv)
+{
+ struct sphor* sphor = NULL;
+ struct sphor_create_args args = SPHOR_CREATE_ARGS_DEFAULT;
+ struct txtrdr* txtrdr = NULL;
+ struct str line;
+ char* token = NULL;
+ char* token_ptr = NULL;
+ char input_filename[] = "input";
+ char output_filename[] = "output";
+ double avg = 0;
+ double std = 0;
+ double nsamples = 0;
+
+ FILE* fp = NULL;
+
+ fp = fopen(input_filename, "w+");
+ CHK(NULL != fp);
+
+ write_input_file(fp);
+ CHK(fclose(fp) == 0);
+
+ args.input_filename = input_filename;
+ args.output_filename = output_filename;
+ args.force = 1;
+
+ (void)argc;
+ (void)argv;
+
+ /* Test sphor_create with NULL logger and allocator */
+ CHK(sphor_create(&args, &sphor) == RES_OK);
+
+ /* Run sphor with the config file */
+ CHK(sphor_run(sphor) == RES_OK);
+
+ CHK(sphor_ref_put(sphor) == RES_OK);
+
+ fp = fopen(output_filename, "r");
+ CHK(NULL != fp);
+
+ str_init(NULL, &line);
+ CHK(txtrdr_stream(NULL, fp, output_filename, '#', &txtrdr) == RES_OK);
+
+ nsamples = (double)args.samples;
+
+ /* Verify first level: total scene MVREA */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "MVREA") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 0, 2 * std / sqrt(nsamples)));
+
+ /* Verify first level: total scene losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 1, 2 * std / sqrt(nsamples)));
+
+ /* Verify second level: per surface losses */
+ CHK(txtrdr_read_line(txtrdr) == RES_OK);
+ CHK(str_set(&line, txtrdr_get_cline(txtrdr)) == RES_OK);
+ token = strtok_r(str_get(&line), ":", &token_ptr);
+ CHK(strcmp(token, "LOSSES") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(strcmp(token, "absorber") == 0);
+ token = strtok_r(NULL, ":", &token_ptr);
+ CHK(cstr_to_double(token, &avg) == RES_OK);
+ CHK(cstr_to_double(token_ptr, &std) == RES_OK);
+ /* Compare to analytical solution */
+ CHK(eq_eps(avg, 1, 2 * std / sqrt(nsamples)));
+
+ CHK(fclose(fp) == 0);
+
+ str_release(&line);
+ txtrdr_ref_put(txtrdr);
+
+ CHK(mem_allocated_size() == 0);
+
+ return 0;
+}
