ssol_mc_receiver.c (6907B)
1 /* Copyright (C) 2018-2026 |Meso|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_estimator_c.h" 19 #include "ssol_object_c.h" 20 21 #include <rsys/double3.h> 22 #include <star/s3d.h> 23 24 #ifdef COMPILER_CL 25 #pragma warning(push) 26 #pragma warning(disable:4706) /* Assignment within a condition */ 27 #endif 28 29 /******************************************************************************* 30 * Exported functions 31 ******************************************************************************/ 32 res_T 33 ssol_estimator_get_mc_receiver 34 (struct ssol_estimator* estimator, 35 const struct ssol_instance* instance, 36 const enum ssol_side_flag side, 37 struct ssol_mc_receiver* rcv) 38 { 39 struct mc_receiver* mc_rcv = NULL; 40 struct mc_receiver_1side* mc_rcv1 = NULL; 41 42 if(!estimator || !instance || !rcv 43 || !(instance->receiver_mask & (int)side)) 44 return RES_BAD_ARG; 45 46 memset(rcv, 0, sizeof(rcv[0])); 47 48 mc_rcv = htable_receiver_find(&estimator->mc_receivers, &instance); 49 if(!mc_rcv) { 50 /* The receiver has no MC estimation */ 51 return RES_OK; 52 } 53 54 mc_rcv1 = side == SSOL_FRONT ? &mc_rcv->front : &mc_rcv->back; 55 #define SETUP_MC_RESULT(Name) { \ 56 const double N = (double)estimator->realisation_count; \ 57 struct mc_data* data = &mc_rcv1->Name; \ 58 double weight, sqr_weight; \ 59 mc_data_get(data, &weight, &sqr_weight); \ 60 rcv->Name.E = weight / N; \ 61 rcv->Name.V = sqr_weight/N - rcv->Name.E*rcv->Name.E; \ 62 rcv->Name.V = rcv->Name.V > 0 ? rcv->Name.V : 0; \ 63 rcv->Name.SE = sqrt(rcv->Name.V / N); \ 64 } (void)0 65 #define MC_SETUP_ALL { \ 66 SETUP_MC_RESULT(incoming_flux); \ 67 SETUP_MC_RESULT(incoming_if_no_atm_loss); \ 68 SETUP_MC_RESULT(incoming_if_no_field_loss); \ 69 SETUP_MC_RESULT(incoming_lost_in_atmosphere); \ 70 SETUP_MC_RESULT(incoming_lost_in_field); \ 71 SETUP_MC_RESULT(absorbed_flux); \ 72 SETUP_MC_RESULT(absorbed_if_no_atm_loss); \ 73 SETUP_MC_RESULT(absorbed_if_no_field_loss); \ 74 SETUP_MC_RESULT(absorbed_lost_in_atmosphere); \ 75 SETUP_MC_RESULT(absorbed_lost_in_field); \ 76 } (void)0 77 MC_SETUP_ALL; 78 #undef SETUP_MC_RESULT 79 rcv->mc__ = mc_rcv1; 80 rcv->N__ = estimator->realisation_count; 81 rcv->instance__ = instance; 82 return RES_OK; 83 } 84 85 res_T 86 ssol_mc_receiver_get_mc_shape 87 (struct ssol_mc_receiver* rcv, 88 const struct ssol_shape* shape, 89 struct ssol_mc_shape* mc) 90 { 91 struct mc_receiver_1side* mc_rcv1; 92 93 if(!rcv || !shape || !mc) return RES_BAD_ARG; 94 if(!object_has_shape(rcv->instance__->object, shape)) return RES_BAD_ARG; 95 mc_rcv1 = rcv->mc__; 96 mc->N__ = rcv->N__; 97 mc->mc__ = htable_shape2mc_find(&mc_rcv1->shape2mc, &shape); 98 mc->shape__ = shape; 99 return RES_OK; 100 } 101 102 res_T 103 ssol_mc_shape_get_mc_primitive 104 (struct ssol_mc_shape* shape, 105 const unsigned i, 106 struct ssol_mc_primitive* prim) 107 { 108 struct mc_shape_1side* mc_shape1; 109 struct mc_primitive_1side* mc_prim1; 110 unsigned ntris; 111 112 if(!shape || !prim) return RES_BAD_ARG; 113 114 SSOL(shape_get_triangles_count(shape->shape__, &ntris)); 115 if(i >= ntris) return RES_BAD_ARG; 116 117 mc_shape1 = shape->mc__; 118 if(!mc_shape1 || !(mc_prim1 = htable_prim2mc_find(&mc_shape1->prim2mc, &i))) { 119 #define SETUP_MC_RESULT(Name) { \ 120 prim->Name.E = 0; \ 121 prim->Name.V = 0; \ 122 prim->Name.SE = 0; \ 123 } (void)0 124 MC_SETUP_ALL; 125 #undef SETUP_MC_RESULT 126 } else { 127 struct s3d_attrib attr; 128 struct s3d_shape* s3d_shape; 129 double v0[3], v1[3], v2[3], E0[3], E1[3], normal[3]; 130 double area; 131 unsigned ids[3]; 132 res_T res = RES_OK; 133 134 s3d_shape = shape->shape__->shape_rt; 135 136 /* Retrieve the primitive indices */ 137 res = s3d_mesh_get_triangle_indices(s3d_shape, i, ids); 138 if(res != RES_OK) return res; 139 140 /* Fetch the primitive vertices */ 141 S3D(mesh_get_vertex_attrib(s3d_shape, ids[0], S3D_POSITION, &attr)); 142 d3_set_f3(v0, attr.value); 143 S3D(mesh_get_vertex_attrib(s3d_shape, ids[1], S3D_POSITION, &attr)); 144 d3_set_f3(v1, attr.value); 145 S3D(mesh_get_vertex_attrib(s3d_shape, ids[2], S3D_POSITION, &attr)); 146 d3_set_f3(v2, attr.value); 147 148 /* Compute the primitive area */ 149 d3_sub(E0, v1, v0); 150 d3_sub(E1, v2, v0); 151 d3_cross(normal, E0, E1); 152 area = d3_len(normal) * 0.5; 153 154 #define SETUP_MC_RESULT(Name) { \ 155 const double N = (double)shape->N__; \ 156 struct mc_data* data = &mc_prim1->Name; \ 157 double weight, sqr_weight; \ 158 mc_data_get(data, &weight, &sqr_weight); \ 159 prim->Name.E = weight / N; \ 160 prim->Name.V = sqr_weight/N - prim->Name.E*prim->Name.E; \ 161 prim->Name.V = prim->Name.V > 0 ? prim->Name.V : 0; \ 162 prim->Name.SE = sqrt(prim->Name.V / N); \ 163 prim->Name.E /= area; \ 164 prim->Name.V /= area*area; \ 165 prim->Name.SE /= area; \ 166 } (void)0 167 MC_SETUP_ALL; 168 #undef SETUP_MC_RESULT 169 #undef MC_SETUP_ALL 170 } 171 172 return RES_OK; 173 } 174 175 #ifdef COMPILER_CL 176 #pragma warning(pop) 177 #endif 178