[Bf-blender-cvs] [b6b78a9] temp-mball-refactor: Patch: T43678 BVH tree for Metaball calculation

Campbell Barton noreply at git.blender.org
Fri Feb 27 16:52:13 CET 2015


Commit: b6b78a9a80c9039e5d3009c85dab1217bb9238fe
Author: Campbell Barton
Date:   Sat Feb 28 01:07:55 2015 +1100
Branches: temp-mball-refactor
https://developer.blender.org/rBb6b78a9a80c9039e5d3009c85dab1217bb9238fe

Patch: T43678 BVH tree for Metaball calculation

T43678#291857 by @chrisr (mball_26_02_fix.diff)

===================================================================

M	source/blender/blenkernel/intern/mball_tessellate.c

===================================================================

diff --git a/source/blender/blenkernel/intern/mball_tessellate.c b/source/blender/blenkernel/intern/mball_tessellate.c
index d57f493..f8b6f73 100644
--- a/source/blender/blenkernel/intern/mball_tessellate.c
+++ b/source/blender/blenkernel/intern/mball_tessellate.c
@@ -44,6 +44,7 @@
 #include "BLI_path_util.h"
 #include "BLI_math.h"
 #include "BLI_utildefines.h"
+#include "BLI_memarena.h"
 
 #include "BKE_global.h"
 
@@ -55,16 +56,6 @@
 
 /* Data types */
 
