test_ssp_ran_triangle.h (4243B)
1 /* Copyright (C) 2015-2025 |Méso|Star> (contact@meso-star.com) 2 * 3 * This program is free software: you can redistribute it and/or modify 4 * it under the terms of the GNU General Public License as published by 5 * the Free Software Foundation, either version 3 of the License, or 6 * (at your option) any later version. 7 * 8 * This program is distributed in the hope that it will be useful, 9 * but WITHOUT ANY WARRANTY; without even the implied warranty of 10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11 * GNU General Public License for more details. 12 * 13 * You should have received a copy of the GNU General Public License 14 * along with this program. If not, see <http://www.gnu.org/licenses/>. */ 15 16 #ifndef TEST_SSP_RAN_TRIANGLE_H 17 #define TEST_SSP_RAN_TRIANGLE_H 18 19 #include "ssp.h" 20 #include "test_ssp_utils.h" 21 22 #include <rsys/float4.h> 23 24 #define NSAMPS 128 25 26 #endif /* TEST_SSP_RAN_TRIANGLE_H */ 27 28 #if TYPE_FLOAT==0 29 #define REAL double 30 #define TEST test_double 31 #define RNG_UNIFORM_R ssp_rng_uniform_double 32 #define RAN_TRIANGLE_UNIFORM ssp_ran_triangle_uniform 33 #define RAN_TRIANGLE_UNIFORM_PDF ssp_ran_triangle_uniform_pdf 34 #define EQ_EPS_R eq_eps 35 #define R3 d3 36 #define R3_DOT d3_dot 37 #define R3_SUB d3_sub 38 #define R3_MINUS d3_minus 39 #define R3_CROSS d3_cross 40 #define R3_LEN d3_len 41 42 #elif TYPE_FLOAT==1 43 #define REAL float 44 #define TEST test_float 45 #define RNG_UNIFORM_R ssp_rng_uniform_float 46 #define RAN_TRIANGLE_UNIFORM ssp_ran_triangle_uniform_float 47 #define RAN_TRIANGLE_UNIFORM_PDF ssp_ran_triangle_uniform_float_pdf 48 #define EQ_EPS_R eq_epsf 49 #define R3 f3 50 #define R3_DOT f3_dot 51 #define R3_SUB f3_sub 52 #define R3_MINUS f3_minus 53 #define R3_CROSS f3_cross 54 #define R3_LEN f3_len 55 56 #else 57 #error "TYPE_FLOAT must be defined either 0 or 1" 58 #endif 59 60 static void 61 TEST() 62 { 63 struct ssp_rng* rng; 64 struct mem_allocator allocator; 65 REAL samps[NSAMPS][3]; 66 REAL A[3], B[3], C[3]; 67 REAL v0[3], v1[3], v2[3], m0[3], m1[3], m2[3]; 68 REAL plane[4]; 69 REAL pdf; 70 size_t counter[2]; 71 size_t nsteps; 72 size_t i; 73 74 mem_init_proxy_allocator(&allocator, &mem_default_allocator); 75 76 CHK(ssp_rng_create(&allocator, SSP_RNG_MT19937_64, &rng) == RES_OK); 77 78 FOR_EACH(i, 0, 3) { 79 A[i] = RNG_UNIFORM_R(rng, 0, 1); 80 B[i] = RNG_UNIFORM_R(rng, 0, 1); 81 C[i] = RNG_UNIFORM_R(rng, 0, 1); 82 } 83 84 R3_SUB(v0, B, A); 85 R3_SUB(v1, C, A); 86 R3_SUB(v2, C, B); 87 R3_MINUS(m0, v0); 88 R3_MINUS(m1, v1); 89 R3_MINUS(m2, v2); 90 R3_CROSS(plane, v0, v1); 91 plane[3] = -R3_DOT(plane, C); 92 93 FOR_EACH(i, 0, NSAMPS) { 94 REAL tmp0[3], tmp1[3]; 95 REAL dot = 0; 96 REAL area = 0; 97 98 CHK(RAN_TRIANGLE_UNIFORM(rng, A, B, C, samps[i], &pdf) == samps[i]); 99 CHK(EQ_EPS_R(R3_DOT(plane, samps[i]), -plane[3], (REAL)1.e-6) == 1); 100 101 R3_SUB(tmp0, samps[i], A); 102 dot = R3_DOT(R3_CROSS(tmp0, tmp0, v0), R3_CROSS(tmp1, v1, v0)); 103 CHK(sign(dot) == 1); 104 R3_SUB(tmp0, samps[i], B); 105 dot = R3_DOT(R3_CROSS(tmp0, tmp0, v2), R3_CROSS(tmp1, m0, v2)); 106 CHK(sign(dot) == 1); 107 R3_SUB(tmp0, samps[i], C); 108 dot = R3_DOT(R3_CROSS(tmp0, tmp0, m1), R3_CROSS(tmp1, m2, m1)); 109 CHK(sign(dot) == 1); 110 111 area = R3_LEN(tmp1) * 0.5f; 112 CHK(EQ_EPS_R(pdf, RAN_TRIANGLE_UNIFORM_PDF(A, B, C), (REAL)1.e-8) == 1); 113 CHK(EQ_EPS_R(1 / area, pdf, (REAL)1.e-6) == 1); 114 } 115 116 nsteps = 10000; 117 counter[0] = counter[1] = 0; 118 R3(A, -1, 0, 0); 119 R3(B, 1, 0, 0); 120 R3(C, 0, 1, 0); 121 FOR_EACH(i, 0, nsteps) { 122 RAN_TRIANGLE_UNIFORM(rng, A, B, C, samps[0], NULL); 123 if(samps[0][0] < 0.0) 124 counter[0] += 1; 125 else 126 counter[1] += 1; 127 } 128 CHK(labs((long)(counter[1] - counter[0])) < 100); 129 130 counter[0] = counter[1] = 0; 131 FOR_EACH(i, 0, nsteps) { 132 RAN_TRIANGLE_UNIFORM(rng, A, B, C, samps[0], NULL); 133 if(samps[0][1] < 1 - 1/sqrt(2)) 134 counter[0] += 1; 135 else 136 counter[1] += 1; 137 } 138 CHK(labs((long)(counter[1] - counter[0])) < 100); 139 140 ssp_rng_ref_put(rng); 141 check_memory_allocator(&allocator); 142 mem_shutdown_proxy_allocator(&allocator); 143 CHK(mem_allocated_size() == 0); 144 } 145 146 #undef REAL 147 #undef TEST 148 #undef RNG_UNIFORM_R 149 #undef RAN_TRIANGLE_UNIFORM 150 #undef RAN_TRIANGLE_UNIFORM_PDF 151 #undef EQ_EPS_R 152 #undef R3 153 #undef R3_DOT 154 #undef R3_SUB 155 #undef R3_MINUS 156 #undef R3_CROSS 157 #undef R3_LEN 158 #undef TYPE_FLOAT