b3e4e896c3e9b7760db27bfbed5ac5dc0b4356d4
[blender.git] / source / blender / bmesh / intern / bmesh_interp.c
1 /*
2  * ***** BEGIN GPL LICENSE BLOCK *****
3  *
4  * This program is free software; you can redistribute it and/or
5  * modify it under the terms of the GNU General Public License
6  * as published by the Free Software Foundation; either version 2
7  * of the License, or (at your option) any later version.
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program; if not, write to the Free Software Foundation,
16  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17  *
18  * The Original Code is Copyright (C) 2007 Blender Foundation.
19  * All rights reserved.
20  *
21  * The Original Code is: all of this file.
22  *
23  * Contributor(s): Geoffrey Bantle.
24  *
25  * ***** END GPL LICENSE BLOCK *****
26  */
27
28 /** \file blender/bmesh/intern/bmesh_interp.c
29  *  \ingroup bmesh
30  *
31  * Functions for interpolating data across the surface of a mesh.
32  */
33
34 #include "MEM_guardedalloc.h"
35
36 #include "DNA_mesh_types.h"
37 #include "DNA_meshdata_types.h"
38
39 #include "BKE_customdata.h"
40 #include "BKE_multires.h"
41
42 #include "BLI_array.h"
43 #include "BLI_math.h"
44
45 #include "bmesh.h"
46 #include "intern/bmesh_private.h"
47
48 /* edge and vertex share, currently theres no need to have different logic */
49 static void bm_data_interp_from_elem(CustomData *data_layer, BMElem *ele1, BMElem *ele2, BMElem *ele_dst, const float fac)
50 {
51         if (ele1->head.data && ele2->head.data) {
52                 /* first see if we can avoid interpolation */
53                 if (fac <= 0.0f) {
54                         if (ele1 == ele_dst) {
55                                 /* do nothing */
56                         }
57                         else {
58                                 CustomData_bmesh_free_block(data_layer, &ele_dst->head.data);
59                                 CustomData_bmesh_copy_data(data_layer, data_layer, ele1->head.data, &ele_dst->head.data);
60                         }
61                 }
62                 else if (fac >= 1.0f) {
63                         if (ele2 == ele_dst) {
64                                 /* do nothing */
65                         }
66                         else {
67                                 CustomData_bmesh_free_block(data_layer, &ele_dst->head.data);
68                                 CustomData_bmesh_copy_data(data_layer, data_layer, ele2->head.data, &ele_dst->head.data);
69                         }
70                 }
71                 else {
72                         void *src[2];
73                         float w[2];
74
75                         src[0] = ele1->head.data;
76                         src[1] = ele2->head.data;
77                         w[0] = 1.0f - fac;
78                         w[1] = fac;
79                         CustomData_bmesh_interp(data_layer, src, w, NULL, 2, ele_dst->head.data);
80                 }
81         }
82 }
83
84 /**
85  * \brief Data, Interp From Verts
86  *
87  * Interpolates per-vertex data from two sources to a target.
88  *
89  * \note This is an exact match to #BM_data_interp_from_edges
90  */
91 void BM_data_interp_from_verts(BMesh *bm, BMVert *v1, BMVert *v2, BMVert *v, const float fac)
92 {
93         bm_data_interp_from_elem(&bm->vdata, (BMElem *)v1, (BMElem *)v2, (BMElem *)v, fac);
94 }
95
96 /**
97  * \brief Data, Interp From Edges
98  *
99  * Interpolates per-edge data from two sources to a target.
100  *
101  * \note This is an exact match to #BM_data_interp_from_verts
102  */
103 void BM_data_interp_from_edges(BMesh *bm, BMEdge *e1, BMEdge *e2, BMEdge *e, const float fac)
104 {
105         bm_data_interp_from_elem(&bm->edata, (BMElem *)e1, (BMElem *)e2, (BMElem *)e, fac);
106 }
107
108 /**
109  * \brief Data Vert Average
110  *
111  * Sets all the customdata (e.g. vert, loop) associated with a vert
112  * to the average of the face regions surrounding it.
113  */
114 static void UNUSED_FUNCTION(BM_Data_Vert_Average)(BMesh *UNUSED(bm), BMFace *UNUSED(f))
115 {
116         // BMIter iter;
117 }
118
119 /**
120  * \brief Data Face-Vert Edge Interp
121  *
122  * Walks around the faces of an edge and interpolates the per-face-edge
123  * data between two sources to a target.
124  */
125 void BM_data_interp_face_vert_edge(BMesh *bm, BMVert *v1, BMVert *UNUSED(v2), BMVert *v, BMEdge *e1, const float fac)
126 {
127         void *src[2];
128         float w[2];
129         BMLoop *v1loop = NULL, *vloop = NULL, *v2loop = NULL;
130         BMLoop *l_iter = NULL;
131
132         if (!e1->l) {
133                 return;
134         }
135
136         w[1] = 1.0f - fac;
137         w[0] = fac;
138
139         l_iter = e1->l;
140         do {
141                 if (l_iter->v == v1) {
142                         v1loop = l_iter;
143                         vloop = v1loop->next;
144                         v2loop = vloop->next;
145                 }
146                 else if (l_iter->v == v) {
147                         v1loop = l_iter->next;
148                         vloop = l_iter;
149                         v2loop = l_iter->prev;
150                 }
151                 
152                 if (!v1loop || !v2loop)
153                         return;
154                 
155                 src[0] = v1loop->head.data;
156                 src[1] = v2loop->head.data;
157
158                 CustomData_bmesh_interp(&bm->ldata, src, w, NULL, 2, vloop->head.data);
159         } while ((l_iter = l_iter->radial_next) != e1->l);
160 }
161
162 /**
163  * \brief Data Interp From Face
164  *
165  * projects target onto source, and pulls interpolated customdata from
166  * source.
167  *
168  * \note Only handles loop customdata. multires is handled.
169  */
170 void BM_face_interp_from_face(BMesh *bm, BMFace *target, BMFace *source)
171 {
172         BMLoop *l_iter;
173         BMLoop *l_first;
174
175         void **blocks = NULL;
176         float (*cos)[3] = NULL, *w = NULL;
177         BLI_array_fixedstack_declare(cos,     BM_NGON_STACK_SIZE, source->len, __func__);
178         BLI_array_fixedstack_declare(w,       BM_NGON_STACK_SIZE, source->len, __func__);
179         BLI_array_fixedstack_declare(blocks,  BM_NGON_STACK_SIZE, source->len, __func__);
180         int i;
181         
182         BM_elem_attrs_copy(bm, bm, source, target);
183
184         i = 0;
185         l_iter = l_first = BM_FACE_FIRST_LOOP(source);
186         do {
187                 copy_v3_v3(cos[i], l_iter->v->co);
188                 blocks[i] = l_iter->head.data;
189                 i++;
190         } while ((l_iter = l_iter->next) != l_first);
191
192         i = 0;
193         l_iter = l_first = BM_FACE_FIRST_LOOP(target);
194         do {
195                 interp_weights_poly_v3(w, cos, source->len, l_iter->v->co);
196                 CustomData_bmesh_interp(&bm->ldata, blocks, w, NULL, source->len, l_iter->head.data);
197                 i++;
198         } while ((l_iter = l_iter->next) != l_first);
199
200         BLI_array_fixedstack_free(cos);
201         BLI_array_fixedstack_free(w);
202         BLI_array_fixedstack_free(blocks);
203 }
204
205 /**
206  * \brief Multires Interpolation
207  *
208  * mdisps is a grid of displacements, ordered thus:
209  * <pre>
210  *      v1/center----v4/next -> x
211  *          |           |
212  *          |           |
213  *       v2/prev------v3/cur
214  *          |
215  *          V
216  *          y
217  * </pre>
218  */
219 static int compute_mdisp_quad(BMLoop *l, float v1[3], float v2[3], float v3[3], float v4[3],
220                               float e1[3], float e2[3])
221 {
222         float cent[3], n[3], p[3];
223
224         /* computer center */
225         BM_face_calc_center_mean(l->f, cent);
226
227         add_v3_v3v3(p, l->prev->v->co, l->v->co);
228         mul_v3_fl(p, 0.5);
229         add_v3_v3v3(n, l->next->v->co, l->v->co);
230         mul_v3_fl(n, 0.5);
231         
232         copy_v3_v3(v1, cent);
233         copy_v3_v3(v2, p);
234         copy_v3_v3(v3, l->v->co);
235         copy_v3_v3(v4, n);
236         
237         sub_v3_v3v3(e1, v2, v1);
238         sub_v3_v3v3(e2, v3, v4);
239         
240         return 1;
241 }
242
243 /* funnily enough, I think this is identical to face_to_crn_interp, heh */
244 static float quad_coord(float aa[3], float bb[3], float cc[3], float dd[3], int a1, int a2)
245 {
246         float x, y, z, f1;
247         
248         x = aa[a1] * cc[a2] - cc[a1] * aa[a2];
249         y = aa[a1] * dd[a2] + bb[a1] * cc[a2] - cc[a1] * bb[a2] - dd[a1] * aa[a2];
250         z = bb[a1] * dd[a2] - dd[a1] * bb[a2];
251         
252         if (fabsf(2.0f * (x - y + z)) > FLT_EPSILON * 10.0f) {
253                 float f2;
254
255                 f1 = ( sqrtf(y * y - 4.0f * x * z) - y + 2.0f * z) / (2.0f * (x - y + z));
256                 f2 = (-sqrtf(y * y - 4.0f * x * z) - y + 2.0f * z) / (2.0f * (x - y + z));
257
258                 f1 = fabsf(f1);
259                 f2 = fabsf(f2);
260                 f1 = minf(f1, f2);
261                 CLAMP(f1, 0.0f, 1.0f + FLT_EPSILON);
262         }
263         else {
264                 f1 = -z / (y - 2 * z);
265                 CLAMP(f1, 0.0f, 1.0f + FLT_EPSILON);
266                 
267                 if (isnan(f1) || f1 > 1.0f || f1 < 0.0f) {
268                         int i;
269                         
270                         for (i = 0; i < 2; i++) {
271                                 if (fabsf(aa[i]) < FLT_EPSILON * 100.0f)
272                                         return aa[(i + 1) % 2] / fabsf(bb[(i + 1) % 2] - aa[(i + 1) % 2]);
273                                 if (fabsf(cc[i]) < FLT_EPSILON * 100.0f)
274                                         return cc[(i + 1) % 2] / fabsf(dd[(i + 1) % 2] - cc[(i + 1) % 2]);
275                         }
276                 }
277         }
278
279         return f1;
280 }
281
282 static int quad_co(float *x, float *y, float v1[3], float v2[3], float v3[3], float v4[3],
283                    float p[3], float n[3])
284 {
285         float projverts[5][3], n2[3];
286         float dprojverts[4][3], origin[3] = {0.0f, 0.0f, 0.0f};
287         int i;
288
289         /* project points into 2d along normal */
290         copy_v3_v3(projverts[0], v1);
291         copy_v3_v3(projverts[1], v2);
292         copy_v3_v3(projverts[2], v3);
293         copy_v3_v3(projverts[3], v4);
294         copy_v3_v3(projverts[4], p);
295
296         normal_quad_v3(n2, projverts[0], projverts[1], projverts[2], projverts[3]);
297
298         if (dot_v3v3(n, n2) < -FLT_EPSILON) {
299                 return 0;
300         }
301
302         /* rotate */
303         poly_rotate_plane(n, projverts, 5);
304         
305         /* flatten */
306         for (i = 0; i < 5; i++) {
307                 projverts[i][2] = 0.0f;
308         }
309         
310         /* subtract origin */
311         for (i = 0; i < 4; i++) {
312                 sub_v3_v3(projverts[i], projverts[4]);
313         }
314         
315         copy_v3_v3(dprojverts[0], projverts[0]);
316         copy_v3_v3(dprojverts[1], projverts[1]);
317         copy_v3_v3(dprojverts[2], projverts[2]);
318         copy_v3_v3(dprojverts[3], projverts[3]);
319
320         if (!isect_point_quad_v2(origin, dprojverts[0], dprojverts[1], dprojverts[2], dprojverts[3])) {
321                 return 0;
322         }
323         
324         *y = quad_coord(dprojverts[1], dprojverts[0], dprojverts[2], dprojverts[3], 0, 1);
325         *x = quad_coord(dprojverts[2], dprojverts[1], dprojverts[3], dprojverts[0], 0, 1);
326
327         return 1;
328 }
329
330 static void mdisp_axis_from_quad(float v1[3], float v2[3], float UNUSED(v3[3]), float v4[3],
331                                 float axis_x[3], float axis_y[3])
332 {
333         sub_v3_v3v3(axis_x, v4, v1);
334         sub_v3_v3v3(axis_y, v2, v1);
335
336         normalize_v3(axis_x);
337         normalize_v3(axis_y);
338 }
339
340 /* tl is loop to project onto, l is loop whose internal displacement, co, is being
341  * projected.  x and y are location in loop's mdisps grid of point co. */
342 static int mdisp_in_mdispquad(BMLoop *l, BMLoop *tl, float p[3], float *x, float *y,
343                               int res, float axis_x[3], float axis_y[3])
344 {
345         float v1[3], v2[3], c[3], v3[3], v4[3], e1[3], e2[3];
346         float eps = FLT_EPSILON * 4000;
347         
348         if (len_v3(l->v->no) == 0.0f)
349                 BM_vert_normal_update_all(l->v);
350         if (len_v3(tl->v->no) == 0.0f)
351                 BM_vert_normal_update_all(tl->v);
352
353         compute_mdisp_quad(tl, v1, v2, v3, v4, e1, e2);
354
355         /* expand quad a bit */
356         cent_quad_v3(c, v1, v2, v3, v4);
357         
358         sub_v3_v3(v1, c); sub_v3_v3(v2, c);
359         sub_v3_v3(v3, c); sub_v3_v3(v4, c);
360         mul_v3_fl(v1, 1.0f + eps); mul_v3_fl(v2, 1.0f + eps);
361         mul_v3_fl(v3, 1.0f + eps); mul_v3_fl(v4, 1.0f + eps);
362         add_v3_v3(v1, c); add_v3_v3(v2, c);
363         add_v3_v3(v3, c); add_v3_v3(v4, c);
364         
365         if (!quad_co(x, y, v1, v2, v3, v4, p, l->v->no))
366                 return 0;
367         
368         *x *= res - 1;
369         *y *= res - 1;
370
371         mdisp_axis_from_quad(v1, v2, v3, v4, axis_x, axis_y);
372
373         return 1;
374 }
375
376 static float bm_loop_flip_equotion(float mat[2][2], float b[2], float target_axis_x[3], float target_axis_y[3],
377                                    float coord[3], int i, int j)
378 {
379         mat[0][0] = target_axis_x[i];
380         mat[0][1] = target_axis_y[i];
381         mat[1][0] = target_axis_x[j];
382         mat[1][1] = target_axis_y[j];
383         b[0] = coord[i];
384         b[1] = coord[j];
385
386         return mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0];
387 }
388
389 static void bm_loop_flip_disp(float source_axis_x[3], float source_axis_y[3],
390                                  float target_axis_x[3], float target_axis_y[3], float disp[3])
391 {
392         float vx[3], vy[3], coord[3];
393         float n[3], vec[3];
394         float b[2], mat[2][2], d;
395
396         mul_v3_v3fl(vx, source_axis_x, disp[0]);
397         mul_v3_v3fl(vy, source_axis_y, disp[1]);
398         add_v3_v3v3(coord, vx, vy);
399
400         /* project displacement from source grid plane onto target grid plane */
401         cross_v3_v3v3(n, target_axis_x, target_axis_y);
402         project_v3_v3v3(vec, coord, n);
403         sub_v3_v3v3(coord, coord, vec);
404
405         d = bm_loop_flip_equotion(mat, b, target_axis_x, target_axis_y, coord, 0, 1);
406
407         if (fabsf(d) < 1e-4f) {
408                 d = bm_loop_flip_equotion(mat, b, target_axis_x, target_axis_y, coord, 0, 2);
409                 if (fabsf(d) < 1e-4f)
410                         d = bm_loop_flip_equotion(mat, b, target_axis_x, target_axis_y, coord, 1, 2);
411         }
412
413         disp[0] = (b[0] * mat[1][1] - mat[0][1] * b[1]) / d;
414         disp[1] = (mat[0][0] * b[1] - b[0] * mat[1][0]) / d;
415 }
416
417 static void bm_loop_interp_mdisps(BMesh *bm, BMLoop *target, BMFace *source)
418 {
419         MDisps *mdisps;
420         BMLoop *l_iter;
421         BMLoop *l_first;
422         float x, y, d, v1[3], v2[3], v3[3], v4[3] = {0.0f, 0.0f, 0.0f}, e1[3], e2[3];
423         int ix, iy, res;
424         float axis_x[3], axis_y[3];
425         
426         /* ignore 2-edged faces */
427         if (UNLIKELY(target->f->len < 3))
428                 return;
429         
430         if (!CustomData_has_layer(&bm->ldata, CD_MDISPS))
431                 return;
432         
433         mdisps = CustomData_bmesh_get(&bm->ldata, target->head.data, CD_MDISPS);
434         compute_mdisp_quad(target, v1, v2, v3, v4, e1, e2);
435         
436         /* if no disps data allocate a new grid, the size of the first grid in source. */
437         if (!mdisps->totdisp) {
438                 MDisps *md2 = CustomData_bmesh_get(&bm->ldata, BM_FACE_FIRST_LOOP(source)->head.data, CD_MDISPS);
439                 
440                 mdisps->totdisp = md2->totdisp;
441                 mdisps->level = md2->level;
442                 if (mdisps->totdisp) {
443                         mdisps->disps = MEM_callocN(sizeof(float) * 3 * mdisps->totdisp,
444                                                     "mdisp->disps in bmesh_loop_intern_mdisps");
445                 }
446                 else {
447                         return;
448                 }
449         }
450         
451         mdisp_axis_from_quad(v1, v2, v3, v4, axis_x, axis_y);
452
453         res = (int)sqrt(mdisps->totdisp);
454         d = 1.0f / (float)(res - 1);
455         for (x = 0.0f, ix = 0; ix < res; x += d, ix++) {
456                 for (y = 0.0f, iy = 0; iy < res; y += d, iy++) {
457                         float co1[3], co2[3], co[3];
458                         
459                         copy_v3_v3(co1, e1);
460                         
461                         mul_v3_fl(co1, y);
462                         add_v3_v3(co1, v1);
463                         
464                         copy_v3_v3(co2, e2);
465                         mul_v3_fl(co2, y);
466                         add_v3_v3(co2, v4);
467                         
468                         sub_v3_v3v3(co, co2, co1);
469                         mul_v3_fl(co, x);
470                         add_v3_v3(co, co1);
471                         
472                         l_iter = l_first = BM_FACE_FIRST_LOOP(source);
473                         do {
474                                 float x2, y2;
475                                 MDisps *md1, *md2;
476                                 float src_axis_x[3], src_axis_y[3];
477
478                                 md1 = CustomData_bmesh_get(&bm->ldata, target->head.data, CD_MDISPS);
479                                 md2 = CustomData_bmesh_get(&bm->ldata, l_iter->head.data, CD_MDISPS);
480                                 
481                                 if (mdisp_in_mdispquad(target, l_iter, co, &x2, &y2, res, src_axis_x, src_axis_y)) {
482                                         old_mdisps_bilinear(md1->disps[iy * res + ix], md2->disps, res, (float)x2, (float)y2);
483                                         bm_loop_flip_disp(src_axis_x, src_axis_y, axis_x, axis_y, md1->disps[iy * res + ix]);
484
485                                         break;
486                                 }
487                         } while ((l_iter = l_iter->next) != l_first);
488                 }
489         }
490 }
491
492 /**
493  * smooths boundaries between multires grids,
494  * including some borders in adjacent faces
495  */
496 void BM_face_multires_bounds_smooth(BMesh *bm, BMFace *f)
497 {
498         BMLoop *l;
499         BMIter liter;
500         
501         if (!CustomData_has_layer(&bm->ldata, CD_MDISPS))
502                 return;
503         
504         BM_ITER_ELEM (l, &liter, f, BM_LOOPS_OF_FACE) {
505                 MDisps *mdp = CustomData_bmesh_get(&bm->ldata, l->prev->head.data, CD_MDISPS);
506                 MDisps *mdl = CustomData_bmesh_get(&bm->ldata, l->head.data, CD_MDISPS);
507                 MDisps *mdn = CustomData_bmesh_get(&bm->ldata, l->next->head.data, CD_MDISPS);
508                 float co1[3];
509                 int sides;
510                 int y;
511                 
512                 /*
513                  *  mdisps is a grid of displacements, ordered thus:
514                  *
515                  *                     v4/next
516                  *                       |
517                  *   |      v1/cent-----mid2 ---> x
518                  *   |         |         |
519                  *   |         |         |
520                  *  v2/prev---mid1-----v3/cur
521                  *             |
522                  *             V
523                  *             y
524                  */
525
526                 sides = (int)sqrt(mdp->totdisp);
527                 for (y = 0; y < sides; y++) {
528                         add_v3_v3v3(co1, mdn->disps[y * sides], mdl->disps[y]);
529                         mul_v3_fl(co1, 0.5);
530
531                         copy_v3_v3(mdn->disps[y * sides], co1);
532                         copy_v3_v3(mdl->disps[y], co1);
533                 }
534         }
535         
536         BM_ITER_ELEM (l, &liter, f, BM_LOOPS_OF_FACE) {
537                 MDisps *mdl1 = CustomData_bmesh_get(&bm->ldata, l->head.data, CD_MDISPS);
538                 MDisps *mdl2;
539                 float co1[3], co2[3], co[3];
540                 int sides;
541                 int y;
542                 
543                 /*
544                  *  mdisps is a grid of displacements, ordered thus:
545                  *
546                  *                     v4/next
547                  *                       |
548                  *   |      v1/cent-----mid2 ---> x
549                  *   |         |         |
550                  *   |         |         |
551                  *  v2/prev---mid1-----v3/cur
552                  *             |
553                  *             V
554                  *             y
555                  */
556
557                 if (l->radial_next == l)
558                         continue;
559
560                 if (l->radial_next->v == l->v)
561                         mdl2 = CustomData_bmesh_get(&bm->ldata, l->radial_next->head.data, CD_MDISPS);
562                 else
563                         mdl2 = CustomData_bmesh_get(&bm->ldata, l->radial_next->next->head.data, CD_MDISPS);
564
565                 sides = (int)sqrt(mdl1->totdisp);
566                 for (y = 0; y < sides; y++) {
567                         int a1, a2, o1, o2;
568                         
569                         if (l->v != l->radial_next->v) {
570                                 a1 = sides * y + sides - 2;
571                                 a2 = (sides - 2) * sides + y;
572                                 
573                                 o1 = sides * y + sides - 1;
574                                 o2 = (sides - 1) * sides + y;
575                         }
576                         else {
577                                 a1 = sides * y + sides - 2;
578                                 a2 = sides * y + sides - 2;
579                                 o1 = sides * y + sides - 1;
580                                 o2 = sides * y + sides - 1;
581                         }
582                         
583                         /* magic blending numbers, hardcoded! */
584                         add_v3_v3v3(co1, mdl1->disps[a1], mdl2->disps[a2]);
585                         mul_v3_fl(co1, 0.18);
586                         
587                         add_v3_v3v3(co2, mdl1->disps[o1], mdl2->disps[o2]);
588                         mul_v3_fl(co2, 0.32);
589                         
590                         add_v3_v3v3(co, co1, co2);
591                         
592                         copy_v3_v3(mdl1->disps[o1], co);
593                         copy_v3_v3(mdl2->disps[o2], co);
594                 }
595         }
596 }
597
598 /**
599  * project the multires grid in target onto source's set of multires grids
600  */
601 void BM_loop_interp_multires(BMesh *bm, BMLoop *target, BMFace *source)
602 {
603         bm_loop_interp_mdisps(bm, target, source);
604 }
605
606 /**
607  * projects a single loop, target, onto source for customdata interpolation. multires is handled.
608  * if do_vertex is true, target's vert data will also get interpolated.
609  */
610 void BM_loop_interp_from_face(BMesh *bm, BMLoop *target, BMFace *source,
611                               int do_vertex, int do_multires)
612 {
613         BMLoop *l_iter;
614         BMLoop *l_first;
615         void **blocks = NULL;
616         void **vblocks = NULL;
617         float (*cos)[3] = NULL, co[3], *w = NULL;
618         float cent[3] = {0.0f, 0.0f, 0.0f};
619         BLI_array_fixedstack_declare(cos,      BM_NGON_STACK_SIZE, source->len, __func__);
620         BLI_array_fixedstack_declare(w,        BM_NGON_STACK_SIZE, source->len, __func__);
621         BLI_array_fixedstack_declare(blocks,   BM_NGON_STACK_SIZE, source->len, __func__);
622         BLI_array_fixedstack_declare(vblocks,  BM_NGON_STACK_SIZE, do_vertex ? source->len : 0, __func__);
623         int i, ax, ay;
624
625         BM_elem_attrs_copy(bm, bm, source, target->f);
626
627         i = 0;
628         l_iter = l_first = BM_FACE_FIRST_LOOP(source);
629         do {
630                 copy_v3_v3(cos[i], l_iter->v->co);
631                 add_v3_v3(cent, cos[i]);
632
633                 w[i] = 0.0f;
634                 blocks[i] = l_iter->head.data;
635
636                 if (do_vertex) {
637                         vblocks[i] = l_iter->v->head.data;
638                 }
639                 i++;
640
641         } while ((l_iter = l_iter->next) != l_first);
642
643         /* find best projection of face XY, XZ or YZ: barycentric weights of
644          * the 2d projected coords are the same and faster to compute */
645
646         axis_dominant_v3(&ax, &ay, source->no);
647
648         /* scale source face coordinates a bit, so points sitting directly on an
649          * edge will work. */
650         mul_v3_fl(cent, 1.0f / (float)source->len);
651         for (i = 0; i < source->len; i++) {
652                 float vec[3], tmp[3];
653                 sub_v3_v3v3(vec, cent, cos[i]);
654                 mul_v3_fl(vec, 0.001f);
655                 add_v3_v3(cos[i], vec);
656
657                 copy_v3_v3(tmp, cos[i]);
658                 cos[i][0] = tmp[ax];
659                 cos[i][1] = tmp[ay];
660                 cos[i][2] = 0.0f;
661         }
662
663
664         /* interpolate */
665         co[0] = target->v->co[ax];
666         co[1] = target->v->co[ay];
667         co[2] = 0.0f;
668
669         interp_weights_poly_v3(w, cos, source->len, co);
670         CustomData_bmesh_interp(&bm->ldata, blocks, w, NULL, source->len, target->head.data);
671         if (do_vertex) {
672                 CustomData_bmesh_interp(&bm->vdata, vblocks, w, NULL, source->len, target->v->head.data);
673                 BLI_array_fixedstack_free(vblocks);
674         }
675
676         BLI_array_fixedstack_free(cos);
677         BLI_array_fixedstack_free(w);
678         BLI_array_fixedstack_free(blocks);
679
680         if (do_multires) {
681                 if (CustomData_has_layer(&bm->ldata, CD_MDISPS)) {
682                         bm_loop_interp_mdisps(bm, target, source);
683                 }
684         }
685 }
686
687
688 void BM_vert_interp_from_face(BMesh *bm, BMVert *v, BMFace *source)
689 {
690         BMLoop *l_iter;
691         BMLoop *l_first;
692         void **blocks = NULL;
693         float (*cos)[3] = NULL, *w = NULL;
694         float cent[3] = {0.0f, 0.0f, 0.0f};
695         BLI_array_fixedstack_declare(cos,      BM_NGON_STACK_SIZE, source->len, __func__);
696         BLI_array_fixedstack_declare(w,        BM_NGON_STACK_SIZE, source->len, __func__);
697         BLI_array_fixedstack_declare(blocks,   BM_NGON_STACK_SIZE, source->len, __func__);
698         int i;
699
700         i = 0;
701         l_iter = l_first = BM_FACE_FIRST_LOOP(source);
702         do {
703                 copy_v3_v3(cos[i], l_iter->v->co);
704                 add_v3_v3(cent, cos[i]);
705
706                 w[i] = 0.0f;
707                 blocks[i] = l_iter->v->head.data;
708                 i++;
709         } while ((l_iter = l_iter->next) != l_first);
710
711         /* scale source face coordinates a bit, so points sitting directly on an
712          * edge will work. */
713         mul_v3_fl(cent, 1.0f / (float)source->len);
714         for (i = 0; i < source->len; i++) {
715                 float vec[3];
716                 sub_v3_v3v3(vec, cent, cos[i]);
717                 mul_v3_fl(vec, 0.01f);
718                 add_v3_v3(cos[i], vec);
719         }
720
721         /* interpolate */
722         interp_weights_poly_v3(w, cos, source->len, v->co);
723         CustomData_bmesh_interp(&bm->vdata, blocks, w, NULL, source->len, v->head.data);
724
725         BLI_array_fixedstack_free(cos);
726         BLI_array_fixedstack_free(w);
727         BLI_array_fixedstack_free(blocks);
728 }
729
730 static void update_data_blocks(BMesh *bm, CustomData *olddata, CustomData *data)
731 {
732         BMIter iter;
733         BLI_mempool *oldpool = olddata->pool;
734         void *block;
735
736         if (data == &bm->vdata) {
737                 BMVert *eve;
738
739                 CustomData_bmesh_init_pool(data, bm->totvert, BM_VERT);
740
741                 BM_ITER_MESH (eve, &iter, bm, BM_VERTS_OF_MESH) {
742                         block = NULL;
743                         CustomData_bmesh_set_default(data, &block);
744                         CustomData_bmesh_copy_data(olddata, data, eve->head.data, &block);
745                         CustomData_bmesh_free_block(olddata, &eve->head.data);
746                         eve->head.data = block;
747                 }
748         }
749         else if (data == &bm->edata) {
750                 BMEdge *eed;
751
752                 CustomData_bmesh_init_pool(data, bm->totedge, BM_EDGE);
753
754                 BM_ITER_MESH (eed, &iter, bm, BM_EDGES_OF_MESH) {
755                         block = NULL;
756                         CustomData_bmesh_set_default(data, &block);
757                         CustomData_bmesh_copy_data(olddata, data, eed->head.data, &block);
758                         CustomData_bmesh_free_block(olddata, &eed->head.data);
759                         eed->head.data = block;
760                 }
761         }
762         else if (data == &bm->ldata) {
763                 BMIter liter;
764                 BMFace *efa;
765                 BMLoop *l;
766
767                 CustomData_bmesh_init_pool(data, bm->totloop, BM_LOOP);
768                 BM_ITER_MESH (efa, &iter, bm, BM_FACES_OF_MESH) {
769                         BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
770                                 block = NULL;
771                                 CustomData_bmesh_set_default(data, &block);
772                                 CustomData_bmesh_copy_data(olddata, data, l->head.data, &block);
773                                 CustomData_bmesh_free_block(olddata, &l->head.data);
774                                 l->head.data = block;
775                         }
776                 }
777         }
778         else if (data == &bm->pdata) {
779                 BMFace *efa;
780
781                 CustomData_bmesh_init_pool(data, bm->totface, BM_FACE);
782
783                 BM_ITER_MESH (efa, &iter, bm, BM_FACES_OF_MESH) {
784                         block = NULL;
785                         CustomData_bmesh_set_default(data, &block);
786                         CustomData_bmesh_copy_data(olddata, data, efa->head.data, &block);
787                         CustomData_bmesh_free_block(olddata, &efa->head.data);
788                         efa->head.data = block;
789                 }
790         }
791         else {
792                 /* should never reach this! */
793                 BLI_assert(0);
794         }
795
796         if (oldpool) {
797                 /* this should never happen but can when dissolve fails - [#28960] */
798                 BLI_assert(data->pool != oldpool);
799
800                 BLI_mempool_destroy(oldpool);
801         }
802 }
803
804 void BM_data_layer_add(BMesh *bm, CustomData *data, int type)
805 {
806         CustomData olddata;
807
808         olddata = *data;
809         olddata.layers = (olddata.layers) ? MEM_dupallocN(olddata.layers): NULL;
810
811         /* the pool is now owned by olddata and must not be shared */
812         data->pool = NULL;
813
814         CustomData_add_layer(data, type, CD_DEFAULT, NULL, 0);
815
816         update_data_blocks(bm, &olddata, data);
817         if (olddata.layers) MEM_freeN(olddata.layers);
818 }
819
820 void BM_data_layer_add_named(BMesh *bm, CustomData *data, int type, const char *name)
821 {
822         CustomData olddata;
823
824         olddata = *data;
825         olddata.layers = (olddata.layers) ? MEM_dupallocN(olddata.layers): NULL;
826
827         /* the pool is now owned by olddata and must not be shared */
828         data->pool = NULL;
829
830         CustomData_add_layer_named(data, type, CD_DEFAULT, NULL, 0, name);
831
832         update_data_blocks(bm, &olddata, data);
833         if (olddata.layers) MEM_freeN(olddata.layers);
834 }
835
836 void BM_data_layer_free(BMesh *bm, CustomData *data, int type)
837 {
838         CustomData olddata;
839
840         olddata = *data;
841         olddata.layers = (olddata.layers) ? MEM_dupallocN(olddata.layers): NULL;
842
843         /* the pool is now owned by olddata and must not be shared */
844         data->pool = NULL;
845
846         CustomData_free_layer_active(data, type, 0);
847
848         update_data_blocks(bm, &olddata, data);
849         if (olddata.layers) MEM_freeN(olddata.layers);
850 }
851
852 void BM_data_layer_free_n(BMesh *bm, CustomData *data, int type, int n)
853 {
854         CustomData olddata;
855
856         olddata = *data;
857         olddata.layers = (olddata.layers) ? MEM_dupallocN(olddata.layers): NULL;
858
859         /* the pool is now owned by olddata and must not be shared */
860         data->pool = NULL;
861
862         CustomData_free_layer(data, type, 0, CustomData_get_layer_index_n(data, type, n));
863         
864         update_data_blocks(bm, &olddata, data);
865         if (olddata.layers) MEM_freeN(olddata.layers);
866 }
867
868 void BM_data_layer_copy(BMesh *bm, CustomData *data, int type, int src_n, int dst_n)
869 {
870         BMIter iter;
871
872         if (&bm->vdata == data) {
873                 BMVert *eve;
874
875                 BM_ITER_MESH (eve, &iter, bm, BM_VERTS_OF_MESH) {
876                         void *ptr = CustomData_bmesh_get_n(data, eve->head.data, type, src_n);
877                         CustomData_bmesh_set_n(data, eve->head.data, type, dst_n, ptr);
878                 }
879         }
880         else if (&bm->edata == data) {
881                 BMEdge *eed;
882
883                 BM_ITER_MESH (eed, &iter, bm, BM_EDGES_OF_MESH) {
884                         void *ptr = CustomData_bmesh_get_n(data, eed->head.data, type, src_n);
885                         CustomData_bmesh_set_n(data, eed->head.data, type, dst_n, ptr);
886                 }
887         }
888         else if (&bm->pdata == data) {
889                 BMFace *efa;
890
891                 BM_ITER_MESH (efa, &iter, bm, BM_FACES_OF_MESH) {
892                         void *ptr = CustomData_bmesh_get_n(data, efa->head.data, type, src_n);
893                         CustomData_bmesh_set_n(data, efa->head.data, type, dst_n, ptr);
894                 }
895         }
896         else if (&bm->ldata == data) {
897                 BMIter liter;
898                 BMFace *efa;
899                 BMLoop *l;
900
901                 BM_ITER_MESH (efa, &iter, bm, BM_FACES_OF_MESH) {
902                         BM_ITER_ELEM (l, &liter, efa, BM_LOOPS_OF_FACE) {
903                                 void *ptr = CustomData_bmesh_get_n(data, l->head.data, type, src_n);
904                                 CustomData_bmesh_set_n(data, l->head.data, type, dst_n, ptr);
905                         }
906                 }
907         }
908         else {
909                 /* should never reach this! */
910                 BLI_assert(0);
911         }
912 }
913
914 float BM_elem_float_data_get(CustomData *cd, void *element, int type)
915 {
916         float *f = CustomData_bmesh_get(cd, ((BMHeader *)element)->data, type);
917         return f ? *f : 0.0f;
918 }
919
920 void BM_elem_float_data_set(CustomData *cd, void *element, int type, const float val)
921 {
922         float *f = CustomData_bmesh_get(cd, ((BMHeader *)element)->data, type);
923         if (f) *f = val;
924 }