-typedef struct vertex {         /* surface vertex */
-	float co[3];  /* position and surface normal */
-	float no[3];
-} VERTEX;
-
-typedef struct vertices {       /* list of vertices in polygonization */
-	int count, max;             /* # vertices, max # allowed */
-	VERTEX *ptr;                /* dynamically allocated */
-} VERTICES;
-
 typedef struct corner {         /* corner of a cube */
 	int i, j, k;                /* (i, j, k) is index within lattice */
 	float co[3], value;       /* location and function value */
@@ -102,70 +93,141 @@ typedef struct intlists {       /* list of list of integers */
 	struct intlists *next;      /* remaining elements */
 } INTLISTS;
 
-/* dividing scene using octal tree makes polygonisation faster */
-typedef struct ml_pointer {
-	struct ml_pointer *next, *prev;
-	struct MetaElem *ml;
-} ml_pointer;
-
-typedef struct octal_node {
-	struct octal_node *nodes[8];/* children of current node */
-	struct octal_node *parent;  /* parent of current node */
-	struct ListBase elems;      /* ListBase of MetaElem pointers (ml_pointer) */
-	float x_min, y_min, z_min;  /* 1st border point */
-	float x_max, y_max, z_max;  /* 7th border point */
-	float x, y, z;              /* center of node */
-	int pos, neg;               /* number of positive and negative MetaElements in the node */
-	int count;                  /* number of MetaElems, which belongs to the node */
-} octal_node;
-
-typedef struct octal_tree {
-	struct octal_node *first;   /* first node */
-	int pos, neg;               /* number of positive and negative MetaElements in the scene */
-	short depth;                /* number of scene subdivision */
-} octal_tree;
-
-struct pgn_elements {
-	struct pgn_elements *next, *prev;
-	char *data;
-};
+typedef struct Box{				/* an AABB with pointer to metalelem */
+	float min[3], max[3];
+	MetaElem *ml;
+} Box;
+
+typedef struct MetaballBVHNode{	/* BVH node */
+	Box bb[2];					/* AABB of children */
+	struct MetaballBVHNode *child[2];
+} MetaballBVHNode;
+
+typedef struct process {        /* parameters, storage */
+	float thresh, size;			/* mball threshold, single cube size */
+	float delta;				/* small delta for calculating normals */
+	int converge_res;			/* converge procedure resolution (more = slower) */
+
+	MetaElem **mainb;			/* array of all metaelems */
+	int totelem, mem;			/* number of metaelems */
 
-typedef struct process {        /* parameters, function, storage */
-	/* ** old G_mb contents ** */
-	float thresh;
-	int totelem;
-	MetaElem **mainb;
-	octal_tree *metaball_tree;
-
-	/* ** old process contents ** */
-
-	/* what happens here? floats, I think. */
-	/*  float (*function)(void);	 */	/* implicit surface function */
-	float (*function)(struct process *, float, float, float);
-	float size, delta;          /* cube size, normal delta */
-	int bounds;                 /* cube range within lattice */
-	CUBES *cubes;               /* active cubes */
-	VERTICES vertices;          /* surface vertices */
+	MetaballBVHNode metaball_bvh; /* The simplest bvh */
+	Box allbb;                   /* Bounding box of all metaelems */
+
+	MetaballBVHNode **bvh_queue; /* Queue used during bvh traversal */
+	int bvh_queue_size;
+
+	CUBES *cubes;               /* stack of cubes waiting for polygonization */
 	CENTERLIST **centers;       /* cube center hash table */
 	CORNER **corners;           /* corner value hash table */
 	EDGELIST **edges;           /* edge and vertex id hash table */
 
-	/* Runtime things */
-	int *indices;
-	int totindex, curindex;
+	int *indices;				/* output indices */
+	int totindex;				/* size of memory allocated for indices */
+	int curindex;				/* number of currently added indices */
 
-	int pgn_offset;
-	struct pgn_elements *pgn_current;
-	ListBase pgn_list;
+	float *co, *no;				/* surface vertices - positions and normals */
+	int totvertex, curvertex;	/* memory size, currently added vertices */
+
+	/* memory allocation from common pool */
+	MemArena *pgn_elements;
 } PROCESS;
 
 /* Forward declarations */
-static int vertid(PROCESS *process, const CORNER *c1, const CORNER *c2, MetaBall *mb);
-static int setcenter(PROCESS *process, CENTERLIST *table[], const int i, const int j, const int k);
-static CORNER *setcorner(PROCESS *process, int i, int j, int k);
-static void converge(PROCESS *process, const float p1[3], const float p2[3], float v1, float v2,
-                     float p[3], MetaBall *mb, int f);
+static int vertid(PROCESS *process, const CORNER *c1, const CORNER *c2);
+static void add_cube(PROCESS *process, int i, int j, int k);
+static void make_face(PROCESS *process, int i1, int i2, int i3, int i4);
+static void converge(PROCESS *process, CORNER c1, CORNER c2, float r_p[3]);
+
+/* ******************* SIMPLE BVH ********************* */
 
+static void make_union(BoundBox *a, Box *b, Box *out)
+{
+	out->min[0] = min_ff(a->vec[0][0], b->min[0]);
+	out->min[1] = min_ff(a->vec[0][1], b->min[1]);
+	out->min[2] = min_ff(a->vec[0][2], b->min[2]);
+
+	out->max[0] = max_ff(a->vec[6][0], b->max[0]);
+	out->max[1] = max_ff(a->vec[6][1], b->max[1]);
+	out->max[2] = max_ff(a->vec[6][2], b->max[2]);
+}
+
+static void make_box_from_ml(Box *r, MetaElem *ml)
+{
+	copy_v3_v3(r->max, ml->bb->vec[6]);
+	copy_v3_v3(r->min, ml->bb->vec[0]);
+	r->ml = ml;
+}
+
+/* Partitions part of mainb array [start, end) along axis s. Returns i,
+* where centroids of elements in the [start, i) segment lie "on the right side" of div,
+* and elements in the [i, end) segment lie "on the left" */
+static int partition_mainb(PROCESS *process, int start, int end, int s, float div)
+{
+	int i = start, j = end - 1;
+	div *= 2.0f;
+
+	while (1) {
+		while (i < j && process->mainb[i]->bb->vec[6][s] + process->mainb[i]->bb->vec[0][s] < div) i++;
+		while (j > i && div < process->mainb[j]->bb->vec[6][s] + process->mainb[j]->bb->vec[0][s]) j--;
+
+		if (i >= j)	break;
+
+		SWAP(MetaElem*, process->mainb[i], process->mainb[j]);
+		i++; j--;
+	}
+
+	if (i == start) i++;
+
+	return i;
+}
+
+/* Recursively builds a BVH, dividing elements along the middle of the longest axis of allbox. */
+static void build_bvh_spatial(PROCESS *process, MetaballBVHNode *node, int start, int end, Box allbox)
+{
+	int part, j, s;
+	float dim[3], div;
+
+	/* Maximum bvh queue size is number of nodes which are made, equals calls to this function. */
+	process->bvh_queue_size++;
+
+	dim[0] = allbox.max[0] - allbox.min[0];
+	dim[1] = allbox.max[1] - allbox.min[1];
+	dim[2] = allbox.max[2] - allbox.min[2];
+
+	s = 0;
+	if (dim[1] > dim[0] && dim[1] > dim[2]) s = 1;
+	else if (dim[2] > dim[1] && dim[2] > dim[0]) s = 2;
+
+	div = allbox.min[s] + (dim[s] / 2.0f);
+
+	part = partition_mainb(process, start, end, s, div);
+
+	make_box_from_ml(&node->bb[0], process->mainb[start]);
+	node->child[0] = NULL;
+
+	if (part > start + 1) {
+		for (j = start; j < part; j++)
+			make_union(process->mainb[j]->bb, &node->bb[0], &node->bb[0]);
+
+		node->child[0] = BLI_memarena_alloc(process->pgn_elements, sizeof(MetaballBVHNode));
+		build_bvh_spatial(process, node->child[0], start, part, node->bb[0]);
+	}
+
+	node->child[1] = NULL;
+	if (part < end) {
+		make_box_from_ml(&node->bb[1], process->mainb[part]);
+
+		if (part < end - 1) {
+			for (j = part; j < end; j++)
+				make_union(process->mainb[j]->bb, &node->bb[1], &node->bb[1]);
+
+			node->child[1] = BLI_memarena_alloc(process->pgn_elements, sizeof(MetaballBVHNode));
+			build_bvh_spatial(process, node->child[1], part, end, node->bb[1]);
+		}
+	}
+	else INIT_MINMAX(node->bb[1].min, node->bb[1].max);
+}
 
 /* ******************** ARITH ************************* */
 
@@ -180,8 +242,6 @@ static void converge(PROCESS *process, const float p1[3], const float p2[3], flo
  * Permission is granted to reproduce, use and distribute this code for
  * any and all purposes, provided that this notice appears in all copies. */
 
-#define RES 12 /* # converge iterations    */
-
 #define L   0  /* left direction:	-x, -i */
 #define R   1  /* right direction:	+x, +i */
 #define B   2  /* bottom direction: -y, -j */
@@ -208,16 +268,14 @@ static void converge(PROCESS *process, const float p1[3], const float p2[3], flo
 #define MB_BIT(i, bit) (((i) >> (bit)) & 1)
 #define FLIP(i, bit) ((i) ^ 1 << (bit)) /* flip the given bit of i */
 
+/* ******************** DENSITY COPMPUTATION ********************* */
 
-/* **************** POLYGONIZATION ************************ */
-
-static void calc_mballco(MetaElem *ml, float vec[3])
-{
-	if (ml->mat) {
-		mul_m4_v3((float (*)[4])ml->mat, vec);
-	}
-}
-
+/* Computes density from given metaball at given position.
+* Metaball equation is: (1 - r^2 / R^2)^3 * s
+* r = distance from center
+* R = metaball radius
+* s - metaball stiffness
+**/
 static float densfunc(MetaElem *ball, float x, float y, float z)
 {
 	float dist2;
@@ -229,37 +287,26 @@ static float densfunc(MetaElem *ball, float x, float y, float z)
 		case MB_BALL:
 			/* do nothing */
 			break;
-		case MB_TUBE:
-			if      (dvec[0] >  ball->expx) dvec[0] -= ball->expx;
-			else if (dvec[0] < -ball->expx) dvec[0] += ball->expx;
-			else                            dvec[0] = 0.0;
-			break;
+		case MB_CUBE:
+			if      (dvec[2] > ball->expz)  dvec[2] -= ball->expz;
+			else if (dvec[2] < -ball->expz) dvec[2] += ball->expz;
+			else                            dvec[2] = 0.0;
+			/* fall through */
 		case MB_PLANE:
-			if      (dvec[0] >  ball->expx) dvec[0] -= ball->expx;
-			else if (dvec[0] < -ball->expx) dvec[0] += ball->expx;
-			else                            dvec[0] = 0.0;
 			if      (dvec[1] >  ball->expy) dvec[1] -= ball->expy;
 			else if (dvec[1] < -ball->expy) dvec[1] += ball->expy;
 			else                            dvec[1] = 0.0;
+			/* fall through */
+		case MB_TUBE:
+			if

@@ Diff output truncated at 10240 characters. @@




More information about the Bf-blender-cvs mailing list