solstice-solver

Solver library of the solstice app
git clone git://git.meso-star.com/solstice-solver.git
Log | Files | Refs | README | LICENSE

ssol_sun.c (8517B)


      1 /* Copyright (C) 2018-2026 |Méso|Star> (contact@meso-star.com)
      2  * Copyright (C) 2016, 2018 CNRS
      3  *
      4  * This program is free software: you can redistribute it and/or modify
      5  * it under the terms of the GNU General Public License as published by
      6  * the Free Software Foundation, either version 3 of the License, or
      7  * (at your option) any later version.
      8  *
      9  * This program is distributed in the hope that it will be useful,
     10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     12  * GNU General Public License for more details.
     13  *
     14  * You should have received a copy of the GNU General Public License
     15  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     16 
     17 #include "ssol.h"
     18 #include "ssol_device_c.h"
     19 #include "ssol_sun_c.h"
     20 #include "ssol_ranst_sun_dir.h"
     21 #include "ssol_ranst_sun_wl.h"
     22 #include "ssol_spectrum_c.h"
     23 
     24 #include <star/scem.h>
     25 
     26 #include <rsys/rsys.h>
     27 #include <rsys/mem_allocator.h>
     28 #include <rsys/ref_count.h>
     29 #include <rsys/math.h>
     30 #include <rsys/double3.h>
     31 
     32 #include <string.h>
     33 
     34 /*******************************************************************************
     35  * Helper functions
     36  ******************************************************************************/
     37 static void
     38 sun_release(ref_T* ref)
     39 {
     40   struct ssol_device* dev;
     41   struct ssol_sun* sun = CONTAINER_OF(ref, struct ssol_sun, ref);
     42   ASSERT(ref);
     43   dev = sun->dev;
     44   ASSERT(dev && dev->allocator);
     45   if(sun->spectrum) SSOL(spectrum_ref_put(sun->spectrum));
     46   MEM_RM(dev->allocator, sun);
     47   SSOL(device_ref_put(dev));
     48 }
     49 
     50 static res_T
     51 sun_create
     52   (struct ssol_device* dev, struct ssol_sun** out_sun, enum sun_type type)
     53 {
     54   struct ssol_sun* sun = NULL;
     55   res_T res = RES_OK;
     56   if(!dev || !out_sun || type >= SUN_TYPES_COUNT__) {
     57     return RES_BAD_ARG;
     58   }
     59 
     60   sun = (struct ssol_sun*)MEM_CALLOC
     61     (dev->allocator, 1, sizeof(struct ssol_sun));
     62   if(!sun) {
     63     res = RES_MEM_ERR;
     64     goto error;
     65   }
     66 
     67   SSOL(device_ref_get(dev));
     68   sun->dev = dev;
     69   sun->type = type;
     70   sun->algorithm = SSOL_SUN_DIRECTION_ALGORITHM_MEEUS;
     71   ref_init(&sun->ref);
     72 
     73 exit:
     74   if(out_sun) *out_sun = sun;
     75   return res;
     76 error:
     77   if(sun) {
     78     SSOL(sun_ref_put(sun));
     79     sun = NULL;
     80   }
     81   goto exit;
     82 }
     83 
     84 /*******************************************************************************
     85  * Exported ssol_image functions
     86  ******************************************************************************/
     87 res_T
     88 ssol_sun_create_directional(struct ssol_device* dev, struct ssol_sun** out_sun)
     89 {
     90   return sun_create(dev, out_sun, SUN_DIRECTIONAL);
     91 }
     92 
     93 res_T
     94 ssol_sun_create_pillbox(struct ssol_device* dev, struct ssol_sun** out_sun)
     95 {
     96   return sun_create(dev, out_sun, SUN_PILLBOX);
     97 }
     98 
     99 res_T
    100 ssol_sun_create_gaussian(struct ssol_device* dev, struct ssol_sun** out_sun)
    101 {
    102   return sun_create(dev, out_sun, SUN_GAUSSIAN);
    103 }
    104 
    105 res_T
    106 ssol_sun_create_buie
    107   (struct ssol_device* dev, struct ssol_sun** out_sun)
    108 {
    109   return sun_create(dev, out_sun, SUN_BUIE);
    110 }
    111 
    112 res_T
    113 ssol_sun_ref_get(struct ssol_sun* sun)
    114 {
    115   if(!sun)
    116     return RES_BAD_ARG;
    117   ref_get(&sun->ref);
    118   return RES_OK;
    119 }
    120 
    121 res_T
    122 ssol_sun_ref_put(struct ssol_sun* sun)
    123 {
    124   if(!sun)
    125     return RES_BAD_ARG;
    126   ref_put(&sun->ref, sun_release);
    127   return RES_OK;
    128 }
    129 
    130 res_T
    131 ssol_sun_set_direction(struct ssol_sun* sun, const double direction[3])
    132 {
    133   if(!sun || !direction)
    134     return RES_BAD_ARG;
    135   if(0 == d3_normalize(sun->direction, direction))
    136     /* zero vector */
    137     return RES_BAD_ARG;
    138   return RES_OK;
    139 }
    140 
    141 res_T
    142 ssol_sun_set_location_and_date
    143   (struct ssol_sun* sun,
    144    const struct ssol_location* loc,
    145    const struct tm* date)
    146 {
    147   res_T res = RES_OK;
    148   struct scem_sun_pos sun_pos;
    149   struct scem_location scem_loc;
    150   enum scem_sun_algo algorithm;
    151 
    152   if(!sun || !loc|| !date) {
    153     res = RES_BAD_ARG;
    154     goto error;
    155   }
    156 
    157   if(loc->latitude < -90 || loc->latitude > 90
    158       || loc->longitude < -180 || loc->longitude > 180)
    159   {
    160     res = RES_BAD_ARG;
    161     goto error;
    162   }
    163 
    164   scem_loc.latitude = loc->latitude;
    165   scem_loc.longitude = loc->longitude;
    166 
    167   switch(sun->algorithm ) {
    168     case SSOL_SUN_DIRECTION_ALGORITHM_PSA :
    169       algorithm = SCEM_SUN_PSA ;
    170       break;
    171     case SSOL_SUN_DIRECTION_ALGORITHM_MEEUS :
    172       algorithm = SCEM_SUN_MEEUS ;
    173       break;
    174     default: FATAL("Unknown sun positioning algorithm.");
    175   }
    176   res = scem_sun_position_from_earth(date, &scem_loc, algorithm, &sun_pos);
    177   if(res != RES_OK) goto error;
    178   res = scem_sun_position_to_sun_vector(&sun_pos, sun->direction);
    179   if(res != RES_OK) goto error;
    180 
    181 exit:
    182   return res;
    183 error:
    184   goto exit;
    185 }
    186 
    187 res_T
    188 ssol_sun_set_algorithm
    189   (struct ssol_sun* sun,
    190    const enum ssol_sun_direction_algorithm algorithm)
    191 {
    192   if(!sun ||
    193       (algorithm != SSOL_SUN_DIRECTION_ALGORITHM_MEEUS
    194        && algorithm != SSOL_SUN_DIRECTION_ALGORITHM_PSA))
    195     return RES_BAD_ARG;
    196   sun->algorithm = algorithm;
    197   return RES_OK;
    198 }
    199 
    200 res_T
    201 ssol_sun_get_algorithm
    202   (const struct ssol_sun* sun,
    203    enum ssol_sun_direction_algorithm* algorithm)
    204 {
    205   if(!sun || !algorithm)
    206     return RES_BAD_ARG;
    207   *algorithm = sun->algorithm;
    208   return RES_OK;
    209 }
    210 
    211 res_T
    212 ssol_sun_get_direction(const struct ssol_sun* sun, double direction[3])
    213 {
    214   if(!sun || !direction)
    215     return RES_BAD_ARG;
    216   d3_set(direction, sun->direction);
    217   return RES_OK;
    218 }
    219 
    220 res_T
    221 ssol_sun_set_dni(struct ssol_sun* sun, const double dni)
    222 {
    223   if(!sun || dni <= 0)
    224     return RES_BAD_ARG;
    225   sun->dni = dni;
    226   return RES_OK;
    227 }
    228 
    229 res_T
    230 ssol_sun_get_dni(const struct ssol_sun* sun, double* dni)
    231 {
    232   if(!sun || !dni)
    233     return RES_BAD_ARG;
    234   *dni = sun->dni;
    235   return RES_OK;
    236 }
    237 
    238 res_T
    239 ssol_sun_set_spectrum(struct ssol_sun* sun, struct ssol_spectrum* spectrum)
    240 {
    241   if(!sun || !spectrum)
    242     return RES_BAD_ARG;
    243   if(spectrum == sun->spectrum) /* no change */
    244     return RES_OK;
    245   if(sun->spectrum)
    246     SSOL(spectrum_ref_put(sun->spectrum));
    247   SSOL(spectrum_ref_get(spectrum));
    248   sun->spectrum = spectrum;
    249   return RES_OK;
    250 }
    251 
    252 res_T
    253 ssol_sun_pillbox_set_half_angle(struct ssol_sun* sun, const double half_angle)
    254 {
    255   if(!sun || half_angle <= 0 || half_angle > PI * 0.5 || sun->type != SUN_PILLBOX)
    256     return RES_BAD_ARG;
    257   sun->data.pillbox.half_angle = half_angle;
    258   return RES_OK;
    259 }
    260 
    261 res_T
    262 ssol_sun_gaussian_set_std_dev(struct ssol_sun* sun, const double std_dev)
    263 {
    264   if(!sun || std_dev <= 0 || sun->type != SUN_GAUSSIAN)
    265     return RES_BAD_ARG;
    266   sun->data.gaussian.std_dev = std_dev;
    267   return RES_OK;
    268 }
    269 
    270 res_T
    271 ssol_sun_set_buie_param
    272   (struct ssol_sun* sun,
    273    const double ratio)
    274 {
    275   if(!sun
    276   || ratio <= 0
    277   || ratio >= 1
    278   || sun->type != SUN_BUIE)
    279     return RES_BAD_ARG;
    280   sun->data.csr.ratio = ratio;
    281   return RES_OK;
    282 }
    283 
    284 /*******************************************************************************
    285  * Local function
    286  ******************************************************************************/
    287 res_T
    288 sun_create_direction_distribution
    289   (struct ssol_sun* sun, struct ranst_sun_dir** out_ran_dir)
    290 {
    291   struct ranst_sun_dir* ran_dir = NULL;
    292   res_T res = RES_OK;
    293   ASSERT(sun && out_ran_dir);
    294 
    295   res = ranst_sun_dir_create(sun->dev->allocator, &ran_dir);
    296   if(res != RES_OK) goto error;
    297   switch(sun->type) {
    298     case SUN_DIRECTIONAL:
    299       res = ranst_sun_dir_dirac_setup(ran_dir, sun->direction);
    300       break;
    301     case SUN_PILLBOX:
    302       res = ranst_sun_dir_pillbox_setup
    303         (ran_dir, sun->data.pillbox.half_angle, sun->direction);
    304       break;
    305     case SUN_GAUSSIAN:
    306       res = ranst_sun_dir_gaussian_setup
    307         (ran_dir, sun->data.gaussian.std_dev, sun->direction);
    308       break;
    309     case SUN_BUIE:
    310       res = ranst_sun_dir_buie_setup
    311         (ran_dir, sun->data.csr.ratio, sun->direction);
    312       break;
    313     default: FATAL("Unreachable code\n"); break;
    314   }
    315 exit:
    316   *out_ran_dir = ran_dir;
    317   return res;
    318 error:
    319   if(ran_dir) {
    320     CHK(ranst_sun_dir_ref_put(ran_dir) == RES_OK);
    321     ran_dir = NULL;
    322   }
    323   goto exit;
    324 }
    325 
    326 res_T
    327 sun_create_wavelength_distribution
    328   (struct ssol_sun* sun, struct ranst_sun_wl** out_ran_wl)
    329 {
    330   struct ranst_sun_wl* ran_wl = NULL;
    331   res_T res = RES_OK;
    332   ASSERT(sun && out_ran_wl);
    333 
    334   res = ranst_sun_wl_create(sun->dev->allocator, &ran_wl);
    335   if(res != RES_OK) goto error;
    336 
    337   res = ranst_sun_wl_setup(ran_wl,
    338     darray_double_cdata_get(&sun->spectrum->wavelengths),
    339     darray_double_cdata_get(&sun->spectrum->intensities),
    340     darray_double_size_get(&sun->spectrum->wavelengths));
    341   if(res != RES_OK) goto error;
    342 
    343 exit:
    344   *out_ran_wl = ran_wl;
    345   return res;
    346 error:
    347   if(ran_wl) {
    348     CHK(ranst_sun_wl_ref_put(ran_wl) == RES_OK);
    349     ran_wl = NULL;
    350   }
    351   goto exit;
    352 }