rnsf

Define and load a phase function data format
git clone git://git.meso-star.fr/rnsf.git
Log | Files | Refs | README | LICENSE

rnsf.c (23016B)


      1 /* Copyright (C) 2022, 2023 Centre National de la Recherche Scientifique
      2  * Copyright (C) 2022, 2023 Institut Pierre-Simon Laplace
      3  * Copyright (C) 2022, 2023 Institut de Physique du Globe de Paris
      4  * Copyright (C) 2022, 2023 |Méso|Star> (contact@meso-star.com)
      5  * Copyright (C) 2022, 2023 Observatoire de Paris
      6  * Copyright (C) 2022, 2023 Université de Reims Champagne-Ardenne
      7  * Copyright (C) 2022, 2023 Université de Versaille Saint-Quentin
      8  * Copyright (C) 2022, 2023 Université Paul Sabatier
      9  *
     10  * This program is free software: you can redistribute it and/or modify
     11  * it under the terms of the GNU General Public License as published by
     12  * the Free Software Foundation, either version 3 of the License, or
     13  * (at your option) any later version.
     14  *
     15  * This program is distributed in the hope that it will be useful,
     16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     18  * GNU General Public License for more details.
     19  *
     20  * You should have received a copy of the GNU General Public License
     21  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     22 
     23 #define _POSIX_C_SOURCE 200112L /* strtok_r */
     24 
     25 #define PI_EPSILON 1.0e-4
     26 
     27 #include "rnsf.h"
     28 #include "rnsf_c.h"
     29 #include "rnsf_log.h"
     30 
     31 #include <rsys/cstr.h>
     32 #include <rsys/math.h>
     33 #include <rsys/text_reader.h>
     34 
     35 #include <math.h>
     36 #include <string.h>
     37 
     38 /*******************************************************************************
     39  * Helper functions
     40  ******************************************************************************/
     41 static INLINE res_T
     42 check_rnsf_create_args(const struct rnsf_create_args* args)
     43 {
     44   /* Nothing to check. Only return RES_BAD_ARG if args is NULL */
     45   return args ? RES_OK : RES_BAD_ARG;
     46 }
     47 
     48 static void
     49 normalize_phase_fn_discret(struct phase_fn_discrete* discrete)
     50 {
     51   struct rnsf_discrete_item* items = NULL;
     52   size_t i, n;
     53   double phi_prev;
     54   double alpha;
     55   ASSERT(discrete);
     56 
     57   n = darray_discrete_item_size_get(&discrete->items);
     58   ASSERT(n >= 2);
     59 
     60   items = darray_discrete_item_data_get(&discrete->items);
     61   ASSERT(items[0].theta == 0 && items[n-1].theta == PI);
     62 
     63   alpha = 0;
     64   phi_prev = 1;
     65   FOR_EACH(i, 1, n) {
     66     const double phi_curr = cos(items[i].theta);
     67     ASSERT(phi_prev > phi_curr);
     68     alpha += (items[i-1].value + items[i].value) * (phi_prev - phi_curr);
     69     phi_prev = phi_curr;
     70   }
     71   alpha *= PI;
     72 
     73   if(!eq_eps(alpha, 1, 1.e-6)) {
     74     FOR_EACH(i, 0, n) items[i].value /= alpha;
     75   }
     76 }
     77 
     78 static res_T
     79 parse_phase_fn_hg
     80   (struct rnsf* rnsf,
     81    struct txtrdr* txtrdr,
     82    char** tk_ctx,
     83    struct phase_fn_hg* phase)
     84 {
     85   char* tk = NULL;
     86   res_T res = RES_OK;
     87   ASSERT(rnsf && txtrdr && tk_ctx && phase);
     88 
     89   tk = strtok_r(NULL, " \t", tk_ctx);
     90   if(!tk) {
     91     log_err(rnsf,
     92       "%s:%lu: missing the Henyey & Greenstein asymmetric parameter.\n",
     93       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
     94     res = RES_BAD_ARG;
     95     goto error;
     96   }
     97 
     98   res = cstr_to_double(tk, &phase->g);
     99   if(res == RES_OK && (phase->g < -1 || phase->g > 1)) res = RES_BAD_ARG;
    100   if(res != RES_OK) {
    101     log_err(rnsf,
    102       "%s:%lu: invalid Henyey & Greenstein asymmetric parameter `%s'.\n",
    103       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    104     goto error;
    105   }
    106 
    107   tk = strtok_r(NULL, " \t", tk_ctx);
    108   if(tk) {
    109     log_warn(rnsf, "%s:%lu: unexpected text `%s'.\n",
    110       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    111   }
    112 
    113 exit:
    114   return res;
    115 error:
    116   goto exit;
    117 }
    118 
    119 static res_T
    120 parse_discrete_item
    121   (struct rnsf* rnsf,
    122    struct txtrdr* txtrdr,
    123    struct rnsf_discrete_item* item)
    124 {
    125   char* tk = NULL;
    126   char* tk_ctx = NULL;
    127   res_T res = RES_OK;
    128   ASSERT(rnsf && txtrdr && item);
    129 
    130   res = txtrdr_read_line(txtrdr);
    131   if(res != RES_OK) {
    132     log_err(rnsf, "%s: could not read the line `%lu' -- %s.\n",
    133       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    134       res_to_cstr(res));
    135     goto error;
    136   }
    137 
    138   if(!txtrdr_get_cline(txtrdr)) {
    139     log_err(rnsf, "%s:%lu: missing an item of the discretized phase function.\n",
    140       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    141     res = RES_BAD_ARG;
    142     goto error;
    143   }
    144 
    145   tk = strtok_r(txtrdr_get_line(txtrdr), " \t", &tk_ctx);
    146   ASSERT(tk);
    147 
    148   res = cstr_to_double(tk, &item->theta);
    149   if(res == RES_OK 
    150   && (  item->theta < 0 
    151      || (!eq_eps(item->theta, PI, PI_EPSILON) && item->theta > PI))) {
    152     res = RES_BAD_ARG;
    153   }
    154   if(res != RES_OK) {
    155     log_err(rnsf, "%s:%lu: invalid phase function angle `%s'.\n",
    156       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    157     goto error;
    158   }
    159 
    160   tk = strtok_r(NULL, " \t", &tk_ctx);
    161   if(!tk) {
    162     log_err(rnsf,
    163       "%s:%lu: the value of the discretized phase function is missing.\n",
    164       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    165     res = RES_BAD_ARG;
    166     goto error;
    167   }
    168 
    169   res = cstr_to_double(tk, &item->value);
    170   if(res == RES_OK && item->value < 0) res = RES_BAD_ARG;
    171   if(res != RES_OK) {
    172     log_err(rnsf, "%s:%lu: invalid value `%s' of the discretized phase function.\n",
    173       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    174     goto error;
    175   }
    176 
    177   tk = strtok_r(NULL, " \t", &tk_ctx);
    178   if(tk) {
    179     log_warn(rnsf, "%s:%lu: unexpected text `%s'.\n",
    180       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    181   }
    182 
    183 exit:
    184   return res;
    185 error:
    186   goto exit;
    187 }
    188 
    189 static res_T
    190 parse_phase_fn_discrete
    191   (struct rnsf* rnsf,
    192    struct txtrdr* txtrdr,
    193    char** tk_ctx,
    194    struct phase_fn_discrete* phase)
    195 {
    196   char* tk = NULL;
    197   unsigned long nangles = 0;
    198   unsigned long i;
    199   res_T res = RES_OK;
    200   ASSERT(rnsf && txtrdr && tk_ctx && phase);
    201 
    202   tk = strtok_r(NULL, " \t", tk_ctx);
    203   if(!tk) {
    204     log_err(rnsf,
    205       "%s:%lu: missing the angle number of the discretized phase function.\n",
    206       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    207     res = RES_BAD_ARG;
    208     goto error;
    209   }
    210 
    211   res = cstr_to_ulong(tk, &nangles);
    212   if(res == RES_OK && nangles < 2) res = RES_BAD_ARG;
    213   if(res != RES_OK) {
    214     log_err(rnsf,
    215       "%s:%lu: invalid angle number of the discretized phase function `%s'.\n",
    216       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    217     goto error;
    218   }
    219 
    220   res = darray_discrete_item_resize(&phase->items, nangles);
    221   if(res != RES_OK) {
    222     log_err(rnsf,
    223       "%s:%lu: could not allocate the list of values of the discretized "
    224       "phase function -- %s.\n",
    225       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    226       res_to_cstr(res));
    227     goto error;
    228   }
    229 
    230   tk = strtok_r(NULL, " \t", tk_ctx);
    231   if(tk) {
    232     log_warn(rnsf, "%s:%lu: unexpected text `%s'.\n",
    233       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    234   }
    235 
    236   FOR_EACH(i, 0, nangles) {
    237     struct rnsf_discrete_item* item = NULL;
    238     item = darray_discrete_item_data_get(&phase->items) + i;
    239 
    240     res = parse_discrete_item(rnsf, txtrdr, item);
    241     if(res != RES_OK) goto error;
    242 
    243     if(i == 0 && item->theta != 0) {
    244       log_err(rnsf,
    245         "%s:%lu: invalid angle value `%g'. The first angle must be 0.\n",
    246         txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    247         item->theta);
    248       res = RES_BAD_ARG;
    249       goto error;
    250     }
    251 
    252     if(i == nangles - 1) {
    253       if(eq_eps(item->theta, PI, PI_EPSILON)) {
    254         item->theta = PI;
    255       } else {
    256         log_err(rnsf,
    257           "%s:%lu: invalid angle value `%g'. The last angle must be 3.14159.\n",
    258           txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    259           item->theta);
    260         res = RES_BAD_ARG;
    261         goto error;
    262       }
    263     }
    264 
    265     if(i > 0 && item[0].theta <= item[-1].theta) {
    266       log_err(rnsf,
    267         "%s:%lu: the discretized phase function must be sorted in ascending "
    268         "order of angles (item %lu at %g rad; item %lu at %g rad).\n",
    269         txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    270         (unsigned long)(i-1), item[-1].theta,
    271         (unsigned long)(i),   item[ 0].theta);
    272       res = RES_BAD_ARG;
    273       goto error;
    274     }
    275   }
    276 
    277 exit:
    278   return res;
    279 error:
    280   goto exit;
    281 }
    282 
    283 static res_T
    284 parse_phase_fn
    285   (struct rnsf* rnsf,
    286    struct txtrdr* txtrdr,
    287    char** tk_ctx,
    288    struct rnsf_phase_fn* phase)
    289 {
    290   char* tk = NULL;
    291   res_T res = RES_OK;
    292   ASSERT(rnsf && txtrdr && tk_ctx && phase);
    293 
    294   tk = strtok_r(NULL, " \t", tk_ctx);
    295   if(!tk) {
    296     log_err(rnsf, "%s:%lu: missing phase function type.\n",
    297       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    298     res = RES_BAD_ARG;
    299     goto error;
    300   }
    301 
    302   if(!strcmp(tk, "HG")) {
    303     phase->type = RNSF_PHASE_FN_HG;
    304     res = parse_phase_fn_hg(rnsf, txtrdr, tk_ctx, &phase->param.hg);
    305     if(res != RES_OK) goto error;
    306 
    307   } else if(!strcmp(tk, "discrete")) {
    308     phase->type = RNSF_PHASE_FN_DISCRETE;
    309     phase_fn_discrete_init(rnsf->allocator, &phase->param.discrete);
    310     res = parse_phase_fn_discrete(rnsf, txtrdr, tk_ctx, &phase->param.discrete);
    311     if(res != RES_OK) goto error;
    312 
    313     normalize_phase_fn_discret(&phase->param.discrete);
    314 
    315   } else {
    316     log_err(rnsf, "%s:%lu: invalid phase function type `%s'.\n",
    317       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    318     res = RES_BAD_ARG;
    319     goto error;
    320   }
    321 
    322 exit:
    323   return res;
    324 error:
    325   goto exit;
    326 }
    327 
    328 
    329 static res_T
    330 parse_wlen_phase_fn
    331   (struct rnsf* rnsf,
    332    struct txtrdr* txtrdr,
    333    struct rnsf_phase_fn* phase)
    334 {
    335   char* tk = NULL;
    336   char* tk_ctx = NULL;
    337   double wlen = 0;
    338   res_T res = RES_OK;
    339   ASSERT(rnsf && txtrdr && phase);
    340 
    341   res = txtrdr_read_line(txtrdr);
    342   if(res != RES_OK) {
    343     log_err(rnsf, "%s: could not read the line `%lu' -- %s.\n",
    344       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    345       res_to_cstr(res));
    346     goto error;
    347   }
    348 
    349   if(!txtrdr_get_cline(txtrdr)) {
    350     const size_t nexpect = darray_phase_fn_size_get(&rnsf->phases);
    351     const size_t nparsed = (size_t)
    352       (phase - darray_phase_fn_cdata_get(&rnsf->phases));
    353     log_err(rnsf,
    354       "%s:%lu: missing phase function. "
    355       "Expecting %lu function%swhile %lu %s parsed.\n",
    356       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    357       (unsigned long)nexpect, nexpect == 1 ? " " : "s ",
    358       (unsigned long)nparsed, nparsed > 1 ? "were" : "was");
    359     res = RES_BAD_ARG;
    360     goto error;
    361   }
    362 
    363   tk = strtok_r(txtrdr_get_line(txtrdr), " \t", &tk_ctx);
    364   ASSERT(tk);
    365 
    366   res = cstr_to_double(tk, &wlen);
    367   if(res == RES_OK && wlen < 0) res = RES_BAD_ARG;
    368   if(res != RES_OK) {
    369     log_err(rnsf, "%s:%lu: invalid wavelength `%s'.\n",
    370       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    371     goto error;
    372   }
    373   phase->wlen[0] = wlen;
    374   phase->wlen[1] = wlen;
    375 
    376   res = parse_phase_fn(rnsf, txtrdr, &tk_ctx, phase);
    377   if(res != RES_OK) goto error;
    378 
    379 exit:
    380   return res;
    381 error:
    382   goto exit;
    383 }
    384 
    385 static res_T
    386 parse_per_wlen_phase_fn
    387   (struct rnsf* rnsf,
    388    struct txtrdr* txtrdr,
    389    char** tk_ctx)
    390 {
    391   char* tk = NULL;
    392   unsigned long count = 0;
    393   unsigned long i;
    394   res_T res = RES_OK;
    395   ASSERT(rnsf && txtrdr && tk_ctx);
    396 
    397   tk = strtok_r(NULL, " \t", tk_ctx);
    398   if(!tk) {
    399     log_err(rnsf, "%s:%lu: the number of wavelengths is missing.\n",
    400       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    401     res = RES_BAD_ARG;
    402     goto error;
    403   }
    404 
    405   res = cstr_to_ulong(tk, &count);
    406   if(res != RES_OK) {
    407     log_err(rnsf, "%s:%lu: invalid number of wavelengths `%s'.\n",
    408       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    409     res = RES_BAD_ARG;
    410     goto error;
    411   }
    412 
    413   tk = strtok_r(NULL, " \t", tk_ctx);
    414   if(tk) {
    415     log_warn(rnsf, "%s:%lu: unexpected text `%s'.\n",
    416       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    417   }
    418 
    419   res = darray_phase_fn_resize(&rnsf->phases, count);
    420   if(res != RES_OK) {
    421     log_err(rnsf,
    422       "%s:%lu: could not allocate the list of phase functions -- %s.\n",
    423       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    424       res_to_cstr(res));
    425     goto error;
    426   }
    427 
    428   FOR_EACH(i, 0, count) {
    429     struct rnsf_phase_fn* phase = NULL;
    430 
    431     phase = darray_phase_fn_data_get(&rnsf->phases) + i;
    432 
    433     res = parse_wlen_phase_fn(rnsf, txtrdr, phase);
    434     if(res != RES_OK) goto error;
    435     ASSERT(phase->wlen[0] == phase->wlen[1]);
    436 
    437     if(i > 0 && phase[0].wlen[0] <= phase[-1].wlen[0]) {
    438       log_err(rnsf,
    439         "%s:%lu: the phase functions must be sorted in ascending order of "
    440         "wavelengths (phase function %lu at %g nm; phase function %lu at %g nm).\n",
    441         txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    442         (unsigned long)(i-1), phase[-1].wlen[0],
    443         (unsigned long)(i),   phase[ 0].wlen[0]);
    444       res = RES_BAD_ARG;
    445       goto error;
    446     }
    447   }
    448 
    449 exit:
    450   return res;
    451 error:
    452   darray_phase_fn_clear(&rnsf->phases);
    453   goto exit;
    454 }
    455 
    456 static res_T
    457 parse_band_phase
    458   (struct rnsf* rnsf,
    459    struct txtrdr* txtrdr,
    460    struct rnsf_phase_fn* phase)
    461 {
    462   char* tk = NULL;
    463   char* tk_ctx = NULL;
    464   res_T res = RES_OK;
    465   ASSERT(rnsf && txtrdr && phase);
    466 
    467   res = txtrdr_read_line(txtrdr);
    468   if(res != RES_OK) {
    469     log_err(rnsf, "%s: could not read the line `%lu' -- %s.\n",
    470       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    471       res_to_cstr(res));
    472     goto error;
    473   }
    474 
    475   if(!txtrdr_get_cline(txtrdr)) {
    476     const size_t nexpect = darray_phase_fn_size_get(&rnsf->phases);
    477     const size_t nparsed = (size_t)
    478       (phase - darray_phase_fn_cdata_get(&rnsf->phases));
    479     log_err(rnsf,
    480       "%s:%lu: missing phase function. "
    481       "Expecting %lu function%swhile %lu %s parsed.\n",
    482       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    483       (unsigned long)nexpect, nexpect == 1 ? " " : "s ",
    484       (unsigned long)nparsed, nparsed > 1 ? "were" : "was");
    485     res = RES_BAD_ARG;
    486     goto error;
    487   }
    488 
    489   tk = strtok_r(txtrdr_get_line(txtrdr), " \t", &tk_ctx);
    490   ASSERT(tk);
    491 
    492   res = cstr_to_double(tk, &phase->wlen[0]);
    493   if(res == RES_OK && phase->wlen[0] < 0) res = RES_BAD_ARG;
    494   if(res != RES_OK) {
    495     log_err(rnsf, "%s:%lu: invalid band min boundary `%s'.\n",
    496       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    497     goto error;
    498   }
    499 
    500   tk = strtok_r(NULL, " \t", &tk_ctx);
    501   if(!tk) {
    502     log_err(rnsf, "%s:%lu: missing band max boundary.\n",
    503       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    504     res = RES_BAD_ARG;
    505     goto error;
    506   }
    507 
    508   res = cstr_to_double(tk, &phase->wlen[1]);
    509   if(res == RES_OK && phase->wlen[1] < 0) res = RES_BAD_ARG;
    510   if(res != RES_OK) {
    511     log_err(rnsf, "%s:%lu: invalid band max boundary `%s'.\n",
    512       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    513     goto error;
    514   }
    515 
    516   if(phase->wlen[0] > phase->wlen[1]) {
    517     log_err(rnsf, "%s:%lu: band boundaries are degenerated -- [%g, %g].\n",
    518       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    519       phase->wlen[0], phase->wlen[1]);
    520     res = RES_BAD_ARG;
    521     goto error;
    522   }
    523 
    524   res = parse_phase_fn(rnsf, txtrdr, &tk_ctx, phase);
    525   if(res != RES_OK) goto error;
    526 
    527 exit:
    528   return res;
    529 error:
    530   goto exit;
    531 }
    532 
    533 static res_T
    534 parse_per_band_phase_fn
    535   (struct rnsf* rnsf,
    536    struct txtrdr* txtrdr,
    537    char** tk_ctx)
    538 {
    539   char* tk = NULL;
    540   unsigned long count = 0;
    541   unsigned long i = 0;
    542   res_T res = RES_OK;
    543   ASSERT(rnsf && txtrdr && tk_ctx);
    544 
    545   tk = strtok_r(NULL, " \t", tk_ctx);
    546   if(!tk) {
    547     log_err(rnsf, "%s:%lu: the number of bands is missing.\n",
    548       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr));
    549     res = RES_BAD_ARG;
    550     goto error;
    551   }
    552 
    553   res = cstr_to_ulong(tk, &count);
    554   if(res != RES_OK) {
    555     log_err(rnsf, "%s:%lu: invalid number of bands `%s'.\n",
    556       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    557     res = RES_BAD_ARG;
    558     goto error;
    559   }
    560 
    561   tk = strtok_r(NULL, " \t", tk_ctx);
    562   if(tk) {
    563     log_warn(rnsf, "%s:%lu: unexpected text `%s'.\n",
    564       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    565   }
    566 
    567   res = darray_phase_fn_resize(&rnsf->phases, count);
    568   if(res != RES_OK) {
    569     log_err(rnsf,
    570       "%s:%lu: could not allocate the list of phase functions -- %s.\n",
    571       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    572       res_to_cstr(res));
    573     goto error;
    574   }
    575 
    576   FOR_EACH(i, 0, count) {
    577     struct rnsf_phase_fn* phase = NULL;
    578 
    579     phase = darray_phase_fn_data_get(&rnsf->phases) + i;
    580 
    581     res = parse_band_phase(rnsf, txtrdr, phase);
    582     if(res != RES_OK) goto error;
    583 
    584     if(i > 0) {
    585       if(phase[0].wlen[0] <  phase[-1].wlen[1]
    586       || phase[0].wlen[0] == phase[-1].wlen[0]) {
    587         log_err(rnsf,
    588           "%s:%lu: the phase functions must be sorted in ascending order of "
    589           "wavelengths and must not overlap (phase function %lu in [%g, %g] nm; "
    590           "phase function %lu in [%g %g] nm).\n",
    591           txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    592           (unsigned long)(i-1), phase[-1].wlen[0], phase[-1].wlen[1],
    593           (unsigned long)(i),   phase[ 0].wlen[0], phase[ 0].wlen[1]);
    594         res = RES_BAD_ARG;
    595         goto error;
    596       }
    597     }
    598   }
    599 
    600 exit:
    601   return res;
    602 error:
    603   darray_phase_fn_clear(&rnsf->phases);
    604   goto exit;
    605 }
    606 
    607 static res_T
    608 load_stream
    609   (struct rnsf* rnsf,
    610    FILE* stream,
    611    const char* stream_name)
    612 {
    613   char* tk = NULL;
    614   char* tk_ctx = NULL;
    615   struct txtrdr* txtrdr = NULL;
    616   res_T res = RES_OK;
    617   ASSERT(rnsf && stream && stream_name);
    618 
    619   darray_phase_fn_clear(&rnsf->phases); /* Clean-up */
    620 
    621   res = txtrdr_stream(rnsf->allocator, stream, stream_name, '#', &txtrdr);
    622   if(res != RES_OK) {
    623     log_err(rnsf, "Could not create the text reader to parse `%s' -- %s.\n",
    624       stream_name, res_to_cstr(res));
    625     goto error;
    626   }
    627 
    628   res = txtrdr_read_line(txtrdr);
    629   if(res != RES_OK) {
    630     log_err(rnsf, "%s: could not read the line `%lu' -- %s.\n",
    631       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr),
    632       res_to_cstr(res));
    633     goto error;
    634   }
    635 
    636   if(!txtrdr_get_cline(txtrdr)) goto exit; /* No line to parse */
    637 
    638   tk = strtok_r(txtrdr_get_line(txtrdr), " \t", &tk_ctx);
    639   ASSERT(tk);
    640 
    641   if(!strcmp(tk, "wavelengths")) {
    642     res = parse_per_wlen_phase_fn(rnsf, txtrdr, &tk_ctx);
    643     if(res != RES_BAD_ARG) goto error;
    644   } else if(!strcmp(tk, "bands")) {
    645     res = parse_per_band_phase_fn(rnsf, txtrdr, &tk_ctx);
    646     if(res != RES_BAD_ARG) goto error;
    647   } else {
    648     log_err(rnsf, "%s:%lu: invalid directive `%s'.\n",
    649       txtrdr_get_name(txtrdr), (unsigned long)txtrdr_get_line_num(txtrdr), tk);
    650     res = RES_BAD_ARG;
    651     goto error;
    652   }
    653 
    654 exit:
    655   if(txtrdr) txtrdr_ref_put(txtrdr);
    656   return res;
    657 error:
    658   goto exit;
    659 }
    660 
    661 static void
    662 release_rnsf(ref_T* ref)
    663 {
    664   struct rnsf* rnsf = NULL;
    665   ASSERT(ref);
    666   rnsf = CONTAINER_OF(ref, struct rnsf, ref);
    667   if(rnsf->logger == &rnsf->logger__) logger_release(&rnsf->logger__);
    668   darray_phase_fn_release(&rnsf->phases);
    669   MEM_RM(rnsf->allocator, rnsf);
    670 }
    671 
    672 /*******************************************************************************
    673  * Exported functions
    674  ******************************************************************************/
    675 res_T
    676 rnsf_create
    677   (const struct rnsf_create_args* args,
    678    struct rnsf** out_rnsf)
    679 {
    680   struct rnsf* rnsf = NULL;
    681   struct mem_allocator* allocator = NULL;
    682   res_T res = RES_OK;
    683 
    684   if(!out_rnsf) { res = RES_BAD_ARG; goto error; }
    685   res = check_rnsf_create_args(args);
    686   if(res != RES_OK) goto error;
    687 
    688   allocator = args->allocator ? args->allocator : &mem_default_allocator;
    689   rnsf = MEM_CALLOC(allocator, 1, sizeof(*rnsf));
    690   if(!rnsf) {
    691     if(args->verbose) {
    692       #define ERR_STR "Could not allocate the device of the "\
    693         "Rad-Net Scattering Function library.\n"
    694       if(args->logger) {
    695         logger_print(args->logger, LOG_ERROR, ERR_STR);
    696       } else {
    697         fprintf(stderr, MSG_ERROR_PREFIX ERR_STR);
    698       }
    699       #undef ERR_STR
    700     }
    701     res = RES_MEM_ERR;
    702     goto error;
    703   }
    704 
    705   ref_init(&rnsf->ref);
    706   rnsf->allocator = allocator;
    707   rnsf->verbose = args->verbose;
    708   darray_phase_fn_init(allocator, &rnsf->phases);
    709   if(args->logger) {
    710     rnsf->logger = args->logger;
    711   } else {
    712     res = setup_log_default(rnsf);
    713     if(res != RES_OK) {
    714       if(args->verbose) {
    715         fprintf(stderr, MSG_ERROR_PREFIX
    716           "Could not setup the default logger "
    717           "of the Rad-Net Scattering Function library.\n");
    718       }
    719       goto error;
    720     }
    721   }
    722 
    723 exit:
    724   if(out_rnsf) *out_rnsf = rnsf ;
    725   return res;
    726 error:
    727   if(rnsf) { RNSF(ref_put(rnsf)); rnsf = NULL; }
    728   goto exit;
    729 }
    730 
    731 res_T
    732 rnsf_ref_get(struct rnsf* rnsf)
    733 {
    734   if(!rnsf) return RES_BAD_ARG;
    735   ref_get(&rnsf->ref);
    736   return RES_OK;
    737 }
    738 
    739 res_T
    740 rnsf_ref_put(struct rnsf* rnsf)
    741 {
    742   if(!rnsf) return RES_BAD_ARG;
    743   ref_put(&rnsf->ref, release_rnsf);
    744   return RES_OK;
    745 }
    746 
    747 res_T
    748 rnsf_load(struct rnsf* rnsf, const char* filename)
    749 {
    750   FILE* fp = NULL;
    751   res_T res = RES_OK;
    752 
    753   if(!rnsf || !filename) {
    754     res = RES_BAD_ARG;
    755     goto error;
    756   }
    757 
    758   fp = fopen(filename, "r");
    759   if(!fp) {
    760     log_err(rnsf, "%s: error opening file `%s' -- %s.\n",
    761       FUNC_NAME, filename, strerror(errno));
    762     res = RES_IO_ERR;
    763     goto error;
    764   }
    765 
    766   res = load_stream(rnsf, fp, filename);
    767   if(res != RES_OK) goto error;
    768 
    769 exit:
    770   if(fp) fclose(fp);
    771   return res;
    772 error:
    773   goto exit;
    774 }
    775 
    776 res_T
    777 rnsf_load_stream
    778   (struct rnsf* rnsf,
    779    FILE* stream,
    780    const char* stream_name)
    781 {
    782   res_T res = RES_OK;
    783 
    784   if(!rnsf || !stream) {
    785     res = RES_BAD_ARG;
    786     goto error;
    787   }
    788 
    789   res = load_stream(rnsf, stream, stream_name ? stream_name : "<stream>");
    790   if(res != RES_OK) goto error;
    791 
    792 exit:
    793   return res;
    794 error:
    795   goto exit;
    796 }
    797 
    798 size_t
    799 rnsf_get_phase_fn_count(const struct rnsf* rnsf)
    800 {
    801   ASSERT(rnsf);
    802   return darray_phase_fn_size_get(&rnsf->phases);
    803 }
    804 
    805 const struct rnsf_phase_fn*
    806 rnsf_get_phase_fn(const struct rnsf* rnsf, const size_t iphase_fn)
    807 {
    808   ASSERT(rnsf && iphase_fn < rnsf_get_phase_fn_count(rnsf));
    809   return darray_phase_fn_cdata_get(&rnsf->phases) + iphase_fn;
    810 }
    811 
    812 enum rnsf_phase_fn_type
    813 rnsf_phase_fn_get_type(const struct rnsf_phase_fn* phase)
    814 {
    815   ASSERT(phase);
    816   return phase->type;
    817 }
    818 
    819 res_T
    820 rnsf_phase_fn_get_hg
    821   (const struct rnsf_phase_fn* phase,
    822    struct rnsf_phase_fn_hg* hg)
    823 {
    824   if(!phase || !hg) return RES_BAD_ARG;
    825   if(phase->type != RNSF_PHASE_FN_HG) return RES_BAD_ARG;
    826   hg->wavelengths[0] = phase->wlen[0];
    827   hg->wavelengths[1] = phase->wlen[1];
    828   hg->g = phase->param.hg.g;
    829   return RES_OK;
    830 }
    831 
    832 res_T
    833 rnsf_phase_fn_get_discrete
    834   (const struct rnsf_phase_fn* phase,
    835    struct rnsf_phase_fn_discrete* discrete)
    836 {
    837   if(!phase || !discrete) return RES_BAD_ARG;
    838   if(phase->type != RNSF_PHASE_FN_DISCRETE) return RES_BAD_ARG;
    839   discrete->wavelengths[0] = phase->wlen[0];
    840   discrete->wavelengths[1] = phase->wlen[1];
    841   discrete->items = darray_discrete_item_cdata_get(&phase->param.discrete.items);
    842   discrete->nitems = darray_discrete_item_size_get(&phase->param.discrete.items);
    843   return RES_OK;
    844 }