htrdr

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

htrdr_rectangle.c (6268B)


      1 /* Copyright (C) 2018-2019, 2022-2025 Centre National de la Recherche Scientifique
      2  * Copyright (C) 2020-2022 Institut Mines Télécom Albi-Carmaux
      3  * Copyright (C) 2022-2025 Institut Pierre-Simon Laplace
      4  * Copyright (C) 2022-2025 Institut de Physique du Globe de Paris
      5  * Copyright (C) 2018-2025 |Méso|Star> (contact@meso-star.com)
      6  * Copyright (C) 2022-2025 Observatoire de Paris
      7  * Copyright (C) 2022-2025 Université de Reims Champagne-Ardenne
      8  * Copyright (C) 2022-2025 Université de Versaille Saint-Quentin
      9  * Copyright (C) 2018-2019, 2022-2025 Université Paul Sabatier
     10  *
     11  * This program is free software: you can redistribute it and/or modify
     12  * it under the terms of the GNU General Public License as published by
     13  * the Free Software Foundation, either version 3 of the License, or
     14  * (at your option) any later version.
     15  *
     16  * This program is distributed in the hope that it will be useful,
     17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
     18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
     19  * GNU General Public License for more details.
     20  *
     21  * You should have received a copy of the GNU General Public License
     22  * along with this program. If not, see <http://www.gnu.org/licenses/>. */
     23 
     24 #include "core/htrdr.h"
     25 #include "core/htrdr_log.h"
     26 #include "core/htrdr_rectangle.h"
     27 
     28 #include <rsys/double2.h>
     29 #include <rsys/double3.h>
     30 #include <rsys/double33.h>
     31 #include <rsys/mem_allocator.h>
     32 #include <rsys/ref_count.h>
     33 
     34 struct htrdr_rectangle {
     35   /* Frame of the rectangle in world space */
     36   double axis_x[3];
     37   double axis_y[3];
     38   double normal[3];
     39 
     40   double size[2];
     41 
     42   double local2world[12]; /* Rectangle to world transformation matrix */
     43   double world2local[12]; /* World to rectangle transformation matrix */
     44 
     45   double position[3]; /* Center of the rectangle */
     46   struct htrdr* htrdr;
     47   ref_T ref;
     48 };
     49 
     50 /*******************************************************************************
     51  * Helper functions
     52  ******************************************************************************/
     53 static void
     54 rectangle_release(ref_T* ref)
     55 {
     56   struct htrdr_rectangle* rect;
     57   struct htrdr* htrdr;
     58   ASSERT(ref);
     59   rect = CONTAINER_OF(ref, struct htrdr_rectangle, ref);
     60   htrdr = rect->htrdr;
     61   MEM_RM(htrdr_get_allocator(htrdr), rect);
     62   htrdr_ref_put(htrdr);
     63 }
     64 
     65 /*******************************************************************************
     66  * Exported functions
     67  ******************************************************************************/
     68 res_T
     69 htrdr_rectangle_create
     70   (struct htrdr* htrdr,
     71    const double sz[2],
     72    const double pos[3],
     73    const double tgt[3],
     74    const double up[3],
     75    struct htrdr_rectangle** out_rect)
     76 {
     77   struct htrdr_rectangle* rect = NULL;
     78   double x[3], y[3], z[3];
     79   double trans[3];
     80   res_T res = RES_OK;
     81   ASSERT(htrdr && pos && tgt && up && sz && out_rect);
     82 
     83   rect = MEM_CALLOC(htrdr_get_allocator(htrdr), 1, sizeof(*rect));
     84   if(!rect) {
     85     htrdr_log_err(htrdr, "Could not allocate the rectangle data structure.\n");
     86     res = RES_MEM_ERR;
     87     goto error;
     88   }
     89   ref_init(&rect->ref);
     90   htrdr_ref_get(htrdr);
     91   rect->htrdr = htrdr;
     92 
     93   if(sz[0] <= 0 || sz[1] <= 0) {
     94     htrdr_log_err(htrdr,
     95       "Invalid rectangle size `%g %g'. It must be strictly positive.\n",
     96       SPLIT2(sz));
     97     res = RES_BAD_ARG;
     98     goto error;
     99   }
    100 
    101   if(d3_normalize(z, d3_sub(z, tgt, pos)) <= 0
    102   || d3_normalize(x, d3_cross(x, z, up)) <= 0
    103   || d3_normalize(y, d3_cross(y, z, x)) <= 0) {
    104     htrdr_log_err(htrdr, "Invalid rectangle frame:\n"
    105       "\tposition = %g %g %g\n"
    106       "\ttarget = %g %g %g\n"
    107       "\tup = %g %g %g\n",
    108       SPLIT3(pos), SPLIT3(tgt), SPLIT3(up));
    109     res = RES_BAD_ARG;
    110     goto error;
    111   }
    112 
    113   /* Setup the local to world transformation matrix */
    114   d3_set(rect->local2world+0, x);
    115   d3_set(rect->local2world+3, y);
    116   d3_set(rect->local2world+6, z);
    117   d3_set(rect->local2world+9, pos);
    118 
    119   /* Inverse the local to world transformation matrix. Note that since the
    120    * represented frame is orthonormal one can transpose the original rotation
    121    * matrix to cheaply compute its inverse */
    122   d33_transpose(rect->world2local, rect->local2world);
    123 
    124   /* Compute the affine inverse transform */
    125   d3_minus(trans, pos);
    126   d33_muld3(rect->world2local+9, rect->world2local, trans);
    127 
    128   d3_muld(rect->axis_x, x, sz[0]*0.5);
    129   d3_muld(rect->axis_y, y, sz[1]*0.5);
    130   d3_set(rect->normal, z);
    131   d3_set(rect->position, pos);
    132 
    133   d2_set(rect->size, sz);
    134 
    135 exit:
    136   *out_rect = rect;
    137   return res;
    138 error:
    139   if(rect) {
    140     htrdr_rectangle_ref_put(rect);
    141     rect = NULL;
    142   }
    143   goto exit;
    144 }
    145 
    146 void
    147 htrdr_rectangle_sample_pos
    148   (const struct htrdr_rectangle* rect,
    149    const double sample[2], /* In [0, 1[ */
    150    double pos[3])
    151 {
    152   double x[3], y[3];
    153   ASSERT(rect && sample && pos);
    154   d3_muld(x, rect->axis_x, sample[0]*2.0 - 1.0);
    155   d3_muld(y, rect->axis_y, sample[1]*2.0 - 1.0);
    156   d3_add(pos, d3_add(pos, rect->position, x), y);
    157 }
    158 
    159 void
    160 htrdr_rectangle_ref_get(struct htrdr_rectangle* rect)
    161 {
    162   ASSERT(rect);
    163   ref_get(&rect->ref);
    164 }
    165 
    166 void
    167 htrdr_rectangle_ref_put(struct htrdr_rectangle* rect)
    168 {
    169   ASSERT(rect);
    170   ref_put(&rect->ref, rectangle_release);
    171 }
    172 
    173 void
    174 htrdr_rectangle_get_normal(const struct htrdr_rectangle* rect, double normal[3])
    175 {
    176   ASSERT(rect && normal);
    177   d3_set(normal, rect->normal);
    178 }
    179 
    180 void
    181 htrdr_rectangle_get_center(const struct htrdr_rectangle* rect, double pos[3])
    182 {
    183   ASSERT(rect && pos);
    184   d3_set(pos, rect->position);
    185 }
    186 
    187 double*
    188 htrdr_rectangle_get_transform
    189   (const struct htrdr_rectangle* rect,
    190    double transform[12])
    191 {
    192   ASSERT(rect && transform);
    193   d3_set(transform+0, rect->local2world+0);
    194   d3_set(transform+3, rect->local2world+3);
    195   d3_set(transform+6, rect->local2world+6);
    196   d3_set(transform+9, rect->local2world+9);
    197   return transform;
    198 }
    199 
    200 double*
    201 htrdr_rectangle_get_transform_inverse
    202   (const struct htrdr_rectangle* rect,
    203    double transform_inverse[12])
    204 {
    205   ASSERT(rect && transform_inverse);
    206   d3_set(transform_inverse+0, rect->world2local+0);
    207   d3_set(transform_inverse+3, rect->world2local+3);
    208   d3_set(transform_inverse+6, rect->world2local+6);
    209   d3_set(transform_inverse+9, rect->world2local+9);
    210   return transform_inverse;
    211 }
    212 
    213 void
    214 htrdr_rectangle_get_size(const struct htrdr_rectangle* rect, double size[2])
    215 {
    216   ASSERT(rect && size);
    217   d2_set(size, rect->size);
    218 }