htrdr

Solving radiative transfer in heterogeneous media
git clone git://git.meso-star.fr/htrdr.git
Log | Files | Refs | README | LICENSE

commit 9103614c15cbb99d9f00765d4d9bc9b8ba1a0ebf
parent b995d4fb967e539c021d49e2b7dce2f0139ad456
Author: Vincent Forest <vincent.forest@meso-star.com>
Date:   Fri,  6 Jul 2018 16:44:39 +0200

Implement the sun data structure

Diffstat:
Mcmake/CMakeLists.txt | 3+++
Asrc/htrdr_c.h | 34++++++++++++++++++++++++++++++++++
Asrc/htrdr_sun.c | 167+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Asrc/htrdr_sun.h | 61+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
4 files changed, 265 insertions(+), 0 deletions(-)

diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt @@ -55,13 +55,16 @@ set(HTRDR_FILES_SRC htrdr_main.c htrdr_rectangle.c htrdr_sky.c + htrdr_sun.c htrdr_transmission.c) set(HTRDR_FILES_INC htrdr.h + htrdr_c.h htrdr_args.h htrdr_buffer.h htrdr_rectangle.h htrdr_sky.h + htrdr_sun.h htrdr_solve.h) set(HTDRD_FILES_DOC COPYING README.md) diff --git a/src/htrdr_c.h b/src/htrdr_c.h @@ -0,0 +1,34 @@ +/* Copyright (C) 2018 Université Paul Sabatier, |Meso|Star> + * + * 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/>. */ + +#ifndef HTRDR_C_H +#define HTRDR_C_H + +/* In nanometer */ +static FINLINE double +wavenumber_to_wavelength(const double nu/*In cm^-1*/) +{ + return 1.e7 / nu; +} + +/* In cm^-1 */ +static FINLINE double +wavelength_to_wavenumber(const double lambda/*In nanometer*/) +{ + return 1.e7 / lambda; +} + +#endif /* HTRDR_C_H */ + diff --git a/src/htrdr_sun.c b/src/htrdr_sun.c @@ -0,0 +1,167 @@ +/* Copyright (C) 2018 Université Paul Sabatier, |Meso|Star> + * + * 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/>. */ + +#include "htrdr.h" +#include "htrdr_sun.h" + +#include <rsys/double3.h> +#include <rsys/dynamic_array_double.h> +#include <rsys/ref_count.h> +#include <rsys/math.h> + +struct htrdr_sun { + /* Short wave radiance in W.m^-2.sr^-1, for each spectral interval */ + struct darray_double radiance_sw; + + /* Short wave spectral interval boundaries, in cm^-1 */ + struct darray_double wavenumbers_sw; + + double half_angle; /* In radian */ + double solid_angle; /* In sr; solid_angle = 2*PI(1 - cos(half_angle)) */ + + double main_dir[3]; /* Direction *toward* the sun */ + + ref_T ref; + struct htrdr* htrdr; +}; + +/******************************************************************************* + * Helper functions + ******************************************************************************/ +static void +release_sun(ref_T* ref) +{ + struct htrdr_sun* sun; + ASSERT(ref); + sun = CONTAINER_OF(ref, struct htrdr_sun, ref); + darray_double_release(&sun->radiance_sw); + darray_double_release(&sun->wavenumbers_sw); + MEM_RM(sun->htrdr->allocator, sun); +} + +/******************************************************************************* + * Local functions + ******************************************************************************/ +res_T +htrdr_sun_create(struct htrdr* htrdr, struct htrdr_sun** out_sun) +{ + const double incoming_flux_sw[] = { /* In W.m^-2 */ + 12.793835026999544, 12.109561093845551, 20.365091338928245, + 23.729742422870157, 22.427697221814142, 55.626612361454150, + 102.93146523363953, 24.293596268358986, 345.73659325842243, + 218.18441435866691, 347.18437832794524, 129.49426803812202, + 50.146977730963876, 3.1197193425713365 + }; + const double wavenumbers_sw[] = { /* In cm^-1 */ + 820.000, 2600.00, 3250.00, 4000.00, 4650.00, + 5150.00, 6150.00, 7700.00, 8050.00, 12850.0, + 16000.0, 22650.0, 29000.0, 38000.0, 49999.0 + }; + + const size_t nspectral_intervals = sizeof(incoming_flux_sw)/sizeof(double); + struct htrdr_sun* sun = NULL; + size_t i; + res_T res = RES_OK; + ASSERT(htrdr && out_sun); + ASSERT(sizeof(wavenumbers_sw)/sizeof(double) == nspectral_intervals+1); + + sun = MEM_CALLOC(htrdr->allocator, 1, sizeof(*sun)); + if(!sun) { + htrdr_log_err(htrdr, "could not allocate the sun data structure.\n"); + res = RES_MEM_ERR; + goto error; + } + ref_init(&sun->ref); + sun->htrdr = htrdr; + darray_double_init(htrdr->allocator, &sun->radiance_sw); + darray_double_init(htrdr->allocator, &sun->wavenumbers_sw); + sun->half_angle = 4.6524e3; + sun->solid_angle = 2*PI*(1-cos(sun->half_angle)); + + res = darray_double_resize(&sun->radiance_sw, nspectral_intervals); + if(res != RES_OK) { + htrdr_log_err(htrdr, + "could not allocate the list of per spectral band radiance of the sun.\n"); + goto error; + } + res = darray_double_resize(&sun->wavenumbers_sw, nspectral_intervals+1); + if(res != RES_OK) { + htrdr_log_err(htrdr, + "could not allocate the list of spectral band boundaries of the sun.\n"); + goto error; + } + + FOR_EACH(i, 0, darray_double_size_get(&sun->radiance_sw)) { + /* Convert the incoming flux in radiance */ + darray_double_data_get(&sun->radiance_sw)[i] = incoming_flux_sw[i] / PI; + } + FOR_EACH(i, 0, darray_double_size_get(&sun->wavenumbers_sw)) { + darray_double_data_get(&sun->wavenumbers_sw)[i] = wavenumbers_sw[i]; + } + + sun->main_dir[0] = 0; + sun->main_dir[1] = 0; + sun->main_dir[2] = 1; + +exit: + *out_sun = sun; + return res; +error: + if(sun) { + htrdr_sun_ref_put(sun); + sun = NULL; + } + goto exit; +} + +void +htrdr_sun_ref_get(struct htrdr_sun* sun) +{ + ASSERT(sun); + ref_get(&sun->ref); +} + +void +htrdr_sun_ref_put(struct htrdr_sun* sun) +{ + ASSERT(sun); + ref_put(&sun->ref, release_sun); +} + +res_T +htrdr_sun_set_direction(struct htrdr_sun* sun, const double dir[3]) +{ + double len; + ASSERT(sun && dir); + len = d3_normalize(sun->main_dir, dir); + return len <= 0 ? RES_BAD_ARG : RES_OK; +} + +double +htrdr_sun_get_solid_angle(const struct htrdr_sun* sun) +{ + ASSERT(sun); + return sun->solid_angle; +} + +int +htrdr_sun_is_dir_in_solar_cone(const struct htrdr_sun* sun, const double dir[3]) +{ + double alpha; + ASSERT(sun && dir && d3_is_normalized(dir)); + ASSERT(d3_is_normalized(sun->main_dir)); + alpha = acos(d3_dot(dir, sun->main_dir)); + return alpha < sun->half_angle; +} diff --git a/src/htrdr_sun.h b/src/htrdr_sun.h @@ -0,0 +1,61 @@ +/* Copyright (C) 2018 Université Paul Sabatier, |Meso|Star> + * + * 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/>. */ + +#ifndef HTRDR_SUN_H +#define HTRDR_SUN_H + +#include <rsys/rsys.h> + +/* Forward declaration */ +struct htrdr; +struct htrdr_sun; +struct ssp_rng; + +extern LOCAL_SYM res_T +htrdr_sun_create + (struct htrdr* htrdr, + struct htrdr_sun** out_sun); + +extern LOCAL_SYM void +htrdr_sun_ref_get + (struct htrdr_sun* sun); + +extern LOCAL_SYM void +htrdr_sun_ref_put + (struct htrdr_sun* sun); + +/* Setup the direction *toward* the sun "center" */ +extern LOCAL_SYM res_T +htrdr_sun_set_direction + (struct htrdr_sun* sun, + const double direction[3]); + +/* Return a direction that points *toward* the sun */ +extern LOCAL_SYM res_T +htrdr_sun_sample_direction + (struct htrdr_sun* sun, + struct ssp_rng* rng, + double dir[3]); + +extern LOCAL_SYM double +htrdr_sun_get_solid_angle + (const struct htrdr_sun* sun); + +extern LOCAL_SYM int +htrdr_sun_is_dir_in_solar_cone + (const struct htrdr_sun* sun, + const double dir[3]); + +#endif /* HTRDR_SUN_H */