rnsf

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

test_rnsf_wlens.c (10892B)


      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 #include "rnsf.h"
     24 
     25 #include <rsys/math.h>
     26 #include <rsys/mem_allocator.h>
     27 
     28 static void
     29 test_load(struct rnsf* rnsf)
     30 {
     31   struct rnsf_phase_fn_discrete discrete = RNSF_PHASE_FN_DISCRETE_NULL;
     32   struct rnsf_phase_fn_hg hg = RNSF_PHASE_FN_HG_NULL;
     33   const struct rnsf_phase_fn* phase = NULL;
     34   FILE* fp = NULL;
     35 
     36   CHK(fp = tmpfile());
     37   fprintf(fp, "wavelengths 0\n");
     38   rewind(fp);
     39   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
     40   CHK(rnsf_get_phase_fn_count(rnsf) == 0);
     41   CHK(fclose(fp) == 0);
     42 
     43   CHK(fp = tmpfile());
     44   fprintf(fp, "# Comment\n");
     45   fprintf(fp, "wavelengths 1\n");
     46   fprintf(fp, "\n");
     47   fprintf(fp, "200.1 HG 0\n");
     48   rewind(fp);
     49 
     50   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
     51 
     52   CHK(rnsf_get_phase_fn_count(rnsf) == 1);
     53   CHK(phase = rnsf_get_phase_fn(rnsf, 0));
     54   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
     55 
     56   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
     57   CHK(hg.wavelengths[0] == 200.1);
     58   CHK(hg.wavelengths[1] == 200.1);
     59   CHK(hg.g == 0);
     60   CHK(fclose(fp) == 0);
     61 
     62   CHK(fp = tmpfile());
     63   fprintf(fp, "wavelengths 10\n");
     64   fprintf(fp, "\n");
     65   fprintf(fp, "# Short waves\n");
     66   fprintf(fp, "430 discrete 8\n");
     67   fprintf(fp, "  0       0.02\n");
     68   fprintf(fp, "  0.23    0.04\n");
     69   fprintf(fp, "  0.5     0.07\n");
     70   fprintf(fp, "  0.7     0.15\n");
     71   fprintf(fp, "  1.54    1.23\n");
     72   fprintf(fp, "  1.8     0.02\n");
     73   fprintf(fp, "  2       1.23\n");
     74   fprintf(fp, "  3.14159 0.79\n");
     75   fprintf(fp, "450 discrete 2\n");
     76   fprintf(fp, "  0        0.5\n");
     77   fprintf(fp, "  3.14159 0.796\n");
     78   fprintf(fp, "750 discrete 5\n");
     79   fprintf(fp, "  0  %.9g\n", 1.0/(4.0*PI));
     80   fprintf(fp, "  %.9g %.9g\n", PI/4.0, 1.0/(4.0*PI));
     81   fprintf(fp, "  %.9g %.9g\n", PI/2.0, 1.0/(4.0*PI));
     82   fprintf(fp, "  %.9g %.9g\n", 3*PI/4.0, 1.0/(4.0*PI));
     83   fprintf(fp, "  %g %.9g\n", PI, 1.0/(4.0*PI));
     84   fprintf(fp, "\n");
     85   fprintf(fp, "# Long waves\n");
     86   fprintf(fp, "1100    HG   -0.1\n");
     87   fprintf(fp, "1300    HG    0.57\n");
     88   fprintf(fp, "1400    HG    0.4\n");
     89   fprintf(fp, "2100    HG    0.3\n");
     90   fprintf(fp, "2500    HG   -0.9\n");
     91   fprintf(fp, "2900    HG   -0.4\n");
     92   fprintf(fp, "100000  HG    0.0\n");
     93   rewind(fp);
     94 
     95   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
     96   CHK(rnsf_get_phase_fn_count(rnsf) == 10);
     97 
     98   CHK(phase = rnsf_get_phase_fn(rnsf, 0));
     99   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_DISCRETE);
    100   CHK(rnsf_phase_fn_get_discrete(phase, &discrete) == RES_OK);
    101   CHK(discrete.wavelengths[0] == 430);
    102   CHK(discrete.wavelengths[1] == 430);
    103   CHK(discrete.nitems == 8);
    104   CHK(discrete.items[0].theta == 0);
    105   CHK(discrete.items[1].theta == 0.23);
    106   CHK(discrete.items[2].theta == 0.5);
    107   CHK(discrete.items[3].theta == 0.7);
    108   CHK(discrete.items[4].theta == 1.54);
    109   CHK(discrete.items[5].theta == 1.8);
    110   CHK(discrete.items[6].theta == 2);
    111   CHK(eq_eps(discrete.items[7].theta, PI, 1.e-6));
    112 
    113   CHK(phase = rnsf_get_phase_fn(rnsf, 1));
    114   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_DISCRETE);
    115   CHK(rnsf_phase_fn_get_discrete(phase, &discrete) == RES_OK);
    116   CHK(discrete.wavelengths[0] == 450);
    117   CHK(discrete.wavelengths[1] == 450);
    118   CHK(discrete.nitems == 2);
    119   CHK(discrete.items[0].theta == 0);
    120   CHK(eq_eps(discrete.items[1].theta, PI, 1.e-6));
    121 
    122   CHK(phase = rnsf_get_phase_fn(rnsf, 2));
    123   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_DISCRETE);
    124   CHK(rnsf_phase_fn_get_discrete(phase, &discrete) == RES_OK);
    125   CHK(discrete.wavelengths[0] == 750);
    126   CHK(discrete.wavelengths[1] == 750);
    127   CHK(discrete.nitems == 5);
    128   CHK(discrete.items[0].theta == 0);
    129   CHK(eq_eps(discrete.items[1].theta, PI/4.0, 1.e-6));
    130   CHK(eq_eps(discrete.items[2].theta, PI/2.0, 1.e-6));
    131   CHK(eq_eps(discrete.items[3].theta, 3*PI/4.0, 1.e-6));
    132   CHK(eq_eps(discrete.items[4].theta, PI, 1.e-6));
    133   CHK(eq_eps(discrete.items[0].value, 1/(4.0*PI), 1.e-6));
    134   CHK(eq_eps(discrete.items[1].value, 1/(4.0*PI), 1.e-6));
    135   CHK(eq_eps(discrete.items[2].value, 1/(4.0*PI), 1.e-6));
    136   CHK(eq_eps(discrete.items[3].value, 1/(4.0*PI), 1.e-6));
    137   CHK(eq_eps(discrete.items[4].value, 1/(4.0*PI), 1.e-6));
    138 
    139   CHK(phase = rnsf_get_phase_fn(rnsf, 3));
    140   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    141   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    142   CHK(hg.wavelengths[0] == 1100);
    143   CHK(hg.wavelengths[1] == 1100);
    144   CHK(hg.g == -0.1);
    145 
    146   CHK(phase = rnsf_get_phase_fn(rnsf, 4));
    147   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    148   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    149   CHK(hg.wavelengths[0] == 1300);
    150   CHK(hg.wavelengths[1] == 1300);
    151   CHK(hg.g == 0.57);
    152 
    153   CHK(phase = rnsf_get_phase_fn(rnsf, 5));
    154   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    155   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    156   CHK(hg.wavelengths[0] == 1400);
    157   CHK(hg.wavelengths[1] == 1400);
    158   CHK(hg.g == 0.4);
    159 
    160   CHK(phase = rnsf_get_phase_fn(rnsf, 6));
    161   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    162   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    163   CHK(hg.wavelengths[0] == 2100);
    164   CHK(hg.wavelengths[1] == 2100);
    165   CHK(hg.g == 0.3);
    166 
    167   CHK(phase = rnsf_get_phase_fn(rnsf, 7));
    168   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    169   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    170   CHK(hg.wavelengths[0] == 2500);
    171   CHK(hg.wavelengths[1] == 2500);
    172   CHK(hg.g == -0.9);
    173 
    174   CHK(phase = rnsf_get_phase_fn(rnsf, 8));
    175   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    176   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    177   CHK(hg.wavelengths[0] == 2900);
    178   CHK(hg.wavelengths[1] == 2900);
    179   CHK(hg.g == -0.4);
    180 
    181   CHK(phase = rnsf_get_phase_fn(rnsf, 9));
    182   CHK(rnsf_phase_fn_get_type(phase) == RNSF_PHASE_FN_HG);
    183   CHK(rnsf_phase_fn_get_hg(phase, &hg) == RES_OK);
    184   CHK(hg.wavelengths[0] == 100000);
    185   CHK(hg.wavelengths[1] == 100000);
    186   CHK(hg.g == 0);
    187 
    188   CHK(fclose(fp) == 0);
    189 }
    190 
    191 static void
    192 test_fetch(struct rnsf* rnsf)
    193 {
    194   FILE* fp = NULL;
    195   size_t iphase_fn = 0;
    196 
    197   CHK(fp = tmpfile());
    198   fprintf(fp, "wavelengths 0\n");
    199   rewind(fp);
    200   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
    201   CHK(fclose(fp) == 0);
    202   CHK(rnsf_fetch_phase_fn(rnsf, 400, 0.5, &iphase_fn) == RES_OK);
    203   CHK(iphase_fn >= rnsf_get_phase_fn_count(rnsf));
    204 
    205   CHK(fp = tmpfile());
    206   fprintf(fp, "wavelengths 2\n");
    207   fprintf(fp, "200 HG 0\n");
    208   fprintf(fp, "300 HG 0\n");
    209   rewind(fp);
    210   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
    211   CHK(fclose(fp) == 0);
    212 
    213   CHK(rnsf_fetch_phase_fn(rnsf, 299, 0, &iphase_fn) == RES_OK);
    214   CHK(iphase_fn == 0);
    215   CHK(rnsf_fetch_phase_fn(rnsf, 301, 0, &iphase_fn) == RES_OK);
    216   CHK(iphase_fn == 1);
    217   CHK(rnsf_fetch_phase_fn(rnsf, 200, 0, &iphase_fn) == RES_OK);
    218   CHK(iphase_fn == 0);
    219   CHK(rnsf_fetch_phase_fn(rnsf, 300, 0, &iphase_fn) == RES_OK);
    220   CHK(iphase_fn == 1);
    221   CHK(rnsf_fetch_phase_fn(rnsf, 280, 0, &iphase_fn) == RES_OK);
    222   CHK(iphase_fn == 0);
    223   CHK(rnsf_fetch_phase_fn(rnsf, 280, 0.19, &iphase_fn) == RES_OK);
    224   CHK(iphase_fn == 0);
    225   CHK(rnsf_fetch_phase_fn(rnsf, 280, 0.20, &iphase_fn) == RES_OK);
    226   CHK(iphase_fn == 1);
    227 
    228   CHK(fp = tmpfile());
    229   fprintf(fp, "wavelengths 4\n");
    230   fprintf(fp, "100 HG 0\n");
    231   fprintf(fp, "200 HG 1\n");
    232   fprintf(fp, "300 HG -1\n");
    233   fprintf(fp, "400 discrete 4\n");
    234   fprintf(fp, " 0 1\n");
    235   fprintf(fp, " 0.5 1\n");
    236   fprintf(fp, " 1.57 1\n");
    237   fprintf(fp, " 3.1416 1\n");
    238   rewind(fp);
    239   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
    240   CHK(fclose(fp) == 0);
    241 
    242   CHK(rnsf_fetch_phase_fn(rnsf, 50, 0, &iphase_fn) == RES_OK);
    243   CHK(iphase_fn == 0);
    244   CHK(rnsf_fetch_phase_fn(rnsf, 120, 0.8, &iphase_fn) == RES_OK);
    245   CHK(iphase_fn == 1);
    246   CHK(rnsf_fetch_phase_fn(rnsf, 270, 0.29, &iphase_fn) == RES_OK);
    247   CHK(iphase_fn == 1);
    248   CHK(rnsf_fetch_phase_fn(rnsf, 270, 0.3, &iphase_fn) == RES_OK);
    249   CHK(iphase_fn == 2);
    250   CHK(rnsf_fetch_phase_fn(rnsf, 350, 0.5, &iphase_fn) == RES_OK);
    251   CHK(iphase_fn == 3);
    252   CHK(rnsf_fetch_phase_fn(rnsf, 450, 0, &iphase_fn) == RES_OK);
    253   CHK(iphase_fn == 3);
    254 }
    255 
    256 static void
    257 test_load_fail(struct rnsf* rnsf)
    258 {
    259   FILE* fp = NULL;
    260 
    261   /* No wavelength count */
    262   CHK(fp = tmpfile());
    263   fprintf(fp, "wavelengths\n");
    264   fprintf(fp, "380 HG 0.5\n");
    265   rewind(fp);
    266   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_BAD_ARG);
    267   CHK(fclose(fp) == 0);
    268 
    269   /* Missing a phase function */
    270   CHK(fp = tmpfile());
    271   fprintf(fp, "wavelengths 2\n");
    272   fprintf(fp, "380 HG 0.5\n");
    273   rewind(fp);
    274   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_BAD_ARG);
    275   CHK(fclose(fp) == 0);
    276 
    277   /* Invalid wavelength */
    278   CHK(fp = tmpfile());
    279   fprintf(fp, "wavelengths 1\n");
    280   fprintf(fp, "-380 HG 0.5\n");
    281   rewind(fp);
    282   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_BAD_ARG);
    283   CHK(fclose(fp) == 0);
    284 
    285   /* Unsorted phase functions */
    286   CHK(fp = tmpfile());
    287   fprintf(fp, "wavelengths 2\n");
    288   fprintf(fp, "380 HG 0.5\n");
    289   fprintf(fp, "280 HG 0.0\n");
    290   rewind(fp);
    291   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_BAD_ARG);
    292   CHK(fclose(fp) == 0);
    293 
    294   /* Phase functions overlap */
    295   CHK(fp = tmpfile());
    296   fprintf(fp, "wavelengths 2\n");
    297   fprintf(fp, "280 HG 0.5\n");
    298   fprintf(fp, "280 HG 0.0\n");
    299   rewind(fp);
    300   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_BAD_ARG);
    301   CHK(fclose(fp) == 0);
    302 
    303   /* Additional text. Print a warning */
    304   CHK(fp = tmpfile());
    305   fprintf(fp, "wavelengths 2 additional_text\n");
    306   fprintf(fp, "280 HG 0.5\n");
    307   fprintf(fp, "380 HG 0.0\n");
    308   rewind(fp);
    309   CHK(rnsf_load_stream(rnsf, fp, NULL) == RES_OK);
    310   CHK(fclose(fp) == 0);
    311 }
    312 
    313 int
    314 main(int argc, char** argv)
    315 {
    316   struct rnsf_create_args args = RNSF_CREATE_ARGS_DEFAULT;
    317   struct rnsf* rnsf = NULL;
    318   (void)argc, (void)argv;
    319 
    320   args.verbose = 1;
    321   CHK(rnsf_create(&args, &rnsf) == RES_OK);
    322   CHK(rnsf_get_phase_fn_count(rnsf) == 0);
    323 
    324   test_load(rnsf);
    325   test_fetch(rnsf);
    326   test_load_fail(rnsf);
    327 
    328   CHK(rnsf_ref_put(rnsf) == RES_OK);
    329   CHK(mem_allocated_size() == 0);
    330   return 0;
    331 }