--- /dev/null
- int i, xn, yn, zn, ax, ay;
+/*
+ * ***** 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) 2007 Blender Foundation.
+ * All rights reserved.
+ *
+ * The Original Code is: all of this file.
+ *
+ * Contributor(s): Geoffrey Bantle.
+ *
+ * ***** END GPL LICENSE BLOCK *****
+ */
+
+/** \file blender/bmesh/intern/bmesh_interp.c
+ * \ingroup bmesh
+ *
+ * Functions for interpolating data across the surface of a mesh.
+ */
+
+#include "MEM_guardedalloc.h"
+
+#include "DNA_mesh_types.h"
+#include "DNA_meshdata_types.h"
+
+#include "BKE_customdata.h"
+#include "BKE_utildefines.h"
+#include "BKE_multires.h"
+
+#include "BLI_array.h"
+#include "BLI_math.h"
+#include "BLI_cellalloc.h"
+
+#include "bmesh.h"
+#include "bmesh_private.h"
+
+/**
+ * bmesh_data_interp_from_verts
+ *
+ * Interpolates per-vertex data from two sources to a target.
+ *
+ * Returns -
+ * Nothing
+ */
+void BM_Data_Interp_From_Verts(BMesh *bm, BMVert *v1, BMVert *v2, BMVert *v, float fac)
+{
+ void *src[2];
+ float w[2];
+ if (v1->head.data && v2->head.data) {
+ src[0]= v1->head.data;
+ src[1]= v2->head.data;
+ w[0] = 1.0f-fac;
+ w[1] = fac;
+ CustomData_bmesh_interp(&bm->vdata, src, w, NULL, 2, v->head.data);
+ }
+}
+
+/*
+ BM Data Vert Average
+
+ Sets all the customdata (e.g. vert, loop) associated with a vert
+ to the average of the face regions surrounding it.
+*/
+
+
+static void UNUSED_FUNCTION(BM_Data_Vert_Average)(BMesh *UNUSED(bm), BMFace *UNUSED(f))
+{
+ // BMIter iter;
+}
+
+/**
+ * bmesh_data_facevert_edgeinterp
+ *
+ * Walks around the faces of an edge and interpolates the per-face-edge
+ * data between two sources to a target.
+ *
+ * Returns -
+ * Nothing
+*/
+
+void BM_Data_Facevert_Edgeinterp(BMesh *bm, BMVert *v1, BMVert *UNUSED(v2), BMVert *v, BMEdge *e1, float fac)
+{
+ void *src[2];
+ float w[2];
+ BMLoop *l=NULL, *v1loop = NULL, *vloop = NULL, *v2loop = NULL;
+
+ w[1] = 1.0f - fac;
+ w[0] = fac;
+
+ if(!e1->l) return;
+ l = e1->l;
+ do{
+ if(l->v == v1){
+ v1loop = l;
+ vloop = v1loop->next;
+ v2loop = vloop->next;
+ }else if(l->v == v){
+ v1loop = l->next;
+ vloop = l;
+ v2loop = l->prev;
+ }
+
+ if (!v1loop || !v2loop)
+ return;
+
+ src[0] = v1loop->head.data;
+ src[1] = v2loop->head.data;
+
+ CustomData_bmesh_interp(&bm->ldata, src,w, NULL, 2, vloop->head.data);
+ l = l->radial_next;
+ }while(l!=e1->l);
+}
+
+void BM_loops_to_corners(BMesh *bm, Mesh *me, int findex,
+ BMFace *f, int numTex, int numCol)
+{
+ BMLoop *l;
+ BMIter iter;
+ MTFace *texface;
+ MTexPoly *texpoly;
+ MCol *mcol;
+ MLoopCol *mloopcol;
+ MLoopUV *mloopuv;
+ int i, j;
+
+ for(i=0; i < numTex; i++){
+ texface = CustomData_get_n(&me->fdata, CD_MTFACE, findex, i);
+ texpoly = CustomData_bmesh_get_n(&bm->pdata, f->head.data, CD_MTEXPOLY, i);
+
+ texface->tpage = texpoly->tpage;
+ texface->flag = texpoly->flag;
+ texface->transp = texpoly->transp;
+ texface->mode = texpoly->mode;
+ texface->tile = texpoly->tile;
+ texface->unwrap = texpoly->unwrap;
+
+ j = 0;
+ BM_ITER(l, &iter, bm, BM_LOOPS_OF_FACE, f) {
+ mloopuv = CustomData_bmesh_get_n(&bm->ldata, l->head.data, CD_MLOOPUV, i);
+ texface->uv[j][0] = mloopuv->uv[0];
+ texface->uv[j][1] = mloopuv->uv[1];
+
+ j++;
+ }
+
+ }
+
+ for(i=0; i < numCol; i++){
+ mcol = CustomData_get_n(&me->fdata, CD_MCOL, findex, i);
+
+ j = 0;
+ BM_ITER(l, &iter, bm, BM_LOOPS_OF_FACE, f) {
+ mloopcol = CustomData_bmesh_get_n(&bm->ldata, l->head.data, CD_MLOOPCOL, i);
+ mcol[j].r = mloopcol->r;
+ mcol[j].g = mloopcol->g;
+ mcol[j].b = mloopcol->b;
+ mcol[j].a = mloopcol->a;
+
+ j++;
+ }
+ }
+}
+
+/**
+ * BM_data_interp_from_face
+ *
+ * projects target onto source, and pulls interpolated customdata from
+ * source.
+ *
+ * Returns -
+ * Nothing
+*/
+void BM_face_interp_from_face(BMesh *bm, BMFace *target, BMFace *source)
+{
+ BMLoop *l1, *l2;
+ void **blocks=NULL;
+ float (*cos)[3]=NULL, *w=NULL;
+ BLI_array_staticdeclare(cos, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(w, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(blocks, BM_NGON_STACK_SIZE);
+
+ BM_Copy_Attributes(bm, bm, source, target);
+
+ l2 = bm_firstfaceloop(source);
+ do {
+ BLI_array_growone(cos);
+ copy_v3_v3(cos[BLI_array_count(cos)-1], l2->v->co);
+ BLI_array_growone(w);
+ BLI_array_append(blocks, l2->head.data);
+ l2 = l2->next;
+ } while (l2 != bm_firstfaceloop(source));
+
+ l1 = bm_firstfaceloop(target);
+ do {
+ interp_weights_poly_v3(w, cos, source->len, l1->v->co);
+ CustomData_bmesh_interp(&bm->ldata, blocks, w, NULL, BLI_array_count(blocks), l1->head.data);
+ l1 = l1->next;
+ } while (l1 != bm_firstfaceloop(target));
+
+ BLI_array_free(cos);
+ BLI_array_free(w);
+ BLI_array_free(blocks);
+}
+
+/****some math stuff for dealing with doubles, put here to
+ avoid merge errors - joeedh ****/
+
+#define VECMUL(a, b) (((a)[0] = (a)[0] * (b)), ((a)[1] = (a)[1] * (b)), ((a)[2] = (a)[2] * (b)))
+#define VECADD2(a, b) (((a)[0] = (a)[0] + (b)[0]), ((a)[1] = (a)[1] + (b)[1]), ((a)[2] = (a)[2] + (b)[2]))
+#define VECSUB2(a, b) (((a)[0] = (a)[0] - (b)[0]), ((a)[1] = (a)[1] - (b)[1]), ((a)[2] = (a)[2] - (b)[2]))
+
+/* find closest point to p on line through l1,l2 and return lambda,
+ * where (0 <= lambda <= 1) when cp is in the line segement l1,l2
+ */
+static double closest_to_line_v3_d(double cp[3], const double p[3], const double l1[3], const double l2[3])
+{
+ double h[3],u[3],lambda;
+ VECSUB(u, l2, l1);
+ VECSUB(h, p, l1);
+ lambda =INPR(u,h)/INPR(u,u);
+ cp[0] = l1[0] + u[0] * lambda;
+ cp[1] = l1[1] + u[1] * lambda;
+ cp[2] = l1[2] + u[2] * lambda;
+ return lambda;
+}
+
+/* point closest to v1 on line v2-v3 in 3D */
+static void UNUSED_FUNCTION(closest_to_line_segment_v3_d)(double *closest, double v1[3], double v2[3], double v3[3])
+{
+ double lambda, cp[3];
+
+ lambda= closest_to_line_v3_d(cp,v1, v2, v3);
+
+ if(lambda <= 0.0) {
+ VECCOPY(closest, v2);
+ } else if(lambda >= 1.0) {
+ VECCOPY(closest, v3);
+ } else {
+ VECCOPY(closest, cp);
+ }
+}
+
+static double UNUSED_FUNCTION(len_v3_d)(const double a[3])
+{
+ return sqrt(INPR(a, a));
+}
+
+static double UNUSED_FUNCTION(len_v3v3_d)(const double a[3], const double b[3])
+{
+ double d[3];
+
+ VECSUB(d, b, a);
+ return sqrt(INPR(d, d));
+}
+
+static void cent_quad_v3_d(double *cent, double *v1, double *v2, double *v3, double *v4)
+{
+ cent[0]= 0.25*(v1[0]+v2[0]+v3[0]+v4[0]);
+ cent[1]= 0.25*(v1[1]+v2[1]+v3[1]+v4[1]);
+ cent[2]= 0.25*(v1[2]+v2[2]+v3[2]+v4[2]);
+}
+
+static void UNUSED_FUNCTION(cent_tri_v3_d)(double *cent, double *v1, double *v2, double *v3)
+{
+ cent[0]= 0.33333*(v1[0]+v2[0]+v3[0]);
+ cent[1]= 0.33333*(v1[1]+v2[1]+v3[1]);
+ cent[2]= 0.33333*(v1[2]+v2[2]+v3[2]);
+}
+
+static void UNUSED_FUNCTION(cross_v3_v3v3_d)(double r[3], const double a[3], const double b[3])
+{
+ r[0]= a[1]*b[2] - a[2]*b[1];
+ r[1]= a[2]*b[0] - a[0]*b[2];
+ r[2]= a[0]*b[1] - a[1]*b[0];
+}
+
+/* distance v1 to line-piece v2-v3 */
+static double UNUSED_FUNCTION(dist_to_line_segment_v2_d)(double v1[3], double v2[3], double v3[3])
+{
+ double labda, rc[2], pt[2], len;
+
+ rc[0]= v3[0]-v2[0];
+ rc[1]= v3[1]-v2[1];
+ len= rc[0]*rc[0]+ rc[1]*rc[1];
+ if(len==0.0) {
+ rc[0]= v1[0]-v2[0];
+ rc[1]= v1[1]-v2[1];
+ return sqrt(rc[0]*rc[0]+ rc[1]*rc[1]);
+ }
+
+ labda= (rc[0]*(v1[0]-v2[0]) + rc[1]*(v1[1]-v2[1]))/len;
+ if(labda<=0.0) {
+ pt[0]= v2[0];
+ pt[1]= v2[1];
+ }
+ else if(labda>=1.0) {
+ pt[0]= v3[0];
+ pt[1]= v3[1];
+ }
+ else {
+ pt[0]= labda*rc[0]+v2[0];
+ pt[1]= labda*rc[1]+v2[1];
+ }
+
+ rc[0]= pt[0]-v1[0];
+ rc[1]= pt[1]-v1[1];
+ return sqrt(rc[0]*rc[0]+ rc[1]*rc[1]);
+}
+
+
+MINLINE double line_point_side_v2_d(const double *l1, const double *l2, const double *pt)
+{
+ return ((l1[0]-pt[0]) * (l2[1]-pt[1])) -
+ ((l2[0]-pt[0]) * (l1[1]-pt[1]));
+}
+
+/* point in quad - only convex quads */
+static int isect_point_quad_v2_d(double pt[2], double v1[2], double v2[2], double v3[2], double v4[2])
+{
+ if (line_point_side_v2_d(v1,v2,pt)>=0.0) {
+ if (line_point_side_v2_d(v2,v3,pt)>=0.0) {
+ if (line_point_side_v2_d(v3,v4,pt)>=0.0) {
+ if (line_point_side_v2_d(v4,v1,pt)>=0.0) {
+ return 1;
+ }
+ }
+ }
+ } else {
+ if (! (line_point_side_v2_d(v2,v3,pt)>=0.0)) {
+ if (! (line_point_side_v2_d(v3,v4,pt)>=0.0)) {
+ if (! (line_point_side_v2_d(v4,v1,pt)>=0.0)) {
+ return -1;
+ }
+ }
+ }
+ }
+
+ return 0;
+}
+
+/***** multires interpolation*****
+
+mdisps is a grid of displacements, ordered thus:
+
+ v1/center----v4/next -> x
+ | |
+ | |
+ v2/prev------v3/cur
+ |
+ V
+ y
+*/
+
+static int compute_mdisp_quad(BMLoop *l, double v1[3], double v2[3], double v3[3], double v4[3], double e1[3], double e2[3])
+{
+ double cent[3] = {0.0, 0.0, 0.0}, n[3], p[3];
+ BMLoop *l2;
+
+ /*computer center*/
+ l2 = bm_firstfaceloop(l->f);
+ do {
+ cent[0] += (double)l2->v->co[0];
+ cent[1] += (double)l2->v->co[1];
+ cent[2] += (double)l2->v->co[2];
+ l2 = l2->next;
+ } while (l2 != bm_firstfaceloop(l->f));
+
+ VECMUL(cent, (1.0/(double)l->f->len));
+
+ VECADD(p, l->prev->v->co, l->v->co);
+ VECMUL(p, 0.5);
+ VECADD(n, l->next->v->co, l->v->co);
+ VECMUL(n, 0.5);
+
+ VECCOPY(v1, cent);
+ VECCOPY(v2, p);
+ VECCOPY(v3, l->v->co);
+ VECCOPY(v4, n);
+
+ VECSUB(e1, v2, v1);
+ VECSUB(e2, v3, v4);
+
+ return 1;
+}
+
+/*funnily enough, I think this is identical to face_to_crn_interp, heh*/
+static double quad_coord(double aa[3], double bb[3], double cc[3], double dd[3], int a1, int a2)
+{
+ double x, y, z, f1;
+
+ x = aa[a1]*cc[a2]-cc[a1]*aa[a2];
+ y = aa[a1]*dd[a2]+bb[a1]*cc[a2]-cc[a1]*bb[a2]-dd[a1]*aa[a2];
+ z = bb[a1]*dd[a2]-dd[a1]*bb[a2];
+
+ if (fabs(2*(x-y+z)) > DBL_EPSILON*10.0) {
+ double f2;
+
+ f1 = (sqrt(y*y-4.0*x*z) - y + 2.0*z) / (2.0*(x-y+z));
+ f2 = (-sqrt(y*y-4.0*x*z) - y + 2.0*z) / (2.0*(x-y+z));
+
+ f1= fabs(f1);
+ f2= fabs(f2);
+ f1 = MIN2(f1, f2);
+ CLAMP(f1, 0.0, 1.0+DBL_EPSILON);
+ } else {
+ f1 = -z/(y - 2*z);
+ CLAMP(f1, 0.0, 1.0+DBL_EPSILON);
+
+ if (isnan(f1) || f1 > 1.0 || f1 < 0.0) {
+ int i;
+
+ for (i=0; i<2; i++) {
+ if (fabsf(aa[i]) < FLT_EPSILON*100)
+ return aa[(i+1)%2] / fabs(bb[(i+1)%2] - aa[(i+1)%2]);
+ if (fabsf(cc[i]) < FLT_EPSILON*100)
+ return cc[(i+1)%2] / fabs(dd[(i+1)%2] - cc[(i+1)%2]);
+ }
+ }
+ }
+
+ return f1;
+}
+
+static int quad_co(double *x, double *y, double v1[3], double v2[3], double v3[3], double v4[3], double p[3], float n[3])
+{
+ float projverts[5][3], n2[3];
+ double dprojverts[4][3], origin[3]={0.0f, 0.0f, 0.0f};
+ int i;
+
+ /*project points into 2d along normal*/
+ VECCOPY(projverts[0], v1);
+ VECCOPY(projverts[1], v2);
+ VECCOPY(projverts[2], v3);
+ VECCOPY(projverts[3], v4);
+ VECCOPY(projverts[4], p);
+
+ normal_quad_v3(n2, projverts[0], projverts[1], projverts[2], projverts[3]);
+
+ if (INPR(n, n2) < -FLT_EPSILON)
+ return 0;
+
+ /*rotate*/
+ poly_rotate_plane(n, projverts, 5);
+
+ /*flatten*/
+ for (i=0; i<5; i++) projverts[i][2] = 0.0f;
+
+ /*subtract origin*/
+ for (i=0; i<4; i++) {
+ VECSUB2(projverts[i], projverts[4]);
+ }
+
+ VECCOPY(dprojverts[0], projverts[0]);
+ VECCOPY(dprojverts[1], projverts[1]);
+ VECCOPY(dprojverts[2], projverts[2]);
+ VECCOPY(dprojverts[3], projverts[3]);
+
+ if (!isect_point_quad_v2_d(origin, dprojverts[0], dprojverts[1], dprojverts[2], dprojverts[3]))
+ return 0;
+
+ *y = quad_coord(dprojverts[1], dprojverts[0], dprojverts[2], dprojverts[3], 0, 1);
+ *x = quad_coord(dprojverts[2], dprojverts[1], dprojverts[3], dprojverts[0], 0, 1);
+
+ return 1;
+}
+
+
+/*tl is loop to project onto, l is loop whose internal displacement, co, is being
+ projected. x and y are location in loop's mdisps grid of point co.*/
+static int mdisp_in_mdispquad(BMesh *bm, BMLoop *l, BMLoop *tl, double p[3], double *x, double *y, int res)
+{
+ double v1[3], v2[3], c[3], v3[3], v4[3], e1[3], e2[3];
+ double eps = FLT_EPSILON*4000;
+
+ if (len_v3(l->v->no) == 0.0f)
+ BM_Vert_UpdateAllNormals(bm, l->v);
+ if (len_v3(tl->v->no) == 0.0f)
+ BM_Vert_UpdateAllNormals(bm, tl->v);
+
+ compute_mdisp_quad(tl, v1, v2, v3, v4, e1, e2);
+
+ /*expand quad a bit*/
+ cent_quad_v3_d(c, v1, v2, v3, v4);
+
+ VECSUB2(v1, c); VECSUB2(v2, c);
+ VECSUB2(v3, c); VECSUB2(v4, c);
+ VECMUL(v1, 1.0+eps); VECMUL(v2, 1.0+eps);
+ VECMUL(v3, 1.0+eps); VECMUL(v4, 1.0+eps);
+ VECADD2(v1, c); VECADD2(v2, c);
+ VECADD2(v3, c); VECADD2(v4, c);
+
+ if (!quad_co(x, y, v1, v2, v3, v4, p, l->v->no))
+ return 0;
+
+ *x *= res-1;
+ *y *= res-1;
+
+ return 1;
+}
+
+static void bmesh_loop_interp_mdisps(BMesh *bm, BMLoop *target, BMFace *source)
+{
+ MDisps *mdisps;
+ BMLoop *l2;
+ double x, y, d, v1[3], v2[3], v3[3], v4[3] = {0.0f, 0.0f, 0.0f}, e1[3], e2[3];
+ int ix, iy, res;
+
+ /*ignore 2-edged faces*/
+ if (target->f->len < 3)
+ return;
+
+ if (!CustomData_has_layer(&bm->ldata, CD_MDISPS))
+ return;
+
+ mdisps = CustomData_bmesh_get(&bm->ldata, target->head.data, CD_MDISPS);
+ compute_mdisp_quad(target, v1, v2, v3, v4, e1, e2);
+
+ /*if no disps data allocate a new grid, the size of the first grid in source.*/
+ if (!mdisps->totdisp) {
+ MDisps *md2 = CustomData_bmesh_get(&bm->ldata, bm_firstfaceloop(source)->head.data, CD_MDISPS);
+
+ mdisps->totdisp = md2->totdisp;
+ if (mdisps->totdisp)
+ mdisps->disps = BLI_cellalloc_calloc(sizeof(float)*3*mdisps->totdisp, "mdisp->disps in bmesh_loop_intern_mdisps");
+ else
+ return;
+ }
+
+ res = (int)sqrt(mdisps->totdisp);
+ d = 1.0/(double)(res-1);
+ for (x=0.0f, ix=0; ix<res; x += d, ix++) {
+ for (y=0.0f, iy=0; iy<res; y+= d, iy++) {
+ double co1[3], co2[3], co[3];
+ /* double xx, yy; */ /* UNUSED */
+
+ VECCOPY(co1, e1);
+
+ /* if (!iy) yy = y + (double)FLT_EPSILON*20; */
+ /* else yy = y - (double)FLT_EPSILON*20; */
+
+ VECMUL(co1, y);
+ VECADD2(co1, v1);
+
+ VECCOPY(co2, e2);
+ VECMUL(co2, y);
+ VECADD2(co2, v4);
+
+ /* if (!ix) xx = x + (double)FLT_EPSILON*20; */ /* UNUSED */
+ /* else xx = x - (double)FLT_EPSILON*20; */ /* UNUSED */
+
+ VECSUB(co, co2, co1);
+ VECMUL(co, x);
+ VECADD2(co, co1);
+
+ l2 = bm_firstfaceloop(source);
+ do {
+ double x2, y2;
+ MDisps *md1, *md2;
+
+ md1 = CustomData_bmesh_get(&bm->ldata, target->head.data, CD_MDISPS);
+ md2 = CustomData_bmesh_get(&bm->ldata, l2->head.data, CD_MDISPS);
+
+ if (mdisp_in_mdispquad(bm, target, l2, co, &x2, &y2, res)) {
+ /* int ix2 = (int)x2; */ /* UNUSED */
+ /* int iy2 = (int)y2; */ /* UNUSED */
+
+ old_mdisps_bilinear(md1->disps[iy*res+ix], md2->disps, res, (float)x2, (float)y2);
+ }
+ l2 = l2->next;
+ } while (l2 != bm_firstfaceloop(source));
+ }
+ }
+}
+
+void BM_multires_smooth_bounds(BMesh *bm, BMFace *f)
+{
+ BMLoop *l;
+ BMIter liter;
+
+ if (!CustomData_has_layer(&bm->ldata, CD_MDISPS))
+ return;
+
+ BM_ITER(l, &liter, bm, BM_LOOPS_OF_FACE, f) {
+ MDisps *mdp = CustomData_bmesh_get(&bm->ldata, l->prev->head.data, CD_MDISPS);
+ MDisps *mdl = CustomData_bmesh_get(&bm->ldata, l->head.data, CD_MDISPS);
+ MDisps *mdn = CustomData_bmesh_get(&bm->ldata, l->next->head.data, CD_MDISPS);
+ float co1[3];
+ int sides;
+ int y;
+
+ /*****
+ mdisps is a grid of displacements, ordered thus:
+
+ v4/next
+ |
+ | v1/cent-----mid2 ---> x
+ | | |
+ | | |
+ v2/prev---mid1-----v3/cur
+ |
+ V
+ y
+ *****/
+
+ sides = (int)sqrt(mdp->totdisp);
+ for (y=0; y<sides; y++) {
+ add_v3_v3v3(co1, mdn->disps[y*sides], mdl->disps[y]);
+ mul_v3_fl(co1, 0.5);
+
+ copy_v3_v3(mdn->disps[y*sides], co1);
+ copy_v3_v3(mdl->disps[y], co1);
+ }
+ }
+
+ BM_ITER(l, &liter, bm, BM_LOOPS_OF_FACE, f) {
+ MDisps *mdl1 = CustomData_bmesh_get(&bm->ldata, l->head.data, CD_MDISPS);
+ MDisps *mdl2;
+ float co1[3], co2[3], co[3];
+ int sides;
+ int y;
+
+ /*****
+ mdisps is a grid of displacements, ordered thus:
+
+ v4/next
+ |
+ | v1/cent-----mid2 ---> x
+ | | |
+ | | |
+ v2/prev---mid1-----v3/cur
+ |
+ V
+ y
+ *****/
+
+ if (l->radial_next == l)
+ continue;
+
+ if (l->radial_next->v == l->v)
+ mdl2 = CustomData_bmesh_get(&bm->ldata, l->radial_next->head.data, CD_MDISPS);
+ else
+ mdl2 = CustomData_bmesh_get(&bm->ldata, l->radial_next->next->head.data, CD_MDISPS);
+
+ sides = (int)sqrt(mdl1->totdisp);
+ for (y=0; y<sides; y++) {
+ int a1, a2, o1, o2;
+
+ if (l->v != l->radial_next->v) {
+ a1 = sides*y + sides-2;
+ a2 = (sides-2)*sides + y;
+
+ o1 = sides*y + sides-1;
+ o2 = (sides-1)*sides + y;
+ } else {
+ a1 = sides*y + sides-2;
+ a2 = sides*y + sides-2;
+ o1 = sides*y + sides-1;
+ o2 = sides*y + sides-1;
+ }
+
+ /*magic blending numbers, hardcoded!*/
+ add_v3_v3v3(co1, mdl1->disps[a1], mdl2->disps[a2]);
+ mul_v3_fl(co1, 0.18);
+
+ add_v3_v3v3(co2, mdl1->disps[o1], mdl2->disps[o2]);
+ mul_v3_fl(co2, 0.32);
+
+ add_v3_v3v3(co, co1, co2);
+
+ copy_v3_v3(mdl1->disps[o1], co);
+ copy_v3_v3(mdl2->disps[o2], co); //*/
+ }
+ }
+}
+
+void BM_loop_interp_multires(BMesh *bm, BMLoop *target, BMFace *source)
+{
+ bmesh_loop_interp_mdisps(bm, target, source);
+}
+
+void BM_loop_interp_from_face(BMesh *bm, BMLoop *target, BMFace *source,
+ int do_vertex, int do_multires)
+{
+ BMLoop *l;
+ void **blocks=NULL;
+ void **vblocks=NULL;
+ float (*cos)[3]=NULL, co[3], *w=NULL, cent[3] = {0.0f, 0.0f, 0.0f};
+ BLI_array_staticdeclare(cos, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(w, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(blocks, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(vblocks, BM_NGON_STACK_SIZE);
- the 2d projected coords are the same and faster to compute */
- xn= fabsf(source->no[0]);
- yn= fabsf(source->no[1]);
- zn= fabsf(source->no[2]);
- if (zn >= xn && zn >= yn) { ax= 0; ay= 1; }
- else if (yn >= xn && yn >= zn) { ax= 0; ay= 2; }
- else { ax= 1; ay= 2; }
-
++ int i, ax, ay;
+
+ BM_Copy_Attributes(bm, bm, source, target->f);
+
+ l = bm_firstfaceloop(source);
+ do {
+ BLI_array_growone(cos);
+ copy_v3_v3(cos[BLI_array_count(cos)-1], l->v->co);
+ add_v3_v3(cent, cos[BLI_array_count(cos)-1]);
+
+ BLI_array_append(w, 0.0f);
+ BLI_array_append(blocks, l->head.data);
+
+ if (do_vertex)
+ BLI_array_append(vblocks, l->v->head.data);
+
+ l = l->next;
+ } while (l != bm_firstfaceloop(source));
+
+ /* find best projection of face XY, XZ or YZ: barycentric weights of
++ * the 2d projected coords are the same and faster to compute */
++
++ axis_dominant_v3(&ax, &ay, source->no);
++
+ /* scale source face coordinates a bit, so points sitting directonly on an
+ edge will work.*/
+ mul_v3_fl(cent, 1.0f/(float)source->len);
+ for (i=0; i<source->len; i++) {
+ float vec[3], tmp[3];
+ sub_v3_v3v3(vec, cent, cos[i]);
+ mul_v3_fl(vec, 0.001);
+ add_v3_v3(cos[i], vec);
+
+ copy_v3_v3(tmp, cos[i]);
+ cos[i][0] = tmp[ax];
+ cos[i][1] = tmp[ay];
+ cos[i][2] = 0.0;
+ }
+
+
+ /*interpolate*/
+ co[0] = target->v->co[ax];
+ co[1] = target->v->co[ay];
+ co[2] = 0.0f;
+
+ interp_weights_poly_v3(w, cos, source->len, co);
+ CustomData_bmesh_interp(&bm->ldata, blocks, w, NULL, source->len, target->head.data);
+ if (do_vertex) {
+ CustomData_bmesh_interp(&bm->vdata, vblocks, w, NULL, source->len, target->v->head.data);
+ BLI_array_free(vblocks);
+ }
+
+ BLI_array_free(cos);
+ BLI_array_free(w);
+ BLI_array_free(blocks);
+
+ if (do_multires) {
+ if (CustomData_has_layer(&bm->ldata, CD_MDISPS)) {
+ bmesh_loop_interp_mdisps(bm, target, source);
+ }
+ }
+}
+
+
+void BM_vert_interp_from_face(BMesh *bm, BMVert *v, BMFace *source)
+{
+ BMLoop *l;
+ void **blocks=NULL;
+ float (*cos)[3]=NULL, *w=NULL, cent[3] = {0.0f, 0.0f, 0.0f};
+ BLI_array_staticdeclare(cos, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(w, BM_NGON_STACK_SIZE);
+ BLI_array_staticdeclare(blocks, BM_NGON_STACK_SIZE);
+ int i;
+
+ l = bm_firstfaceloop(source);
+ do {
+ BLI_array_growone(cos);
+ copy_v3_v3(cos[BLI_array_count(cos)-1], l->v->co);
+ add_v3_v3(cent, cos[BLI_array_count(cos)-1]);
+
+ BLI_array_append(w, 0.0f);
+ BLI_array_append(blocks, l->v->head.data);
+ l = l->next;
+ } while (l != bm_firstfaceloop(source));
+
+ /*scale source face coordinates a bit, so points sitting directonly on an
+ edge will work.*/
+ mul_v3_fl(cent, 1.0f/(float)source->len);
+ for (i=0; i<source->len; i++) {
+ float vec[3];
+ sub_v3_v3v3(vec, cent, cos[i]);
+ mul_v3_fl(vec, 0.01);
+ add_v3_v3(cos[i], vec);
+ }
+
+ /*interpolate*/
+ interp_weights_poly_v3(w, cos, source->len, v->co);
+ CustomData_bmesh_interp(&bm->vdata, blocks, w, NULL, source->len, v->head.data);
+
+ BLI_array_free(cos);
+ BLI_array_free(w);
+ BLI_array_free(blocks);
+}
+
+static void update_data_blocks(BMesh *bm, CustomData *olddata, CustomData *data)
+{
+ BMIter iter;
+ BLI_mempool *oldpool = olddata->pool;
+ void *block;
+
+ CustomData_bmesh_init_pool(data, data==&bm->ldata ? 2048 : 512);
+
+ if (data == &bm->vdata) {
+ BMVert *eve;
+
+ BM_ITER(eve, &iter, bm, BM_VERTS_OF_MESH, NULL) {
+ block = NULL;
+ CustomData_bmesh_set_default(data, &block);
+ CustomData_bmesh_copy_data(olddata, data, eve->head.data, &block);
+ CustomData_bmesh_free_block(olddata, &eve->head.data);
+ eve->head.data= block;
+ }
+ }
+ else if (data == &bm->edata) {
+ BMEdge *eed;
+
+ BM_ITER(eed, &iter, bm, BM_EDGES_OF_MESH, NULL) {
+ block = NULL;
+ CustomData_bmesh_set_default(data, &block);
+ CustomData_bmesh_copy_data(olddata, data, eed->head.data, &block);
+ CustomData_bmesh_free_block(olddata, &eed->head.data);
+ eed->head.data= block;
+ }
+ }
+ else if (data == &bm->pdata || data == &bm->ldata) {
+ BMIter liter;
+ BMFace *efa;
+ BMLoop *l;
+
+ BM_ITER(efa, &iter, bm, BM_FACES_OF_MESH, NULL) {
+ if (data == &bm->pdata) {
+ block = NULL;
+ CustomData_bmesh_set_default(data, &block);
+ CustomData_bmesh_copy_data(olddata, data, efa->head.data, &block);
+ CustomData_bmesh_free_block(olddata, &efa->head.data);
+ efa->head.data= block;
+ }
+
+ if (data == &bm->ldata) {
+ BM_ITER(l, &liter, bm, BM_LOOPS_OF_FACE, efa) {
+ block = NULL;
+ CustomData_bmesh_set_default(data, &block);
+ CustomData_bmesh_copy_data(olddata, data, l->head.data, &block);
+ CustomData_bmesh_free_block(olddata, &l->head.data);
+ l->head.data= block;
+ }
+ }
+ }
+ }
+
+ if (oldpool) {
+ /* this should never happen but can when dissolve fails - [#28960] */
+ BLI_assert(data->pool != oldpool);
+
+ BLI_mempool_destroy(oldpool);
+ }
+}
+
+void BM_add_data_layer(BMesh *bm, CustomData *data, int type)
+{
+ CustomData olddata;
+
+ olddata= *data;
+ olddata.layers= (olddata.layers)? MEM_dupallocN(olddata.layers): NULL;
+
+ /* the pool is now owned by olddata and must not be shared */
+ data->pool = NULL;
+
+ CustomData_add_layer(data, type, CD_DEFAULT, NULL, 0);
+
+ update_data_blocks(bm, &olddata, data);
+ if (olddata.layers) MEM_freeN(olddata.layers);
+}
+
+void BM_add_data_layer_named(BMesh *bm, CustomData *data, int type, const char *name)
+{
+ CustomData olddata;
+
+ olddata= *data;
+ olddata.layers= (olddata.layers)? MEM_dupallocN(olddata.layers): NULL;
+
+ /* the pool is now owned by olddata and must not be shared */
+ data->pool = NULL;
+
+ CustomData_add_layer_named(data, type, CD_DEFAULT, NULL, 0, name);
+
+ update_data_blocks(bm, &olddata, data);
+ if (olddata.layers) MEM_freeN(olddata.layers);
+}
+
+void BM_free_data_layer(BMesh *bm, CustomData *data, int type)
+{
+ CustomData olddata;
+
+ olddata= *data;
+ olddata.layers= (olddata.layers)? MEM_dupallocN(olddata.layers): NULL;
+
+ /* the pool is now owned by olddata and must not be shared */
+ data->pool = NULL;
+
+ CustomData_free_layer_active(data, type, 0);
+
+ update_data_blocks(bm, &olddata, data);
+ if (olddata.layers) MEM_freeN(olddata.layers);
+}
+
+void BM_free_data_layer_n(BMesh *bm, CustomData *data, int type, int n)
+{
+ CustomData olddata;
+
+ olddata= *data;
+ olddata.layers= (olddata.layers)? MEM_dupallocN(olddata.layers): NULL;
+
+ /* the pool is now owned by olddata and must not be shared */
+ data->pool = NULL;
+
+ CustomData_free_layer(data, type, 0, CustomData_get_layer_index_n(data, type, n));
+
+ update_data_blocks(bm, &olddata, data);
+ if (olddata.layers) MEM_freeN(olddata.layers);
+}
+
+float BM_GetCDf(CustomData *cd, void *element, int type)
+{
+ if (CustomData_has_layer(cd, type)) {
+ float *f = CustomData_bmesh_get(cd, ((BMHeader*)element)->data, type);
+ return *f;
+ }
+
+ return 0.0;
+}
+
+void BM_SetCDf(CustomData *cd, void *element, int type, float val)
+{
+ if (CustomData_has_layer(cd, type)) {
+ float *f = CustomData_bmesh_get(cd, ((BMHeader*)element)->data, type);
+ *f = val;
+ }
+
+ return;
+}
--- /dev/null
- int xn, yn, zn, ax, ay;
+/*
+ * ***** 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.
+ *
+ * Contributor(s): Joseph Eagar, Geoffrey Bantle, Campbell Barton
+ *
+ * ***** END GPL LICENSE BLOCK *****
+ */
+
+/** \file blender/bmesh/intern/bmesh_polygon.c
+ * \ingroup bmesh
+ *
+ * This file contains code for dealing
+ * with polygons (normal/area calculation,
+ * tesselation, etc)
+ *
+ * BMESH_TODO:
+ * - Add in Tesselator frontend that creates
+ * BMTriangles from copied faces
+ *
+ * - Add in Function that checks for and flags
+ * degenerate faces.
+ */
+
+#include <string.h>
+#include <math.h>
+#include <stdlib.h>
+
+#include "BKE_utildefines.h"
+
+#include "BLI_math.h"
+#include "BLI_blenlib.h"
+#include "BLI_array.h"
+#include "BLI_utildefines.h"
+
+#include "MEM_guardedalloc.h"
+
+#include "bmesh.h"
+#include "bmesh_private.h"
+
+/*
+ * TEST EDGE SIDE and POINT IN TRIANGLE
+ *
+ * Point in triangle tests stolen from scanfill.c.
+ * Used for tesselator
+ *
+*/
+
+static short testedgeside(double *v1, double *v2, double *v3)
+/* is v3 to the right of v1-v2 ? With exception: v3==v1 || v3==v2 */
+{
+ double inp;
+
+ //inp= (v2[cox]-v1[cox])*(v1[coy]-v3[coy]) +(v1[coy]-v2[coy])*(v1[cox]-v3[cox]);
+ inp= (v2[0]-v1[0])*(v1[1]-v3[1]) +(v1[1]-v2[1])*(v1[0]-v3[0]);
+
+ if(inp<0.0) return 0;
+ else if(inp==0) {
+ if(v1[0]==v3[0] && v1[1]==v3[1]) return 0;
+ if(v2[0]==v3[0] && v2[1]==v3[1]) return 0;
+ }
+ return 1;
+}
+
+static short testedgesidef(float *v1, float *v2, float *v3)
+/* is v3 to the right of v1-v2 ? With exception: v3==v1 || v3==v2 */
+{
+ double inp;
+
+ //inp= (v2[cox]-v1[cox])*(v1[coy]-v3[coy]) +(v1[coy]-v2[coy])*(v1[cox]-v3[cox]);
+ inp= (v2[0]-v1[0])*(v1[1]-v3[1]) +(v1[1]-v2[1])*(v1[0]-v3[0]);
+
+ if(inp<0.0) return 0;
+ else if(inp==0) {
+ if(v1[0]==v3[0] && v1[1]==v3[1]) return 0;
+ if(v2[0]==v3[0] && v2[1]==v3[1]) return 0;
+ }
+ return 1;
+}
+
+static int point_in_triangle(double *v1, double *v2, double *v3, double *pt)
+{
+ if(testedgeside(v1,v2,pt) && testedgeside(v2,v3,pt) && testedgeside(v3,v1,pt))
+ return 1;
+ return 0;
+}
+
+/*
+ * COMPUTE POLY NORMAL
+ *
+ * Computes the normal of a planar
+ * polygon See Graphics Gems for
+ * computing newell normal.
+ *
+*/
+
+static void compute_poly_normal(float normal[3], float (*verts)[3], int nverts)
+{
+
+ float u[3], v[3], w[3];/*, *w, v1[3], v2[3];*/
+ float n[3]= {0.0f, 0.0f, 0.0f} /*, l, v1[3], v2[3] */;
+ /* double l2; */
+ int i /*, s=0 */;
+
+ /*this fixes some weird numerical error*/
+ add_v3_fl(verts[0], 0.0001f);
+
+ for(i = 0; i < nverts; i++){
+ copy_v3_v3(u, verts[i]);
+ copy_v3_v3(v, verts[(i+1) % nverts]);
+ copy_v3_v3(w, verts[(i+2) % nverts]);
+
+#if 0
+ sub_v3_v3v3(v1, w, v);
+ sub_v3_v3v3(v2, v, u);
+ normalize_v3(v1);
+ normalize_v3(v2);
+
+ l = dot_v3v3(v1, v2);
+ if (fabsf(l-1.0) < 0.1f) {
+ continue;
+ }
+#endif
+ /* newell's method
+
+ so thats?:
+ (a[1] - b[1]) * (a[2] + b[2]);
+ a[1]*b[2] - b[1]*a[2] - b[1]*b[2] + a[1]*a[2]
+
+ odd. half of that is the cross product. . .what's the
+ other half?
+
+ also could be like a[1]*(b[2] + a[2]) - b[1]*(a[2] - b[2])
+ */
+
+ n[0] += (u[1] - v[1]) * (u[2] + v[2]);
+ n[1] += (u[2] - v[2]) * (u[0] + v[0]);
+ n[2] += (u[0] - v[0]) * (u[1] + v[1]);
+ }
+
+ if(normalize_v3_v3(normal, n) == 0.0f) {
+ normal[2] = 1.0f; /* other axis set to 0.0 */
+ }
+
+#if 0
+ l = len_v3(n);
+ /*fast square root, newton/babylonian method:
+ l2 = l*0.1;
+
+ l2 = (l/l2 + l2)*0.5;
+ l2 = (l/l2 + l2)*0.5;
+ l2 = (l/l2 + l2)*0.5;
+ */
+
+ if (l == 0.0) {
+ normal[0] = 0.0f;
+ normal[1] = 0.0f;
+ normal[2] = 1.0f;
+
+ return;
+ } else l = 1.0f / l;
+
+ mul_v3_fl(n, l);
+
+ copy_v3_v3(normal, n);
+#endif
+}
+
+/*
+ * COMPUTE POLY CENTER
+ *
+ * Computes the centroid and
+ * area of a polygon in the X/Y
+ * plane.
+ *
+*/
+
+static int compute_poly_center(float center[3], float *area, float (*verts)[3], int nverts)
+{
+ int i, j;
+ float atmp = 0.0, xtmp = 0.0, ytmp = 0.0, ai;
+
+ zero_v3(center);
+
+ if(nverts < 3)
+ return 0;
+
+ i = nverts-1;
+ j = 0;
+
+ while(j < nverts){
+ ai = verts[i][0] * verts[j][1] - verts[j][0] * verts[i][1];
+ atmp += ai;
+ xtmp += (verts[j][0] + verts[i][0]) * ai;
+ ytmp += (verts[j][1] + verts[i][1]) * ai;
+ i = j;
+ j += 1;
+ }
+
+ if(area)
+ *area = atmp / 2.0f;
+
+ if (atmp != 0){
+ center[0] = xtmp / (3.0f * atmp);
+ center[1] = xtmp / (3.0f * atmp);
+ return 1;
+ }
+ return 0;
+}
+
+float BM_face_area(BMFace *f)
+{
+ BMLoop *l;
+ BMIter iter;
+ float (*verts)[3];
+ float area, center[3];
+ int i;
+
+ BLI_array_fixedstack_declare(verts, BM_NGON_STACK_SIZE, f->len, __func__);
+
+ i = 0;
+ BM_ITER(l, &iter, NULL, BM_LOOPS_OF_FACE, f) {
+ copy_v3_v3(verts[i], l->v->co);
+ i++;
+ }
+
+ compute_poly_center(center, &area, verts, f->len);
+
+ BLI_array_fixedstack_free(verts);
+
+ return area;
+}
+/*
+computes center of face in 3d. uses center of bounding box.
+*/
+
+void BM_Compute_Face_CenterBounds(BMesh *bm, BMFace *f, float r_cent[3])
+{
+ BMIter iter;
+ BMLoop *l;
+ float min[3], max[3];
+ int i;
+
+ INIT_MINMAX(min, max);
+ l = BMIter_New(&iter, bm, BM_LOOPS_OF_FACE, f);
+ for (i=0; l; l=BMIter_Step(&iter), i++) {
+ DO_MINMAX(l->v->co, min, max);
+ }
+
+ mid_v3_v3v3(r_cent, min, max);
+}
+
+void BM_Compute_Face_CenterMean(BMesh *bm, BMFace *f, float r_cent[3])
+{
+ BMIter iter;
+ BMLoop *l;
+ int i;
+
+ zero_v3(r_cent);
+
+ l = BMIter_New(&iter, bm, BM_LOOPS_OF_FACE, f);
+ for (i=0; l; l=BMIter_Step(&iter), i++) {
+ add_v3_v3(r_cent, l->v->co);
+ }
+
+ mul_v3_fl(r_cent, 1.0f / (float)f->len);
+}
+
+/*
+ * COMPUTE POLY PLANE
+ *
+ * Projects a set polygon's vertices to
+ * a plane defined by the average
+ * of its edges cross products
+ *
+*/
+
+void compute_poly_plane(float (*verts)[3], int nverts)
+{
+
+ float avgc[3], norm[3], mag, avgn[3];
+ float *v1, *v2, *v3;
+ int i;
+
+ if(nverts < 3)
+ return;
+
+ zero_v3(avgn);
+ zero_v3(avgc);
+
+ for(i = 0; i < nverts; i++){
+ v1 = verts[i];
+ v2 = verts[(i+1) % nverts];
+ v3 = verts[(i+2) % nverts];
+ normal_tri_v3( norm,v1, v2, v3);
+
+ add_v3_v3(avgn, norm);
+ }
+
+ /*what was this bit for?*/
+ if (is_zero_v3(avgn)) {
+ avgn[0] = 0.0f;
+ avgn[1] = 0.0f;
+ avgn[2] = 1.0f;
+ } else {
+ /* XXX, why is this being divided and _then_ normalized
+ * division could be removed? - campbell */
+ avgn[0] /= nverts;
+ avgn[1] /= nverts;
+ avgn[2] /= nverts;
+ normalize_v3(avgn);
+ }
+
+ for(i = 0; i < nverts; i++){
+ v1 = verts[i];
+ mag = dot_v3v3(v1, avgn);
+ madd_v3_v3fl(v1, avgn, -mag);
+ }
+}
+
+/*
+ BM LEGAL EDGES
+
+ takes in a face and a list of edges, and sets to NULL any edge in
+ the list that bridges a concave region of the face or intersects
+ any of the faces's edges.
+*/
+#if 0 /* needs BLI math double versions of these functions */
+static void shrink_edged(double *v1, double *v2, double fac)
+{
+ double mid[3];
+
+ mid_v3_v3v3(mid, v1, v2);
+
+ sub_v3_v3v3(v1, v1, mid);
+ sub_v3_v3v3(v2, v2, mid);
+
+ mul_v3_fl(v1, fac);
+ mul_v3_fl(v2, fac);
+
+ add_v3_v3v3(v1, v1, mid);
+ add_v3_v3v3(v2, v2, mid);
+}
+#endif
+
+static void shrink_edgef(float *v1, float *v2, float fac)
+{
+ float mid[3];
+
+ mid_v3_v3v3(mid, v1, v2);
+
+ sub_v3_v3v3(v1, v1, mid);
+ sub_v3_v3v3(v2, v2, mid);
+
+ mul_v3_fl(v1, fac);
+ mul_v3_fl(v2, fac);
+
+ add_v3_v3v3(v1, v1, mid);
+ add_v3_v3v3(v2, v2, mid);
+}
+
+
+/*
+ * POLY ROTATE PLANE
+ *
+ * Rotates a polygon so that it's
+ * normal is pointing towards the mesh Z axis
+ *
+*/
+
+void poly_rotate_plane(const float normal[3], float (*verts)[3], const int nverts)
+{
+
+ float up[3] = {0.0f,0.0f,1.0f}, axis[3], q[4];
+ float mat[3][3];
+ double angle;
+ int i;
+
+ cross_v3_v3v3(axis, normal, up);
+
+ angle = saacos(dot_v3v3(normal, up));
+
+ if (angle == 0.0) return;
+
+ axis_angle_to_quat(q, axis, (float)angle);
+ quat_to_mat3(mat, q);
+
+ for(i = 0; i < nverts; i++)
+ mul_m3_v3(mat, verts[i]);
+}
+
+/*
+ * BMESH UPDATE FACE NORMAL
+ *
+ * Updates the stored normal for the
+ * given face. Requires that a buffer
+ * of sufficient length to store projected
+ * coordinates for all of the face's vertices
+ * is passed in as well.
+ *
+*/
+
+void BM_Face_UpdateNormal(BMesh *bm, BMFace *f)
+{
+ if (f->len >= 3) {
+ float (*proj)[3];
+
+ BLI_array_fixedstack_declare(proj, BM_NGON_STACK_SIZE, f->len, __func__);
+
+ bmesh_update_face_normal(bm, f, proj);
+
+ BLI_array_fixedstack_free(proj);
+ }
+}
+
+void BM_Edge_UpdateNormals(BMesh *bm, BMEdge *e)
+{
+ BMIter iter;
+ BMFace *f;
+
+ f = BMIter_New(&iter, bm, BM_FACES_OF_EDGE, e);
+ for (; f; f=BMIter_Step(&iter)) {
+ BM_Face_UpdateNormal(bm, f);
+ }
+
+ BM_Vert_UpdateNormal(bm, e->v1);
+ BM_Vert_UpdateNormal(bm, e->v2);
+}
+
+void BM_Vert_UpdateNormal(BMesh *bm, BMVert *v)
+{
+ BMIter eiter, liter;
+ BMEdge *e;
+ BMLoop *l;
+ float vec1[3], vec2[3], fac;
+ int len=0;
+
+ zero_v3(v->no);
+
+ BM_ITER(e, &eiter, bm, BM_EDGES_OF_VERT, v) {
+ BM_ITER(l, &liter, bm, BM_LOOPS_OF_EDGE, e) {
+ if (l->v == v) {
+ /* Same calculation used in BM_Compute_Normals */
+ sub_v3_v3v3(vec1, l->v->co, l->prev->v->co);
+ sub_v3_v3v3(vec2, l->next->v->co, l->v->co);
+ normalize_v3(vec1);
+ normalize_v3(vec2);
+
+ fac = saacos(-dot_v3v3(vec1, vec2));
+
+ madd_v3_v3fl(v->no, l->f->no, fac);
+
+ len++;
+ }
+ }
+ }
+
+ if (len) {
+ normalize_v3(v->no);
+ }
+}
+
+void BM_Vert_UpdateAllNormals(BMesh *bm, BMVert *v)
+{
+ BMIter iter;
+ BMFace *f;
+ int len=0;
+
+ f = BMIter_New(&iter, bm, BM_FACES_OF_VERT, v);
+ for (; f; f=BMIter_Step(&iter), len++) {
+ BM_Face_UpdateNormal(bm, f);
+ }
+
+ BM_Vert_UpdateNormal(bm, v);
+}
+
+void bmesh_update_face_normal(BMesh *bm, BMFace *f, float (*projectverts)[3])
+{
+ BMIter iter;
+ BMLoop *l;
+
+ if(f->len > 4) {
+ int i = 0;
+ BM_ITER(l, &iter, bm, BM_LOOPS_OF_FACE, f) {
+ copy_v3_v3(projectverts[i], l->v->co);
+ i += 1;
+ }
+
+ compute_poly_normal(f->no, projectverts, f->len);
+ }
+ else if(f->len == 3){
+ BMVert *v1, *v2, *v3;
+ v1 = bm_firstfaceloop(f)->v;
+ v2 = bm_firstfaceloop(f)->next->v;
+ v3 = bm_firstfaceloop(f)->next->next->v;
+ normal_tri_v3( f->no,v1->co, v2->co, v3->co);
+ }
+ else if(f->len == 4){
+ BMVert *v1, *v2, *v3, *v4;
+ v1 = bm_firstfaceloop(f)->v;
+ v2 = bm_firstfaceloop(f)->next->v;
+ v3 = bm_firstfaceloop(f)->next->next->v;
+ v4 = bm_firstfaceloop(f)->prev->v;
+ normal_quad_v3( f->no,v1->co, v2->co, v3->co, v4->co);
+ }
+ else{ /*horrible, two sided face!*/
+ zero_v3(f->no);
+ }
+
+}
+
+
+/*
+ * BMESH FLIP NORMAL
+ *
+ * Reverses the winding of a face.
+ * Note that this updates the calculated
+ * normal.
+*/
+void BM_flip_normal(BMesh *bm, BMFace *f)
+{
+ bmesh_loop_reverse(bm, f);
+ negate_v3(f->no);
+}
+
+/* detects if two line segments cross each other (intersects).
+ note, there could be more winding cases then there needs to be. */
+static int UNUSED_FUNCTION(linecrosses)(double *v1, double *v2, double *v3, double *v4)
+{
+ int w1, w2, w3, w4, w5;
+
+ /*w1 = winding(v1, v3, v4);
+ w2 = winding(v2, v3, v4);
+ w3 = winding(v3, v1, v2);
+ w4 = winding(v4, v1, v2);
+
+ return (w1 == w2) && (w3 == w4);*/
+
+ w1 = testedgeside(v1, v3, v2);
+ w2 = testedgeside(v2, v4, v1);
+ w3 = !testedgeside(v1, v2, v3);
+ w4 = testedgeside(v3, v2, v4);
+ w5 = !testedgeside(v3, v1, v4);
+ return w1 == w2 && w2 == w3 && w3 == w4 && w4==w5;
+}
+
+/* detects if two line segments cross each other (intersects).
+ note, there could be more winding cases then there needs to be. */
+static int linecrossesf(float *v1, float *v2, float *v3, float *v4)
+{
+ int w1, w2, w3, w4, w5 /*, ret*/;
+ float mv1[2], mv2[2], mv3[2], mv4[2];
+
+ /*now test winding*/
+ w1 = testedgesidef(v1, v3, v2);
+ w2 = testedgesidef(v2, v4, v1);
+ w3 = !testedgesidef(v1, v2, v3);
+ w4 = testedgesidef(v3, v2, v4);
+ w5 = !testedgesidef(v3, v1, v4);
+
+ if (w1 == w2 && w2 == w3 && w3 == w4 && w4==w5)
+ return 1;
+
+#define GETMIN2_AXIS(a, b, ma, mb, axis) ma[axis] = MIN2(a[axis], b[axis]), mb[axis] = MAX2(a[axis], b[axis])
+#define GETMIN2(a, b, ma, mb) GETMIN2_AXIS(a, b, ma, mb, 0); GETMIN2_AXIS(a, b, ma, mb, 1);
+
+ GETMIN2(v1, v2, mv1, mv2);
+ GETMIN2(v3, v4, mv3, mv4);
+
+ /*do an interval test on the x and y axes*/
+ /*first do x axis*/
+ #define T FLT_EPSILON*15
+ if (ABS(v1[1]-v2[1]) < T && ABS(v3[1]-v4[1]) < T &&
+ ABS(v1[1]-v3[1]) < T)
+ {
+ return (mv4[0] >= mv1[0] && mv3[0] <= mv2[0]);
+ }
+
+ /*now do y axis*/
+ if (ABS(v1[0]-v2[0]) < T && ABS(v3[0]-v4[0]) < T &&
+ ABS(v1[0]-v3[0]) < T)
+ {
+ return (mv4[1] >= mv1[1] && mv3[1] <= mv2[1]);
+ }
+
+ return 0;
+}
+
+/*
+ BM POINT IN FACE
+
+ Projects co onto face f, and returns true if it is inside
+ the face bounds. Note that this uses a best-axis projection
+ test, instead of projecting co directly into f's orientation
+ space, so there might be accuracy issues.
+*/
+int BM_Point_In_Face(BMesh *bm, BMFace *f, const float co[3])
+{
- the 2d projected coords are the same and faster to compute
-
- this probably isn't all that accurate, but it has the advantage of
- being fast (especially compared to projecting into the face orientation)
- */
- xn= fabsf(f->no[0]);
- yn= fabsf(f->no[1]);
- zn= fabsf(f->no[2]);
-
- if (zn >= xn && zn >= yn) { ax= 0; ay= 1; }
- else if (yn >= xn && yn >= zn) { ax= 0; ay= 2; }
- else { ax= 1; ay= 2; }
++ int ax, ay;
+ float co2[3], cent[3] = {0.0f, 0.0f, 0.0f}, out[3] = {FLT_MAX*0.5f, FLT_MAX*0.5f, 0};
+ BMLoop *l;
+ int crosses = 0;
+ float eps = 1.0f+(float)FLT_EPSILON*150.0f;
+
+ if (dot_v3v3(f->no, f->no) <= FLT_EPSILON*10)
+ BM_Face_UpdateNormal(bm, f);
+
+ /* find best projection of face XY, XZ or YZ: barycentric weights of
++ * the 2d projected coords are the same and faster to compute
++ *
++ * this probably isn't all that accurate, but it has the advantage of
++ * being fast (especially compared to projecting into the face orientation)
++ */
++ axis_dominant_v3(&ax, &ay, f->no);
+
+ co2[0] = co[ax];
+ co2[1] = co[ay];
+ co2[2] = 0;
+
+ l = bm_firstfaceloop(f);
+ do {
+ cent[0] += l->v->co[ax];
+ cent[1] += l->v->co[ay];
+ l = l->next;
+ } while (l != bm_firstfaceloop(f));
+
+ mul_v2_fl(cent, 1.0f/(float)f->len);
+
+ l = bm_firstfaceloop(f);
+ do {
+ float v1[3], v2[3];
+
+ v1[0] = (l->prev->v->co[ax] - cent[ax])*eps + cent[ax];
+ v1[1] = (l->prev->v->co[ay] - cent[ay])*eps + cent[ay];
+ v1[2] = 0.0f;
+
+ v2[0] = (l->v->co[ax] - cent[ax])*eps + cent[ax];
+ v2[1] = (l->v->co[ay] - cent[ay])*eps + cent[ay];
+ v2[2] = 0.0f;
+
+ crosses += linecrossesf(v1, v2, co2, out) != 0;
+
+ l = l->next;
+ } while (l != bm_firstfaceloop(f));
+
+ return crosses%2 != 0;
+}
+
+static int goodline(float (*projectverts)[3], BMFace *f, int v1i,
+ int v2i, int v3i, int UNUSED(nvert))
+{
+ BMLoop *l = bm_firstfaceloop(f);
+ double v1[3], v2[3], v3[3], pv1[3], pv2[3];
+ int i;
+
+ VECCOPY(v1, projectverts[v1i]);
+ VECCOPY(v2, projectverts[v2i]);
+ VECCOPY(v3, projectverts[v3i]);
+
+ if (testedgeside(v1, v2, v3)) return 0;
+
+ //for (i=0; i<nvert; i++) {
+ do {
+ i = BM_GetIndex(l->v);
+ if (i == v1i || i == v2i || i == v3i) {
+ l = l->next;
+ continue;
+ }
+
+ VECCOPY(pv1, projectverts[BM_GetIndex(l->v)]);
+ VECCOPY(pv2, projectverts[BM_GetIndex(l->next->v)]);
+
+ //if (linecrosses(pv1, pv2, v1, v3)) return 0;
+ if (point_in_triangle(v1, v2, v3, pv1)) return 0;
+ if (point_in_triangle(v3, v2, v1, pv1)) return 0;
+
+ l = l->next;
+ } while (l != bm_firstfaceloop(f));
+ return 1;
+}
+/*
+ * FIND EAR
+ *
+ * Used by tesselator to find
+ * the next triangle to 'clip off'
+ * of a polygon while tesselating.
+ *
+*/
+
+static BMLoop *find_ear(BMesh *UNUSED(bm), BMFace *f, float (*verts)[3],
+ int nvert)
+{
+ BMVert *v1, *v2, *v3;
+ BMLoop *bestear = NULL, *l;
+ /*float angle, bestangle = 180.0f;*/
+ int isear /*, i=0*/;
+
+ l = bm_firstfaceloop(f);
+ do {
+ isear = 1;
+
+ v1 = l->prev->v;
+ v2 = l->v;
+ v3 = l->next->v;
+
+ if (BM_Edge_Exist(v1, v3)) isear = 0;
+
+ if (isear && !goodline(verts, f, BM_GetIndex(v1), BM_GetIndex(v2),
+ BM_GetIndex(v3), nvert))
+ isear = 0;
+
+ if(isear) {
+ /*angle = angle_v3v3v3(verts[v1->head.eflag2], verts[v2->head.eflag2], verts[v3->head.eflag2]);
+ if(!bestear || ABS(angle-45.0f) < bestangle) {
+ bestear = l;
+ bestangle = ABS(45.0f-angle);
+ }
+
+ if (angle > 20 && angle < 90) break;
+ if (angle < 100 && i > 5) break;
+ i += 1;*/
+ bestear = l;
+ break;
+ }
+ l = l->next;
+ }
+ while(l != bm_firstfaceloop(f));
+
+ return bestear;
+}
+
+/*
+ * BMESH TRIANGULATE FACE
+ *
+ * Triangulates a face using a
+ * simple 'ear clipping' algorithm
+ * that tries to favor non-skinny
+ * triangles (angles less than
+ * 90 degrees). If the triangulator
+ * has bits left over (or cannot
+ * triangulate at all) it uses a
+ * simple fan triangulation
+ *
+ * newfaces, if non-null, must be an array of BMFace pointers,
+ * with a length equal to f->len. it will be filled with the new
+ * triangles, and will be NULL-terminated.
+*/
+void BM_Triangulate_Face(BMesh *bm, BMFace *f, float (*projectverts)[3],
+ int newedgeflag, int newfaceflag, BMFace **newfaces)
+{
+ int i, done, nvert, nf_i = 0;
+ BMLoop *l, *newl, *nextloop;
+ /* BMVert *v; */ /* UNUSED */
+
+ /*copy vertex coordinates to vertspace array*/
+ i = 0;
+ l = bm_firstfaceloop(f);
+ do{
+ copy_v3_v3(projectverts[i], l->v->co);
+ BM_SetIndex(l->v, i); /* set dirty! */
+ i++;
+ l = l->next;
+ }while(l != bm_firstfaceloop(f));
+
+ bm->elem_index_dirty |= BM_VERT; /* see above */
+
+ ///bmesh_update_face_normal(bm, f, projectverts);
+
+ compute_poly_normal(f->no, projectverts, f->len);
+ poly_rotate_plane(f->no, projectverts, i);
+
+ nvert = f->len;
+
+ //compute_poly_plane(projectverts, i);
+ for (i=0; i<nvert; i++) {
+ projectverts[i][2] = 0.0f;
+ }
+
+ done = 0;
+ while(!done && f->len > 3){
+ done = 1;
+ l = find_ear(bm, f, projectverts, nvert);
+ if(l) {
+ done = 0;
+ /* v = l->v; */ /* UNUSED */
+ f = BM_Split_Face(bm, l->f, l->prev->v,
+ l->next->v,
+ &newl, NULL);
+ copy_v3_v3(f->no, l->f->no);
+
+ if (!f) {
+ fprintf(stderr, "%s: triangulator failed to split face! (bmesh internal error)\n", __func__);
+ break;
+ }
+
+ BMO_SetFlag(bm, newl->e, newedgeflag);
+ BMO_SetFlag(bm, f, newfaceflag);
+
+ if (newfaces) newfaces[nf_i++] = f;
+
+ /*l = f->loopbase;
+ do {
+ if (l->v == v) {
+ f->loopbase = l;
+ break;
+ }
+ l = l->next;
+ } while (l != f->loopbase);*/
+ }
+ }
+
+ if (f->len > 3){
+ l = bm_firstfaceloop(f);
+ while (l->f->len > 3){
+ nextloop = l->next->next;
+ f = BM_Split_Face(bm, l->f, l->v, nextloop->v,
+ &newl, NULL);
+ if (!f) {
+ printf("triangle fan step of triangulator failed.\n");
+
+ /*NULL-terminate*/
+ if (newfaces) newfaces[nf_i] = NULL;
+ return;
+ }
+
+ if (newfaces) newfaces[nf_i++] = f;
+
+ BMO_SetFlag(bm, newl->e, newedgeflag);
+ BMO_SetFlag(bm, f, newfaceflag);
+ l = nextloop;
+ }
+ }
+
+ /*NULL-terminate*/
+ if (newfaces) newfaces[nf_i] = NULL;
+}
+
+/*each pair of loops defines a new edge, a split. this function goes
+ through and sets pairs that are geometrically invalid to null. a
+ split is invalid, if it forms a concave angle or it intersects other
+ edges in the face, or it intersects another split. in the case of
+ intersecting splits, only the first of the set of intersecting
+ splits survives.*/
+void BM_LegalSplits(BMesh *bm, BMFace *f, BMLoop *(*loops)[2], int len)
+{
+ BMIter iter;
+ BMLoop *l;
+ float v1[3], v2[3], v3[3]/*, v4[3]*/, no[3], mid[3], *p1, *p2, *p3, *p4;
+ float out[3] = {-234324.0f, -234324.0f, 0.0f};
+ float (*projverts)[3];
+ float (*edgeverts)[3];
+ float fac1 = 1.0000001f, fac2 = 0.9f; //9999f; //0.999f;
+ int i, j, a=0, clen;
+
+ BLI_array_fixedstack_declare(projverts, BM_NGON_STACK_SIZE, f->len, "projvertsb");
+ BLI_array_fixedstack_declare(edgeverts, BM_NGON_STACK_SIZE * 2, len * 2, "edgevertsb");
+
+ i = 0;
+ l = BMIter_New(&iter, bm, BM_LOOPS_OF_FACE, f);
+ for (; l; l=BMIter_Step(&iter)) {
+ BM_SetIndex(l, i); /* set_loop */
+ copy_v3_v3(projverts[i], l->v->co);
+ i++;
+ }
+
+ for (i=0; i<len; i++) {
+ copy_v3_v3(v1, loops[i][0]->v->co);
+ copy_v3_v3(v2, loops[i][1]->v->co);
+
+ shrink_edgef(v1, v2, fac2);
+
+ copy_v3_v3(edgeverts[a], v1);
+ a++;
+ copy_v3_v3(edgeverts[a], v2);
+ a++;
+ }
+
+ compute_poly_normal(no, projverts, f->len);
+ poly_rotate_plane(no, projverts, f->len);
+ poly_rotate_plane(no, edgeverts, len*2);
+
+ l = bm_firstfaceloop(f);
+ for (i=0; i<f->len; i++) {
+ p1 = projverts[i];
+ out[0] = MAX2(out[0], p1[0]) + 0.01f;
+ out[1] = MAX2(out[1], p1[1]) + 0.01f;
+ out[2] = 0.0f;
+ p1[2] = 0.0f;
+
+ //copy_v3_v3(l->v->co, p1);
+
+ l = l->next;
+ }
+
+ for (i=0; i<len; i++) {
+ edgeverts[i*2][2] = 0.0f;
+ edgeverts[i*2+1][2] = 0.0f;
+ }
+
+ /*do convexity test*/
+ for (i=0; i<len; i++) {
+ copy_v3_v3(v2, edgeverts[i*2]);
+ copy_v3_v3(v3, edgeverts[i*2+1]);
+
+ mid_v3_v3v3(mid, v2, v3);
+
+ clen = 0;
+ for (j=0; j<f->len; j++) {
+ p1 = projverts[j];
+ p2 = projverts[(j+1)%f->len];
+
+ copy_v3_v3(v1, p1);
+ copy_v3_v3(v2, p2);
+
+ shrink_edgef(v1, v2, fac1);
+
+ if (linecrossesf(p1, p2, mid, out)) clen++;
+ }
+
+ if (clen%2 == 0) {
+ loops[i][0] = NULL;
+ }
+ }
+
+ /*do line crossing tests*/
+ for (i=0; i<f->len; i++) {
+ p1 = projverts[i];
+ p2 = projverts[(i+1)%f->len];
+
+ copy_v3_v3(v1, p1);
+ copy_v3_v3(v2, p2);
+
+ shrink_edgef(v1, v2, fac1);
+
+ for (j=0; j<len; j++) {
+ if (!loops[j][0]) continue;
+
+ p3 = edgeverts[j*2];
+ p4 = edgeverts[j*2+1];
+
+ if (linecrossesf(v1, v2, p3, p4))
+ {
+ loops[j][0] = NULL;
+ }
+ }
+ }
+
+ for (i=0; i<len; i++) {
+ for (j=0; j<len; j++) {
+ if (j == i) continue;
+ if (!loops[i][0]) continue;
+ if (!loops[j][0]) continue;
+
+ p1 = edgeverts[i*2];
+ p2 = edgeverts[i*2+1];
+ p3 = edgeverts[j*2];
+ p4 = edgeverts[j*2+1];
+
+ copy_v3_v3(v1, p1);
+ copy_v3_v3(v2, p2);
+
+ shrink_edgef(v1, v2, fac1);
+
+ if (linecrossesf(v1, v2, p3, p4)) {
+ loops[i][0]=NULL;
+ }
+ }
+ }
+
+ BLI_array_fixedstack_free(projverts);
+ BLI_array_fixedstack_free(edgeverts);
+}