rnsf

commit cef3a9f18d00951cf86c132b8ab03535923d996f
parent 94a163739882a8d96de39a2949be44bcd8337935
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Tue, 21 Jun 2022 16:10:30 +0200

Normalize the discretized phase function

Diffstat:
M
cmake/CMakeLists.txt
 | 
2
+-
M
src/rnsf.c
 | 
42
++++++++++++++++++++++++++++++++++++++----

2 files changed, 39 insertions(+), 5 deletions(-)
diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt
@@ -67,7 +67,7 @@ rcmake_prepend_path(RNSF_FILES_INC_API ${RNSF_SOURCE_DIR})
 rcmake_prepend_path(RNSF_FILES_DOC ${PROJECT_SOURCE_DIR}/../)
 
 add_library(rnsf SHARED ${RNSF_FILES_SRC} ${RNSF_FILES_INC} ${RNSF_FILES_INC_API})
-target_link_libraries(rnsf RSys)
+target_link_libraries(rnsf RSys m)
 
 set_target_properties(rnsf PROPERTIES
   DEFINE_SYMBOL RNSF_SHARED_BUILD
diff --git a/src/rnsf.c b/src/rnsf.c
@@ -29,12 +29,13 @@
 #include "rnsf_c.h"
 #include "rnsf_log.h"
 
-#include <string.h>
-
 #include <rsys/cstr.h>
 #include <rsys/math.h>
 #include <rsys/text_reader.h>
 
+#include <math.h>
+#include <string.h>
+
 /*******************************************************************************
  * Helper functions
  ******************************************************************************/
@@ -45,6 +46,36 @@ check_rnsf_create_args(const struct rnsf_create_args* args)
   return args ? RES_OK : RES_BAD_ARG;
 }
 
+static void
+normalize_phase_fn_discret(struct phase_fn_discrete* discrete)
+{
+  struct rnsf_discrete_item* items = NULL;
+  size_t i, n;
+  double phi_prev;
+  double alpha;
+  ASSERT(discrete);
+
+  n = darray_discrete_item_size_get(&discrete->items);
+  ASSERT(n >= 2);
+
+  items = darray_discrete_item_data_get(&discrete->items);
+  ASSERT(items[0].theta == 0 && items[n-1].theta == PI);
+
+  alpha = 0;
+  phi_prev = 1;
+  FOR_EACH(i, 1, n) {
+    const double phi_curr = cos(items[i].theta);
+    ASSERT(phi_prev > phi_curr);
+    alpha += (items[i-1].value + items[i].value) * (phi_prev - phi_curr);
+    phi_prev = phi_curr;
+  }
+  alpha *= PI;
+
+  if(!eq_eps(alpha, 1, 1.e-6)) {
+    FOR_EACH(i, 0, n) items[i].value /= alpha;
+  }
+}
+
 static res_T
 parse_phase_fn_hg
   (struct rnsf* rnsf,
@@ -106,7 +137,7 @@ parse_discrete_item
   }
 
   if(!txtrdr_get_cline(txtrdr)) {
-    log_err(rnsf, "%s:%lu: missing an item of a discretized phase function.\n",
+    log_err(rnsf, "%s:%lu: missing an item of the discretized phase function.\n",
       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
     res = RES_BAD_ARG;
     goto error;
@@ -221,7 +252,8 @@ parse_phase_fn_discrete
 
     if(i == nangles - 1 && item->theta != PI) {
       log_err(rnsf,
-        "%s:%lu: invalid angle value `%g'. The last angle must be PI.\n",
+        "%s:%lu: invalid angle value `%g'. The last angle must be the "
+        "litteral PI.\n",
         txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
         item->theta);
       res = RES_BAD_ARG;
@@ -274,6 +306,8 @@ parse_phase_fn
     res = parse_phase_fn_discrete(rnsf, txtrdr, tk_ctx, &phase->param.discrete);
     if(res != RES_OK) goto error;
 
+    normalize_phase_fn_discret(&phase->param.discrete);
+
   } else {
     log_err(rnsf, "%s:%lu: invalid phase function type `%s'.\n",
       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);