1 module ray_aabb_intersection; 2 3 import Math = math; 4 5 import std.math.traits: isNaN; 6 7 /* 8 * The MIT License 9 * 10 * Copyright (c) 2015-2021 Kai Burjack 11 %#%# Translated by jordan4ibanez 12 * 13 * Permission is hereby granted, free of charge, to any person obtaining a copy 14 * of this software and associated documentation files (the "Software"), to deal 15 * in the Software without restriction, including without limitation the rights 16 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell 17 * copies of the Software, and to permit persons to whom the Software is 18 * furnished to do so, subject to the following conditions: 19 * 20 * The above copyright notice and this permission notice shall be included in 21 * all copies or substantial portions of the Software. 22 * 23 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 24 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 25 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE 26 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 27 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, 28 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN 29 * THE SOFTWARE. 30 */ 31 32 /** 33 * This is an implementation of the <a 34 * href="http://www.cg.cs.tu-bs.de/media/publications/fast-rayaxis-aligned-bounding-box-overlap-tests-using-ray-slopes.pdf">Fast Ray/Axis-Aligned Bounding Box 35 * Overlap Tests using Ray Slopes</a> paper. 36 * <p> 37 * It is an efficient implementation when testing many axis-aligned boxes against the same ray. 38 * <p> 39 * This class is thread-safe and can be used in a multithreaded environment when testing many axis-aligned boxes against the same ray concurrently. 40 * 41 * @author Kai Burjack 42 */ 43 struct RayAabIntersection { 44 private float originX, originY, originZ; 45 private float dirX, dirY, dirZ; 46 47 /* Needed for ray slope intersection method */ 48 private float c_xy, c_yx, c_zy, c_yz, c_xz, c_zx; 49 private float s_xy, s_yx, s_zy, s_yz, s_xz, s_zx; 50 private byte classification; 51 52 /** 53 * Create a new {@link RayAabIntersection} and initialize it with a ray with origin <code>(originX, originY, originZ)</code> 54 * and direction <code>(dirX, dirY, dirZ)</code>. 55 * <p> 56 * In order to change the direction and/or origin of the ray later, use {@link #set(float, float, float, float, float, float) set()}. 57 * 58 * @see #set(float, float, float, float, float, float) 59 * 60 * @param originX 61 * the x coordinate of the origin 62 * @param originY 63 * the y coordinate of the origin 64 * @param originZ 65 * the z coordinate of the origin 66 * @param dirX 67 * the x coordinate of the direction 68 * @param dirY 69 * the y coordinate of the direction 70 * @param dirZ 71 * the z coordinate of the direction 72 */ 73 this(float originX, float originY, float originZ, float dirX, float dirY, float dirZ) { 74 set(originX, originY, originZ, dirX, dirY, dirZ); 75 } 76 77 /** 78 * Update the ray stored by this {@link RayAabIntersection} with the new origin <code>(originX, originY, originZ)</code> 79 * and direction <code>(dirX, dirY, dirZ)</code>. 80 * 81 * @param originX 82 * the x coordinate of the ray origin 83 * @param originY 84 * the y coordinate of the ray origin 85 * @param originZ 86 * the z coordinate of the ray origin 87 * @param dirX 88 * the x coordinate of the ray direction 89 * @param dirY 90 * the y coordinate of the ray direction 91 * @param dirZ 92 * the z coordinate of the ray direction 93 */ 94 public void set(float originX, float originY, float originZ, float dirX, float dirY, float dirZ) { 95 this.originX = originX; 96 this.originY = originY; 97 this.originZ = originZ; 98 this.dirX = dirX; 99 this.dirY = dirY; 100 this.dirZ = dirZ; 101 precomputeSlope(); 102 } 103 /** 104 * Precompute the values necessary for the ray slope algorithm. 105 */ 106 private void precomputeSlope() { 107 float invDirX = 1.0f / dirX; 108 float invDirY = 1.0f / dirY; 109 float invDirZ = 1.0f / dirZ; 110 s_yx = dirX * invDirY; 111 s_xy = dirY * invDirX; 112 s_zy = dirY * invDirZ; 113 s_yz = dirZ * invDirY; 114 s_xz = dirZ * invDirX; 115 s_zx = dirX * invDirZ; 116 c_xy = originY - s_xy * originX; 117 c_yx = originX - s_yx * originY; 118 c_zy = originY - s_zy * originZ; 119 c_yz = originZ - s_yz * originY; 120 c_xz = originZ - s_xz * originX; // <- original paper had a bug here. It switched originZ/originX 121 c_zx = originX - s_zx * originZ; // <- original paper had a bug here. It switched originZ/originX 122 int sgnX = cast(int)Math.math_signum(dirX); 123 int sgnY = cast(int)Math.math_signum(dirY); 124 int sgnZ = cast(int)Math.math_signum(dirZ); 125 classification = cast(byte) ((sgnZ+1) << 4 | (sgnY+1) << 2 | (sgnX+1)); 126 } 127 128 /** 129 * Test whether the ray stored in this {@link RayAabIntersection} intersect the axis-aligned box 130 * given via its minimum corner <code>(minX, minY, minZ)</code> and its maximum corner <code>(maxX, maxY, maxZ)</code>. 131 * <p> 132 * This implementation uses a tableswitch to dispatch to the correct intersection method. 133 * <p> 134 * This method is thread-safe and can be used to test many axis-aligned boxes concurrently. 135 * 136 * @param minX 137 * the x coordinate of the minimum corner 138 * @param minY 139 * the y coordinate of the minimum corner 140 * @param minZ 141 * the z coordinate of the minimum corner 142 * @param maxX 143 * the x coordinate of the maximum corner 144 * @param maxY 145 * the y coordinate of the maximum corner 146 * @param maxZ 147 * the z coordinate of the maximum corner 148 * @return <code>true</code> iff the ray intersects the given axis-aligned box; <code>false</code> otherwise 149 */ 150 public bool test(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 151 // tableswitch with dense and consecutive cases (will be a simple jump based on the switch argument) 152 switch (classification) { 153 case 0: // 0b000000: // MMM 154 return MMM(minX, minY, minZ, maxX, maxY, maxZ); 155 case 1: // 0b000001: // OMM 156 return OMM(minX, minY, minZ, maxX, maxY, maxZ); 157 case 2: // 0b000010: // PMM 158 return PMM(minX, minY, minZ, maxX, maxY, maxZ); 159 case 3: // 0b000011: // not used 160 return false; 161 case 4: // 0b000100: // MOM 162 return MOM(minX, minY, minZ, maxX, maxY, maxZ); 163 case 5: // 0b000101: // OOM 164 return OOM(minX, minY, minZ, maxX, maxY); 165 case 6: // 0b000110: // POM 166 return POM(minX, minY, minZ, maxX, maxY, maxZ); 167 case 7: // 0b000111: // not used 168 return false; 169 case 8: // 0b001000: // MPM 170 return MPM(minX, minY, minZ, maxX, maxY, maxZ); 171 case 9: // 0b001001: // OPM 172 return OPM(minX, minY, minZ, maxX, maxY, maxZ); 173 case 10: // 0b001010: // PPM 174 return PPM(minX, minY, minZ, maxX, maxY, maxZ); 175 case 11: // 0b001011: // not used 176 case 12: // 0b001100: // not used 177 case 13: // 0b001101: // not used 178 case 14: // 0b001110: // not used 179 case 15: // 0b001111: // not used 180 return false; 181 case 16: // 0b010000: // MMO 182 return MMO(minX, minY, minZ, maxX, maxY, maxZ); 183 case 17: // 0b010001: // OMO 184 return OMO(minX, minY, minZ, maxX, maxZ); 185 case 18: // 0b010010: // PMO 186 return PMO(minX, minY, minZ, maxX, maxY, maxZ); 187 case 19: // 0b010011: // not used 188 return false; 189 case 20: // 0b010100: // MOO 190 return MOO(minX, minY, minZ, maxY, maxZ); 191 case 21: // 0b010101: // OOO 192 return false; // <- degenerate case 193 case 22: // 0b010110: // POO 194 return POO(minY, minZ, maxX, maxY, maxZ); 195 case 23: // 0b010111: // not used 196 return false; 197 case 24: // 0b011000: // MPO 198 return MPO(minX, minY, minZ, maxX, maxY, maxZ); 199 case 25: // 0b011001: // OPO 200 return OPO(minX, minZ, maxX, maxY, maxZ); 201 case 26: // 0b011010: // PPO 202 return PPO(minX, minY, minZ, maxX, maxY, maxZ); 203 case 27: // 0b011011: // not used 204 case 28: // 0b011100: // not used 205 case 29: // 0b011101: // not used 206 case 30: // 0b011110: // not used 207 case 31: // 0b011111: // not used 208 return false; 209 case 32: // 0b100000: // MMP 210 return MMP(minX, minY, minZ, maxX, maxY, maxZ); 211 case 33: // 0b100001: // OMP 212 return OMP(minX, minY, minZ, maxX, maxY, maxZ); 213 case 34: // 0b100010: // PMP 214 return PMP(minX, minY, minZ, maxX, maxY, maxZ); 215 case 35: // 0b100011: // not used 216 return false; 217 case 36: // 0b100100: // MOP 218 return MOP(minX, minY, minZ, maxX, maxY, maxZ); 219 case 37: // 0b100101: // OOP 220 return OOP(minX, minY, maxX, maxY, maxZ); 221 case 38: // 0b100110: // POP 222 return POP(minX, minY, minZ, maxX, maxY, maxZ); 223 case 39: // 0b100111: // not used 224 return false; 225 case 40: // 0b101000: // MPP 226 return MPP(minX, minY, minZ, maxX, maxY, maxZ); 227 case 41: // 0b101001: // OPP 228 return OPP(minX, minY, minZ, maxX, maxY, maxZ); 229 case 42: // 0b101010: // PPP 230 return PPP(minX, minY, minZ, maxX, maxY, maxZ); 231 default: 232 return false; 233 } 234 } 235 236 /* Intersection tests for all possible ray direction cases */ 237 238 private bool MMM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 239 return originX >= minX && originY >= minY && originZ >= minZ 240 && s_xy * minX - maxY + c_xy <= 0.0f 241 && s_yx * minY - maxX + c_yx <= 0.0f 242 && s_zy * minZ - maxY + c_zy <= 0.0f 243 && s_yz * minY - maxZ + c_yz <= 0.0f 244 && s_xz * minX - maxZ + c_xz <= 0.0f 245 && s_zx * minZ - maxX + c_zx <= 0.0f; 246 } 247 private bool OMM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 248 return originX >= minX && originX <= maxX && originY >= minY && originZ >= minZ 249 && s_zy * minZ - maxY + c_zy <= 0.0f 250 && s_yz * minY - maxZ + c_yz <= 0.0f; 251 } 252 private bool PMM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 253 return originX <= maxX && originY >= minY && originZ >= minZ 254 && s_xy * maxX - maxY + c_xy <= 0.0f 255 && s_yx * minY - minX + c_yx >= 0.0f 256 && s_zy * minZ - maxY + c_zy <= 0.0f 257 && s_yz * minY - maxZ + c_yz <= 0.0f 258 && s_xz * maxX - maxZ + c_xz <= 0.0f 259 && s_zx * minZ - minX + c_zx >= 0.0f; 260 } 261 private bool MOM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 262 return originY >= minY && originY <= maxY && originX >= minX && originZ >= minZ 263 && s_xz * minX - maxZ + c_xz <= 0.0f 264 && s_zx * minZ - maxX + c_zx <= 0.0f; 265 } 266 private bool OOM(float minX, float minY, float minZ, float maxX, float maxY) { 267 return originZ >= minZ && originX >= minX && originX <= maxX && originY >= minY && originY <= maxY; 268 } 269 private bool POM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 270 return originY >= minY && originY <= maxY && originX <= maxX && originZ >= minZ 271 && s_xz * maxX - maxZ + c_xz <= 0.0f 272 && s_zx * minZ - minX + c_zx >= 0.0f; 273 } 274 private bool MPM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 275 return originX >= minX && originY <= maxY && originZ >= minZ 276 && s_xy * minX - minY + c_xy >= 0.0f 277 && s_yx * maxY - maxX + c_yx <= 0.0f 278 && s_zy * minZ - minY + c_zy >= 0.0f 279 && s_yz * maxY - maxZ + c_yz <= 0.0f 280 && s_xz * minX - maxZ + c_xz <= 0.0f 281 && s_zx * minZ - maxX + c_zx <= 0.0f; 282 } 283 private bool OPM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 284 return originX >= minX && originX <= maxX && originY <= maxY && originZ >= minZ 285 && s_zy * minZ - minY + c_zy >= 0.0f 286 && s_yz * maxY - maxZ + c_yz <= 0.0f; 287 } 288 private bool PPM(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 289 return originX <= maxX && originY <= maxY && originZ >= minZ 290 && s_xy * maxX - minY + c_xy >= 0.0f 291 && s_yx * maxY - minX + c_yx >= 0.0f 292 && s_zy * minZ - minY + c_zy >= 0.0f 293 && s_yz * maxY - maxZ + c_yz <= 0.0f 294 && s_xz * maxX - maxZ + c_xz <= 0.0f 295 && s_zx * minZ - minX + c_zx >= 0.0f; 296 } 297 private bool MMO(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 298 return originZ >= minZ && originZ <= maxZ && originX >= minX && originY >= minY 299 && s_xy * minX - maxY + c_xy <= 0.0f 300 && s_yx * minY - maxX + c_yx <= 0.0f; 301 } 302 private bool OMO(float minX, float minY, float minZ, float maxX, float maxZ) { 303 return originY >= minY && originX >= minX && originX <= maxX && originZ >= minZ && originZ <= maxZ; 304 } 305 private bool PMO(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 306 return originZ >= minZ && originZ <= maxZ && originX <= maxX && originY >= minY 307 && s_xy * maxX - maxY + c_xy <= 0.0f 308 && s_yx * minY - minX + c_yx >= 0.0f; 309 } 310 private bool MOO(float minX, float minY, float minZ, float maxY, float maxZ) { 311 return originX >= minX && originY >= minY && originY <= maxY && originZ >= minZ && originZ <= maxZ; 312 } 313 private bool POO(float minY, float minZ, float maxX, float maxY, float maxZ) { 314 return originX <= maxX && originY >= minY && originY <= maxY && originZ >= minZ && originZ <= maxZ; 315 } 316 private bool MPO(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 317 return originZ >= minZ && originZ <= maxZ && originX >= minX && originY <= maxY 318 && s_xy * minX - minY + c_xy >= 0.0f 319 && s_yx * maxY - maxX + c_yx <= 0.0f; 320 } 321 private bool OPO(float minX, float minZ, float maxX, float maxY, float maxZ) { 322 return originY <= maxY && originX >= minX && originX <= maxX && originZ >= minZ && originZ <= maxZ; 323 } 324 private bool PPO(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 325 return originZ >= minZ && originZ <= maxZ && originX <= maxX && originY <= maxY 326 && s_xy * maxX - minY + c_xy >= 0.0f 327 && s_yx * maxY - minX + c_yx >= 0.0f; 328 } 329 private bool MMP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 330 return originX >= minX && originY >= minY && originZ <= maxZ 331 && s_xy * minX - maxY + c_xy <= 0.0f 332 && s_yx * minY - maxX + c_yx <= 0.0f 333 && s_zy * maxZ - maxY + c_zy <= 0.0f 334 && s_yz * minY - minZ + c_yz >= 0.0f 335 && s_xz * minX - minZ + c_xz >= 0.0f 336 && s_zx * maxZ - maxX + c_zx <= 0.0f; 337 } 338 private bool OMP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 339 return originX >= minX && originX <= maxX && originY >= minY && originZ <= maxZ 340 && s_zy * maxZ - maxY + c_zy <= 0.0f 341 && s_yz * minY - minZ + c_yz >= 0.0f; 342 } 343 private bool PMP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 344 return originX <= maxX && originY >= minY && originZ <= maxZ 345 && s_xy * maxX - maxY + c_xy <= 0.0f 346 && s_yx * minY - minX + c_yx >= 0.0f 347 && s_zy * maxZ - maxY + c_zy <= 0.0f 348 && s_yz * minY - minZ + c_yz >= 0.0f 349 && s_xz * maxX - minZ + c_xz >= 0.0f 350 && s_zx * maxZ - minX + c_zx >= 0.0f; 351 } 352 private bool MOP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 353 return originY >= minY && originY <= maxY && originX >= minX && originZ <= maxZ 354 && s_xz * minX - minZ + c_xz >= 0.0f 355 && s_zx * maxZ - maxX + c_zx <= 0.0f; 356 } 357 private bool OOP(float minX, float minY, float maxX, float maxY, float maxZ) { 358 return originZ <= maxZ && originX >= minX && originX <= maxX && originY >= minY && originY <= maxY; 359 } 360 private bool POP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 361 return originY >= minY && originY <= maxY && originX <= maxX && originZ <= maxZ 362 && s_xz * maxX - minZ + c_xz >= 0.0f 363 && s_zx * maxZ - minX + c_zx <= 0.0f; 364 } 365 private bool MPP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 366 return originX >= minX && originY <= maxY && originZ <= maxZ 367 && s_xy * minX - minY + c_xy >= 0.0f 368 && s_yx * maxY - maxX + c_yx <= 0.0f 369 && s_zy * maxZ - minY + c_zy >= 0.0f 370 && s_yz * maxY - minZ + c_yz >= 0.0f 371 && s_xz * minX - minZ + c_xz >= 0.0f 372 && s_zx * maxZ - maxX + c_zx <= 0.0f; 373 } 374 private bool OPP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 375 return originX >= minX && originX <= maxX && originY <= maxY && originZ <= maxZ 376 && s_zy * maxZ - minY + c_zy <= 0.0f 377 && s_yz * maxY - minZ + c_yz <= 0.0f; 378 } 379 private bool PPP(float minX, float minY, float minZ, float maxX, float maxY, float maxZ) { 380 return originX <= maxX && originY <= maxY && originZ <= maxZ 381 && s_xy * maxX - minY + c_xy >= 0.0f 382 && s_yx * maxY - minX + c_yx >= 0.0f 383 && s_zy * maxZ - minY + c_zy >= 0.0f 384 && s_yz * maxY - minZ + c_yz >= 0.0f 385 && s_xz * maxX - minZ + c_xz >= 0.0f 386 && s_zx * maxZ - minX + c_zx >= 0.0f; 387 } 388 }