diff options
Diffstat (limited to 'source/blender/render/intern/raytrace/rayobject_octree.cpp')
| -rw-r--r-- | source/blender/render/intern/raytrace/rayobject_octree.cpp | 1101 |
1 files changed, 1101 insertions, 0 deletions
diff --git a/source/blender/render/intern/raytrace/rayobject_octree.cpp b/source/blender/render/intern/raytrace/rayobject_octree.cpp new file mode 100644 index 00000000000..4b73e64ca45 --- /dev/null +++ b/source/blender/render/intern/raytrace/rayobject_octree.cpp @@ -0,0 +1,1101 @@ +/* + * ***** BEGIN GPL LICENSE BLOCK ***** + * + * 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 2 + * 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, write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + * + * The Original Code is Copyright (C) 1990-1998 NeoGeo BV. + * All rights reserved. + * + * Contributors: 2004/2005 Blender Foundation, full recode + * + * ***** END GPL LICENSE BLOCK ***** + */ + +/** \file blender/render/intern/raytrace/rayobject_octree.cpp + * \ingroup render + */ + + +/* IMPORTANT NOTE: this code must be independent of any other render code + * to use it outside the renderer! */ + +#include <math.h> +#include <string.h> +#include <stdlib.h> +#include <float.h> +#include <assert.h> + +#include "MEM_guardedalloc.h" + +#include "DNA_material_types.h" + +#include "BLI_math.h" +#include "BLI_utildefines.h" + +#include "rayintersection.h" +#include "rayobject.h" + +/* ********** structs *************** */ +#define BRANCH_ARRAY 1024 +#define NODE_ARRAY 4096 + +typedef struct Branch { + struct Branch *b[8]; +} Branch; + +typedef struct OcVal { + short ocx, ocy, ocz; +} OcVal; + +typedef struct Node { + struct RayFace *v[8]; + struct OcVal ov[8]; + struct Node *next; +} Node; + +typedef struct Octree { + RayObject rayobj; + + struct Branch **adrbranch; + struct Node **adrnode; + float ocsize; /* ocsize: mult factor, max size octree */ + float ocfacx, ocfacy, ocfacz; + float min[3], max[3]; + int ocres; + int branchcount, nodecount; + + /* during building only */ + char *ocface; + + RayFace **ro_nodes; + int ro_nodes_size, ro_nodes_used; + +} Octree; + +static int RE_rayobject_octree_intersect(RayObject *o, Isect *isec); +static void RE_rayobject_octree_add(RayObject *o, RayObject *ob); +static void RE_rayobject_octree_done(RayObject *o); +static void RE_rayobject_octree_free(RayObject *o); +static void RE_rayobject_octree_bb(RayObject *o, float *min, float *max); + +/* + * This function is not expected to be called by current code state. + */ +static float RE_rayobject_octree_cost(RayObject *UNUSED(o)) +{ + return 1.0; +} + +static void RE_rayobject_octree_hint_bb(RayObject *UNUSED(o), RayHint *UNUSED(hint), + float *UNUSED(min), float *UNUSED(max)) +{ + return; +} + +static RayObjectAPI octree_api = +{ + RE_rayobject_octree_intersect, + RE_rayobject_octree_add, + RE_rayobject_octree_done, + RE_rayobject_octree_free, + RE_rayobject_octree_bb, + RE_rayobject_octree_cost, + RE_rayobject_octree_hint_bb +}; + +/* **************** ocval method ******************* */ +/* within one octree node, a set of 3x15 bits defines a 'boundbox' to OR with */ + +#define OCVALRES 15 +#define BROW16(min, max) \ + (((max) >= OCVALRES ? 0xFFFF : (1 << ((max) + 1)) - 1) - (((min) > 0) ? ((1 << (min)) - 1) : 0)) + +static void calc_ocval_face(float *v1, float *v2, float *v3, float *v4, short x, short y, short z, OcVal *ov) +{ + float min[3], max[3]; + int ocmin, ocmax; + + copy_v3_v3(min, v1); + copy_v3_v3(max, v1); + DO_MINMAX(v2, min, max); + DO_MINMAX(v3, min, max); + if (v4) { + DO_MINMAX(v4, min, max); + } + + ocmin = OCVALRES * (min[0] - x); + ocmax = OCVALRES * (max[0] - x); + ov->ocx = BROW16(ocmin, ocmax); + + ocmin = OCVALRES * (min[1] - y); + ocmax = OCVALRES * (max[1] - y); + ov->ocy = BROW16(ocmin, ocmax); + + ocmin = OCVALRES * (min[2] - z); + ocmax = OCVALRES * (max[2] - z); + ov->ocz = BROW16(ocmin, ocmax); + +} + +static void calc_ocval_ray(OcVal *ov, float xo, float yo, float zo, float *vec1, float *vec2) +{ + int ocmin, ocmax; + + if (vec1[0] < vec2[0]) { + ocmin = OCVALRES * (vec1[0] - xo); + ocmax = OCVALRES * (vec2[0] - xo); + } + else { + ocmin = OCVALRES * (vec2[0] - xo); + ocmax = OCVALRES * (vec1[0] - xo); + } + ov->ocx = BROW16(ocmin, ocmax); + + if (vec1[1] < vec2[1]) { + ocmin = OCVALRES * (vec1[1] - yo); + ocmax = OCVALRES * (vec2[1] - yo); + } + else { + ocmin = OCVALRES * (vec2[1] - yo); + ocmax = OCVALRES * (vec1[1] - yo); + } + ov->ocy = BROW16(ocmin, ocmax); + + if (vec1[2] < vec2[2]) { + ocmin = OCVALRES * (vec1[2] - zo); + ocmax = OCVALRES * (vec2[2] - zo); + } + else { + ocmin = OCVALRES * (vec2[2] - zo); + ocmax = OCVALRES * (vec1[2] - zo); + } + ov->ocz = BROW16(ocmin, ocmax); +} + +/* ************* octree ************** */ + +static Branch *addbranch(Octree *oc, Branch *br, short ocb) +{ + int index; + + if (br->b[ocb]) return br->b[ocb]; + + oc->branchcount++; + index = oc->branchcount >> 12; + + if (oc->adrbranch[index] == NULL) + oc->adrbranch[index] = (Branch *)MEM_callocN(4096 * sizeof(Branch), "new oc branch"); + + if (oc->branchcount >= BRANCH_ARRAY * 4096) { + printf("error; octree branches full\n"); + oc->branchcount = 0; + } + + return br->b[ocb] = oc->adrbranch[index] + (oc->branchcount & 4095); +} + +static Node *addnode(Octree *oc) +{ + int index; + + oc->nodecount++; + index = oc->nodecount >> 12; + + if (oc->adrnode[index] == NULL) + oc->adrnode[index] = (Node *)MEM_callocN(4096 * sizeof(Node), "addnode"); + + if (oc->nodecount > NODE_ARRAY * NODE_ARRAY) { + printf("error; octree nodes full\n"); + oc->nodecount = 0; + } + + return oc->adrnode[index] + (oc->nodecount & 4095); +} + +static bool face_in_node(RayFace *face, short x, short y, short z, float rtf[4][3]) +{ + static float nor[3], d; + float fx, fy, fz; + + // init static vars + if (face) { + normal_tri_v3(nor, rtf[0], rtf[1], rtf[2]); + d = -nor[0] * rtf[0][0] - nor[1] * rtf[0][1] - nor[2] * rtf[0][2]; + return 0; + } + + fx = x; + fy = y; + fz = z; + + if ((fx) * nor[0] + (fy) * nor[1] + (fz) * nor[2] + d > 0.0f) { + if ((fx + 1) * nor[0] + (fy ) * nor[1] + (fz ) * nor[2] + d < 0.0f) return 1; + if ((fx ) * nor[0] + (fy + 1) * nor[1] + (fz ) * nor[2] + d < 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy + 1) * nor[1] + (fz ) * nor[2] + d < 0.0f) return 1; + + if ((fx ) * nor[0] + (fy ) * nor[1] + (fz + 1) * nor[2] + d < 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy ) * nor[1] + (fz + 1) * nor[2] + d < 0.0f) return 1; + if ((fx ) * nor[0] + (fy + 1) * nor[1] + (fz + 1) * nor[2] + d < 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy + 1) * nor[1] + (fz + 1) * nor[2] + d < 0.0f) return 1; + } + else { + if ((fx + 1) * nor[0] + (fy ) * nor[1] + (fz ) * nor[2] + d > 0.0f) return 1; + if ((fx ) * nor[0] + (fy + 1) * nor[1] + (fz ) * nor[2] + d > 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy + 1) * nor[1] + (fz ) * nor[2] + d > 0.0f) return 1; + + if ((fx ) * nor[0] + (fy ) * nor[1] + (fz + 1) * nor[2] + d > 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy ) * nor[1] + (fz + 1) * nor[2] + d > 0.0f) return 1; + if ((fx ) * nor[0] + (fy + 1) * nor[1] + (fz + 1) * nor[2] + d > 0.0f) return 1; + if ((fx + 1) * nor[0] + (fy + 1) * nor[1] + (fz + 1) * nor[2] + d > 0.0f) return 1; + } + + return 0; +} + +static void ocwrite(Octree *oc, RayFace *face, int quad, short x, short y, short z, float rtf[4][3]) +{ + Branch *br; + Node *no; + short a, oc0, oc1, oc2, oc3, oc4, oc5; + + x <<= 2; + y <<= 1; + + br = oc->adrbranch[0]; + + if (oc->ocres == 512) { + oc0 = ((x & 1024) + (y & 512) + (z & 256)) >> 8; + br = addbranch(oc, br, oc0); + } + if (oc->ocres >= 256) { + oc0 = ((x & 512) + (y & 256) + (z & 128)) >> 7; + br = addbranch(oc, br, oc0); + } + if (oc->ocres >= 128) { + oc0 = ((x & 256) + (y & 128) + (z & 64)) >> 6; + br = addbranch(oc, br, oc0); + } + + oc0 = ((x & 128) + (y & 64) + (z & 32)) >> 5; + oc1 = ((x & 64) + (y & 32) + (z & 16)) >> 4; + oc2 = ((x & 32) + (y & 16) + (z & 8)) >> 3; + oc3 = ((x & 16) + (y & 8) + (z & 4)) >> 2; + oc4 = ((x & 8) + (y & 4) + (z & 2)) >> 1; + oc5 = ((x & 4) + (y & 2) + (z & 1)); + + br = addbranch(oc, br, oc0); + br = addbranch(oc, br, oc1); + br = addbranch(oc, br, oc2); + br = addbranch(oc, br, oc3); + br = addbranch(oc, br, oc4); + no = (Node *)br->b[oc5]; + if (no == NULL) br->b[oc5] = (Branch *)(no = addnode(oc)); + + while (no->next) no = no->next; + + a = 0; + if (no->v[7]) { /* node full */ + no->next = addnode(oc); + no = no->next; + } + else { + while (no->v[a] != NULL) a++; + } + + no->v[a] = (RayFace *) RE_rayobject_align(face); + + if (quad) + calc_ocval_face(rtf[0], rtf[1], rtf[2], rtf[3], x >> 2, y >> 1, z, &no->ov[a]); + else + calc_ocval_face(rtf[0], rtf[1], rtf[2], NULL, x >> 2, y >> 1, z, &no->ov[a]); +} + +static void d2dda(Octree *oc, short b1, short b2, short c1, short c2, char *ocface, short rts[4][3], float rtf[4][3]) +{ + int ocx1, ocx2, ocy1, ocy2; + int x, y, dx = 0, dy = 0; + float ox1, ox2, oy1, oy2; + float lambda, lambda_o, lambda_x, lambda_y, ldx, ldy; + + ocx1 = rts[b1][c1]; + ocy1 = rts[b1][c2]; + ocx2 = rts[b2][c1]; + ocy2 = rts[b2][c2]; + + if (ocx1 == ocx2 && ocy1 == ocy2) { + ocface[oc->ocres * ocx1 + ocy1] = 1; + return; + } + + ox1 = rtf[b1][c1]; + oy1 = rtf[b1][c2]; + ox2 = rtf[b2][c1]; + oy2 = rtf[b2][c2]; + + if (ox1 != ox2) { + if (ox2 - ox1 > 0.0f) { + lambda_x = (ox1 - ocx1 - 1.0f) / (ox1 - ox2); + ldx = -1.0f / (ox1 - ox2); + dx = 1; + } + else { + lambda_x = (ox1 - ocx1) / (ox1 - ox2); + ldx = 1.0f / (ox1 - ox2); + dx = -1; + } + } + else { + lambda_x = 1.0f; + ldx = 0; + } + + if (oy1 != oy2) { + if (oy2 - oy1 > 0.0f) { + lambda_y = (oy1 - ocy1 - 1.0f) / (oy1 - oy2); + ldy = -1.0f / (oy1 - oy2); + dy = 1; + } + else { + lambda_y = (oy1 - ocy1) / (oy1 - oy2); + ldy = 1.0f / (oy1 - oy2); + dy = -1; + } + } + else { + lambda_y = 1.0f; + ldy = 0; + } + + x = ocx1; y = ocy1; + lambda = MIN2(lambda_x, lambda_y); + + while (true) { + + if (x < 0 || y < 0 || x >= oc->ocres || y >= oc->ocres) { + /* pass*/ + } + else { + ocface[oc->ocres * x + y] = 1; + } + + lambda_o = lambda; + if (lambda_x == lambda_y) { + lambda_x += ldx; + x += dx; + lambda_y += ldy; + y += dy; + } + else { + if (lambda_x < lambda_y) { + lambda_x += ldx; + x += dx; + } + else { + lambda_y += ldy; + y += dy; + } + } + lambda = MIN2(lambda_x, lambda_y); + if (lambda == lambda_o) break; + if (lambda >= 1.0f) break; + } + ocface[oc->ocres * ocx2 + ocy2] = 1; +} + +static void filltriangle(Octree *oc, short c1, short c2, char *ocface, short *ocmin, short *ocmax) +{ + int a, x, y, y1, y2; + + for (x = ocmin[c1]; x <= ocmax[c1]; x++) { + a = oc->ocres * x; + for (y = ocmin[c2]; y <= ocmax[c2]; y++) { + if (ocface[a + y]) { + y++; + while (ocface[a + y] && y != ocmax[c2]) y++; + for (y1 = ocmax[c2]; y1 > y; y1--) { + if (ocface[a + y1]) { + for (y2 = y; y2 <= y1; y2++) ocface[a + y2] = 1; + y1 = 0; + } + } + y = ocmax[c2]; + } + } + } +} + +static void RE_rayobject_octree_free(RayObject *tree) +{ + Octree *oc = (Octree *)tree; + +#if 0 + printf("branches %d nodes %d\n", oc->branchcount, oc->nodecount); + printf("raycount %d\n", raycount); + printf("ray coherent %d\n", coherent_ray); + printf("accepted %d rejected %d\n", accepted, rejected); +#endif + if (oc->ocface) + MEM_freeN(oc->ocface); + + if (oc->adrbranch) { + int a = 0; + while (oc->adrbranch[a]) { + MEM_freeN(oc->adrbranch[a]); + oc->adrbranch[a] = NULL; + a++; + } + MEM_freeN(oc->adrbranch); + oc->adrbranch = NULL; + } + oc->branchcount = 0; + + if (oc->adrnode) { + int a = 0; + while (oc->adrnode[a]) { + MEM_freeN(oc->adrnode[a]); + oc->adrnode[a] = NULL; + a++; + } + MEM_freeN(oc->adrnode); + oc->adrnode = NULL; + } + oc->nodecount = 0; + + MEM_freeN(oc); +} + + +RayObject *RE_rayobject_octree_create(int ocres, int size) +{ + Octree *oc = (Octree *)MEM_callocN(sizeof(Octree), "Octree"); + assert(RE_rayobject_isAligned(oc) ); /* RayObject API assumes real data to be 4-byte aligned */ + + oc->rayobj.api = &octree_api; + + oc->ocres = ocres; + + oc->ro_nodes = (RayFace **)MEM_callocN(sizeof(RayFace *) * size, "octree rayobject nodes"); + oc->ro_nodes_size = size; + oc->ro_nodes_used = 0; + + + return RE_rayobject_unalignRayAPI((RayObject *) oc); +} + + +static void RE_rayobject_octree_add(RayObject *tree, RayObject *node) +{ + Octree *oc = (Octree *)tree; + + assert(RE_rayobject_isRayFace(node) ); + assert(oc->ro_nodes_used < oc->ro_nodes_size); + oc->ro_nodes[oc->ro_nodes_used++] = (RayFace *)RE_rayobject_align(node); +} + +static void octree_fill_rayface(Octree *oc, RayFace *face) +{ + float ocfac[3], rtf[4][3]; + float co1[3], co2[3], co3[3], co4[3]; + short rts[4][3]; + short ocmin[3], ocmax[3]; + char *ocface = oc->ocface; // front, top, size view of face, to fill in + int a, b, c, oc1, oc2, oc3, oc4, x, y, z, ocres2; + + ocfac[0] = oc->ocfacx; + ocfac[1] = oc->ocfacy; + ocfac[2] = oc->ocfacz; + + ocres2 = oc->ocres * oc->ocres; + + copy_v3_v3(co1, face->v1); + copy_v3_v3(co2, face->v2); + copy_v3_v3(co3, face->v3); + if (RE_rayface_isQuad(face)) + copy_v3_v3(co4, face->v4); + + for (c = 0; c < 3; c++) { + rtf[0][c] = (co1[c] - oc->min[c]) * ocfac[c]; + rts[0][c] = (short)rtf[0][c]; + rtf[1][c] = (co2[c] - oc->min[c]) * ocfac[c]; + rts[1][c] = (short)rtf[1][c]; + rtf[2][c] = (co3[c] - oc->min[c]) * ocfac[c]; + rts[2][c] = (short)rtf[2][c]; + if (RE_rayface_isQuad(face)) { + rtf[3][c] = (co4[c] - oc->min[c]) * ocfac[c]; + rts[3][c] = (short)rtf[3][c]; + } + } + + for (c = 0; c < 3; c++) { + oc1 = rts[0][c]; + oc2 = rts[1][c]; + oc3 = rts[2][c]; + if (!RE_rayface_isQuad(face)) { + ocmin[c] = min_iii(oc1, oc2, oc3); + ocmax[c] = max_iii(oc1, oc2, oc3); + } + else { + oc4 = rts[3][c]; + ocmin[c] = min_iiii(oc1, oc2, oc3, oc4); + ocmax[c] = max_iiii(oc1, oc2, oc3, oc4); + } + if (ocmax[c] > oc->ocres - 1) ocmax[c] = oc->ocres - 1; + if (ocmin[c] < 0) ocmin[c] = 0; + } + + if (ocmin[0] == ocmax[0] && ocmin[1] == ocmax[1] && ocmin[2] == ocmax[2]) { + ocwrite(oc, face, RE_rayface_isQuad(face), ocmin[0], ocmin[1], ocmin[2], rtf); + } + else { + + d2dda(oc, 0, 1, 0, 1, ocface + ocres2, rts, rtf); + d2dda(oc, 0, 1, 0, 2, ocface, rts, rtf); + d2dda(oc, 0, 1, 1, 2, ocface + 2 * ocres2, rts, rtf); + d2dda(oc, 1, 2, 0, 1, ocface + ocres2, rts, rtf); + d2dda(oc, 1, 2, 0, 2, ocface, rts, rtf); + d2dda(oc, 1, 2, 1, 2, ocface + 2 * ocres2, rts, rtf); + if (!RE_rayface_isQuad(face)) { + d2dda(oc, 2, 0, 0, 1, ocface + ocres2, rts, rtf); + d2dda(oc, 2, 0, 0, 2, ocface, rts, rtf); + d2dda(oc, 2, 0, 1, 2, ocface + 2 * ocres2, rts, rtf); + } + else { + d2dda(oc, 2, 3, 0, 1, ocface + ocres2, rts, rtf); + d2dda(oc, 2, 3, 0, 2, ocface, rts, rtf); + d2dda(oc, 2, 3, 1, 2, ocface + 2 * ocres2, rts, rtf); + d2dda(oc, 3, 0, 0, 1, ocface + ocres2, rts, rtf); + d2dda(oc, 3, 0, 0, 2, ocface, rts, rtf); + d2dda(oc, 3, 0, 1, 2, ocface + 2 * ocres2, rts, rtf); + } + /* nothing todo with triangle..., just fills :) */ + filltriangle(oc, 0, 1, ocface + ocres2, ocmin, ocmax); + filltriangle(oc, 0, 2, ocface, ocmin, ocmax); + filltriangle(oc, 1, 2, ocface + 2 * ocres2, ocmin, ocmax); + + /* init static vars here */ + face_in_node(face, 0, 0, 0, rtf); + + for (x = ocmin[0]; x <= ocmax[0]; x++) { + a = oc->ocres * x; + for (y = ocmin[1]; y <= ocmax[1]; y++) { + if (ocface[a + y + ocres2]) { + b = oc->ocres * y + 2 * ocres2; + for (z = ocmin[2]; z <= ocmax[2]; z++) { + if (ocface[b + z] && ocface[a + z]) { + if (face_in_node(NULL, x, y, z, rtf)) + ocwrite(oc, face, RE_rayface_isQuad(face), x, y, z, rtf); + } + } + } + } + } + + /* same loops to clear octree, doubt it can be done smarter */ + for (x = ocmin[0]; x <= ocmax[0]; x++) { + a = oc->ocres * x; + for (y = ocmin[1]; y <= ocmax[1]; y++) { + /* x-y */ + ocface[a + y + ocres2] = 0; + + b = oc->ocres * y + 2 * ocres2; + for (z = ocmin[2]; z <= ocmax[2]; z++) { + /* y-z */ + ocface[b + z] = 0; + /* x-z */ + ocface[a + z] = 0; + } + } + } + } +} + +static void RE_rayobject_octree_done(RayObject *tree) +{ + Octree *oc = (Octree *)tree; + int c; + float t00, t01, t02; + int ocres2 = oc->ocres * oc->ocres; + + INIT_MINMAX(oc->min, oc->max); + + /* Calculate Bounding Box */ + for (c = 0; c < oc->ro_nodes_used; c++) + RE_rayobject_merge_bb(RE_rayobject_unalignRayFace(oc->ro_nodes[c]), oc->min, oc->max); + + /* Alloc memory */ + oc->adrbranch = (Branch **)MEM_callocN(sizeof(void *) * BRANCH_ARRAY, "octree branches"); + oc->adrnode = (Node **)MEM_callocN(sizeof(void *) * NODE_ARRAY, "octree nodes"); + + oc->adrbranch[0] = (Branch *)MEM_callocN(4096 * sizeof(Branch), "makeoctree"); + + /* the lookup table, per face, for which nodes to fill in */ + oc->ocface = (char *)MEM_callocN(3 * ocres2 + 8, "ocface"); + memset(oc->ocface, 0, 3 * ocres2); + + for (c = 0; c < 3; c++) { /* octree enlarge, still needed? */ + oc->min[c] -= 0.01f; + oc->max[c] += 0.01f; + } + + t00 = oc->max[0] - oc->min[0]; + t01 = oc->max[1] - oc->min[1]; + t02 = oc->max[2] - oc->min[2]; + + /* this minus 0.1 is old safety... seems to be needed? */ + oc->ocfacx = (oc->ocres - 0.1f) / t00; + oc->ocfacy = (oc->ocres - 0.1f) / t01; + oc->ocfacz = (oc->ocres - 0.1f) / t02; + + oc->ocsize = sqrtf(t00 * t00 + t01 * t01 + t02 * t02); /* global, max size octree */ + + for (c = 0; c < oc->ro_nodes_used; c++) { + octree_fill_rayface(oc, oc->ro_nodes[c]); + } + + MEM_freeN(oc->ocface); + oc->ocface = NULL; + MEM_freeN(oc->ro_nodes); + oc->ro_nodes = NULL; + +#if 0 + printf("%f %f - %f\n", oc->min[0], oc->max[0], oc->ocfacx); + printf("%f %f - %f\n", oc->min[1], oc->max[1], oc->ocfacy); + printf("%f %f - %f\n", oc->min[2], oc->max[2], oc->ocfacz); +#endif +} + +static void RE_rayobject_octree_bb(RayObject *tree, float *min, float *max) +{ + Octree *oc = (Octree *)tree; + DO_MINMAX(oc->min, min, max); + DO_MINMAX(oc->max, min, max); +} + +/* check all faces in this node */ +static int testnode(Octree *UNUSED(oc), Isect *is, Node *no, OcVal ocval) +{ + short nr = 0; + + /* return on any first hit */ + if (is->mode == RE_RAY_SHADOW) { + + for (; no; no = no->next) { + for (nr = 0; nr < 8; nr++) { + RayFace *face = no->v[nr]; + OcVal *ov = no->ov + nr; + + if (!face) break; + + if ( (ov->ocx & ocval.ocx) && (ov->ocy & ocval.ocy) && (ov->ocz & ocval.ocz) ) { + if (RE_rayobject_intersect(RE_rayobject_unalignRayFace(face), is) ) + return 1; + } + } + } + } + else { + /* else mirror or glass or shadowtra, return closest face */ + int found = 0; + + for (; no; no = no->next) { + for (nr = 0; nr < 8; nr++) { + RayFace *face = no->v[nr]; + OcVal *ov = no->ov + nr; + + if (!face) break; + + if ( (ov->ocx & ocval.ocx) && (ov->ocy & ocval.ocy) && (ov->ocz & ocval.ocz) ) { + if (RE_rayobject_intersect(RE_rayobject_unalignRayFace(face), is) ) { + found = 1; + } + } + } + } + + return found; + } + + return 0; +} + +/* find the Node for the octree coord x y z */ +static Node *ocread(Octree *oc, int x, int y, int z) +{ + Branch *br; + int oc1; + + x <<= 2; + y <<= 1; + + br = oc->adrbranch[0]; + + if (oc->ocres == 512) { + oc1 = ((x & 1024) + (y & 512) + (z & 256)) >> 8; + br = br->b[oc1]; + if (br == NULL) { + return NULL; + } + } + if (oc->ocres >= 256) { + oc1 = ((x & 512) + (y & 256) + (z & 128)) >> 7; + br = br->b[oc1]; + if (br == NULL) { + return NULL; + } + } + if (oc->ocres >= 128) { + oc1 = ((x & 256) + (y & 128) + (z & 64)) >> 6; + br = br->b[oc1]; + if (br == NULL) { + return NULL; + } + } + + oc1 = ((x & 128) + (y & 64) + (z & 32)) >> 5; + br = br->b[oc1]; + if (br) { + oc1 = ((x & 64) + (y & 32) + (z & 16)) >> 4; + br = br->b[oc1]; + if (br) { + oc1 = ((x & 32) + (y & 16) + (z & 8)) >> 3; + br = br->b[oc1]; + if (br) { + oc1 = ((x & 16) + (y & 8) + (z & 4)) >> 2; + br = br->b[oc1]; + if (br) { + oc1 = ((x & 8) + (y & 4) + (z & 2)) >> 1; + br = br->b[oc1]; + if (br) { + oc1 = ((x & 4) + (y & 2) + (z & 1)); + return (Node *)br->b[oc1]; + } + } + } + } + } + + return NULL; +} + +static int cliptest(float p, float q, float *u1, float *u2) +{ + float r; + + if (p < 0.0f) { + if (q < p) return 0; + else if (q < 0.0f) { + r = q / p; + if (r > *u2) return 0; + else if (r > *u1) *u1 = r; + } + } + else { + if (p > 0.0f) { + if (q < 0.0f) return 0; + else if (q < p) { + r = q / p; + if (r < *u1) return 0; + else if (r < *u2) *u2 = r; + } + } + else if (q < 0.0f) return 0; + } + return 1; +} + +/* extensive coherence checks/storage cancels out the benefit of it, and gives errors... we + * need better methods, sample code commented out below (ton) */ + +#if 0 + +in top : static int coh_nodes[16 * 16 * 16][6]; +in makeoctree : memset(coh_nodes, 0, sizeof(coh_nodes)); + +static void add_coherence_test(int ocx1, int ocx2, int ocy1, int ocy2, int ocz1, int ocz2) +{ + short *sp; + + sp = coh_nodes[(ocx2 & 15) + 16 * (ocy2 & 15) + 256 * (ocz2 & 15)]; + sp[0] = ocx1; sp[1] = ocy1; sp[2] = ocz1; + sp[3] = ocx2; sp[4] = ocy2; sp[5] = ocz2; + +} + +static int do_coherence_test(int ocx1, int ocx2, int ocy1, int ocy2, int ocz1, int ocz2) +{ + short *sp; + + sp = coh_nodes[(ocx2 & 15) + 16 * (ocy2 & 15) + 256 * (ocz2 & 15)]; + if (sp[0] == ocx1 && sp[1] == ocy1 && sp[2] == ocz1 && + sp[3] == ocx2 && sp[4] == ocy2 && sp[5] == ocz2) return 1; + return 0; +} + +#endif + +/* return 1: found valid intersection */ +/* starts with is->orig.face */ +static int RE_rayobject_octree_intersect(RayObject *tree, Isect *is) +{ + Octree *oc = (Octree *)tree; + Node *no; + OcVal ocval; + float vec1[3], vec2[3], start[3], end[3]; + float u1, u2, ox1, ox2, oy1, oy2, oz1, oz2; + float lambda_o, lambda_x, ldx, lambda_y, ldy, lambda_z, ldz, dda_lambda; + float o_lambda = 0; + int dx, dy, dz; + int xo, yo, zo, c1 = 0; + int ocx1, ocx2, ocy1, ocy2, ocz1, ocz2; + + /* clip with octree */ + if (oc->branchcount == 0) return 0; + + /* do this before intersect calls */ +#if 0 + is->facecontr = NULL; /* to check shared edge */ + is->obcontr = 0; + is->faceisect = is->isect = 0; /* shared edge, quad half flag */ + is->userdata = oc->userdata; +#endif + + copy_v3_v3(start, is->start); + madd_v3_v3v3fl(end, is->start, is->dir, is->dist); + ldx = is->dir[0] * is->dist; + o_lambda = is->dist; + u1 = 0.0f; + u2 = 1.0f; + + /* clip with octree cube */ + if (cliptest(-ldx, start[0] - oc->min[0], &u1, &u2)) { + if (cliptest(ldx, oc->max[0] - start[0], &u1, &u2)) { + ldy = is->dir[1] * is->dist; + if (cliptest(-ldy, start[1] - oc->min[1], &u1, &u2)) { + if (cliptest(ldy, oc->max[1] - start[1], &u1, &u2)) { + ldz = is->dir[2] * is->dist; + if (cliptest(-ldz, start[2] - oc->min[2], &u1, &u2)) { + if (cliptest(ldz, oc->max[2] - start[2], &u1, &u2)) { + c1 = 1; + if (u2 < 1.0f) { + end[0] = start[0] + u2 * ldx; + end[1] = start[1] + u2 * ldy; + end[2] = start[2] + u2 * ldz; + } + + if (u1 > 0.0f) { + start[0] += u1 * ldx; + start[1] += u1 * ldy; + start[2] += u1 * ldz; + } + } + } + } + } + } + } + + if (c1 == 0) return 0; + + /* reset static variables in ocread */ + //ocread(oc, oc->ocres, 0, 0); + + /* setup 3dda to traverse octree */ + ox1 = (start[0] - oc->min[0]) * oc->ocfacx; + oy1 = (start[1] - oc->min[1]) * oc->ocfacy; + oz1 = (start[2] - oc->min[2]) * oc->ocfacz; + ox2 = (end[0] - oc->min[0]) * oc->ocfacx; + oy2 = (end[1] - oc->min[1]) * oc->ocfacy; + oz2 = (end[2] - oc->min[2]) * oc->ocfacz; + + ocx1 = (int)ox1; + ocy1 = (int)oy1; + ocz1 = (int)oz1; + ocx2 = (int)ox2; + ocy2 = (int)oy2; + ocz2 = (int)oz2; + + if (ocx1 == ocx2 && ocy1 == ocy2 && ocz1 == ocz2) { + no = ocread(oc, ocx1, ocy1, ocz1); + if (no) { + /* exact intersection with node */ + vec1[0] = ox1; vec1[1] = oy1; vec1[2] = oz1; + vec2[0] = ox2; vec2[1] = oy2; vec2[2] = oz2; + calc_ocval_ray(&ocval, (float)ocx1, (float)ocy1, (float)ocz1, vec1, vec2); + if (testnode(oc, is, no, ocval) ) return 1; + } + } + else { + int found = 0; + //static int coh_ocx1, coh_ocx2, coh_ocy1, coh_ocy2, coh_ocz1, coh_ocz2; + float dox, doy, doz; + int eqval; + + /* calc lambda en ld */ + dox = ox1 - ox2; + doy = oy1 - oy2; + doz = oz1 - oz2; + + if (dox < -FLT_EPSILON) { + ldx = -1.0f / dox; + lambda_x = (ocx1 - ox1 + 1.0f) * ldx; + dx = 1; + } + else if (dox > FLT_EPSILON) { + ldx = 1.0f / dox; + lambda_x = (ox1 - ocx1) * ldx; + dx = -1; + } + else { + lambda_x = 1.0f; + ldx = 0; + dx = 0; + } + + if (doy < -FLT_EPSILON) { + ldy = -1.0f / doy; + lambda_y = (ocy1 - oy1 + 1.0f) * ldy; + dy = 1; + } + else if (doy > FLT_EPSILON) { + ldy = 1.0f / doy; + lambda_y = (oy1 - ocy1) * ldy; + dy = -1; + } + else { + lambda_y = 1.0f; + ldy = 0; + dy = 0; + } + + if (doz < -FLT_EPSILON) { + ldz = -1.0f / doz; + lambda_z = (ocz1 - oz1 + 1.0f) * ldz; + dz = 1; + } + else if (doz > FLT_EPSILON) { + ldz = 1.0f / doz; + lambda_z = (oz1 - ocz1) * ldz; + dz = -1; + } + else { + lambda_z = 1.0f; + ldz = 0; + dz = 0; + } + + xo = ocx1; yo = ocy1; zo = ocz1; + dda_lambda = min_fff(lambda_x, lambda_y, lambda_z); + + vec2[0] = ox1; + vec2[1] = oy1; + vec2[2] = oz1; + + /* this loop has been constructed to make sure the first and last node of ray + * are always included, even when dda_lambda==1.0f or larger */ + + while (true) { + + no = ocread(oc, xo, yo, zo); + if (no) { + + /* calculate ray intersection with octree node */ + copy_v3_v3(vec1, vec2); + // dox, y, z is negative + vec2[0] = ox1 - dda_lambda * dox; + vec2[1] = oy1 - dda_lambda * doy; + vec2[2] = oz1 - dda_lambda * doz; + calc_ocval_ray(&ocval, (float)xo, (float)yo, (float)zo, vec1, vec2); + + //is->dist = (u1 + dda_lambda * (u2 - u1)) * o_lambda; + if (testnode(oc, is, no, ocval) ) + found = 1; + + if (is->dist < (u1 + dda_lambda * (u2 - u1)) * o_lambda) + return found; + } + + + lambda_o = dda_lambda; + + /* traversing octree nodes need careful detection of smallest values, with proper + * exceptions for equal lambdas */ + eqval = (lambda_x == lambda_y); + if (lambda_y == lambda_z) eqval += 2; + if (lambda_x == lambda_z) eqval += 4; + + if (eqval) { // only 4 cases exist! + if (eqval == 7) { // x=y=z + xo += dx; lambda_x += ldx; + yo += dy; lambda_y += ldy; + zo += dz; lambda_z += ldz; + } + else if (eqval == 1) { // x=y + if (lambda_y < lambda_z) { + xo += dx; lambda_x += ldx; + yo += dy; lambda_y += ldy; + } + else { + zo += dz; lambda_z += ldz; + } + } + else if (eqval == 2) { // y=z + if (lambda_x < lambda_y) { + xo += dx; lambda_x += ldx; + } + else { + yo += dy; lambda_y += ldy; + zo += dz; lambda_z += ldz; + } + } + else { // x=z + if (lambda_y < lambda_x) { + yo += dy; lambda_y += ldy; + } + else { + xo += dx; lambda_x += ldx; + zo += dz; lambda_z += ldz; + } + } + } + else { // all three different, just three cases exist + eqval = (lambda_x < lambda_y); + if (lambda_y < lambda_z) eqval += 2; + if (lambda_x < lambda_z) eqval += 4; + + if (eqval == 7 || eqval == 5) { // x smallest + xo += dx; lambda_x += ldx; + } + else if (eqval == 2 || eqval == 6) { // y smallest + yo += dy; lambda_y += ldy; + } + else { // z smallest + zo += dz; lambda_z += ldz; + } + + } + + dda_lambda = min_fff(lambda_x, lambda_y, lambda_z); + if (dda_lambda == lambda_o) break; + /* to make sure the last node is always checked */ + if (lambda_o >= 1.0f) break; + } + } + + /* reached end, no intersections found */ + return 0; +} + + + |