Cleanup: warnings & style
[blender.git] / source / blender / blenkernel / intern / mesh_evaluate.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) 2001-2002 by NaN Holding BV.
19  * All rights reserved.
20  *
21  * Contributor(s): Blender Foundation
22  *
23  * ***** END GPL LICENSE BLOCK *****
24  */
25
26 /** \file blender/blenkernel/intern/mesh_evaluate.c
27  *  \ingroup bke
28  *
29  * Functions to evaluate mesh data.
30  */
31
32 #include <limits.h>
33
34 #include "MEM_guardedalloc.h"
35
36 #include "DNA_object_types.h"
37 #include "DNA_mesh_types.h"
38 #include "DNA_meshdata_types.h"
39
40 #include "BLI_utildefines.h"
41 #include "BLI_memarena.h"
42 #include "BLI_mempool.h"
43 #include "BLI_math.h"
44 #include "BLI_edgehash.h"
45 #include "BLI_bitmap.h"
46 #include "BLI_polyfill2d.h"
47 #include "BLI_linklist.h"
48 #include "BLI_linklist_stack.h"
49 #include "BLI_alloca.h"
50 #include "BLI_stack.h"
51 #include "BLI_task.h"
52
53 #include "BKE_customdata.h"
54 #include "BKE_global.h"
55 #include "BKE_mesh.h"
56 #include "BKE_multires.h"
57 #include "BKE_report.h"
58
59 #include "BLI_strict_flags.h"
60
61 #include "mikktspace.h"
62
63 // #define DEBUG_TIME
64
65 #include "PIL_time.h"
66 #ifdef DEBUG_TIME
67 #  include "PIL_time_utildefines.h"
68 #endif
69
70 /* -------------------------------------------------------------------- */
71
72 /** \name Mesh Normal Calculation
73  * \{ */
74
75 /**
76  * Call when there are no polygons.
77  */
78 static void mesh_calc_normals_vert_fallback(MVert *mverts, int numVerts)
79 {
80         int i;
81         for (i = 0; i < numVerts; i++) {
82                 MVert *mv = &mverts[i];
83                 float no[3];
84
85                 normalize_v3_v3(no, mv->co);
86                 normal_float_to_short_v3(mv->no, no);
87         }
88 }
89
90 /* Calculate vertex and face normals, face normals are returned in *r_faceNors if non-NULL
91  * and vertex normals are stored in actual mverts.
92  */
93 void BKE_mesh_calc_normals_mapping(
94         MVert *mverts, int numVerts,
95         const MLoop *mloop, const MPoly *mpolys, int numLoops, int numPolys, float (*r_polyNors)[3],
96         const MFace *mfaces, int numFaces, const int *origIndexFace, float (*r_faceNors)[3])
97 {
98         BKE_mesh_calc_normals_mapping_ex(
99                 mverts, numVerts, mloop, mpolys,
100                 numLoops, numPolys, r_polyNors, mfaces, numFaces,
101                 origIndexFace, r_faceNors, false);
102 }
103 /* extended version of 'BKE_mesh_calc_normals_poly' with option not to calc vertex normals */
104 void BKE_mesh_calc_normals_mapping_ex(
105         MVert *mverts, int numVerts,
106         const MLoop *mloop, const MPoly *mpolys,
107         int numLoops, int numPolys, float (*r_polyNors)[3],
108         const MFace *mfaces, int numFaces, const int *origIndexFace, float (*r_faceNors)[3],
109         const bool only_face_normals)
110 {
111         float (*pnors)[3] = r_polyNors, (*fnors)[3] = r_faceNors;
112         int i;
113         const MFace *mf;
114         const MPoly *mp;
115
116         if (numPolys == 0) {
117                 if (only_face_normals == false) {
118                         mesh_calc_normals_vert_fallback(mverts, numVerts);
119                 }
120                 return;
121         }
122
123         /* if we are not calculating verts and no verts were passes then we have nothing to do */
124         if ((only_face_normals == true) && (r_polyNors == NULL) && (r_faceNors == NULL)) {
125                 printf("%s: called with nothing to do\n", __func__);
126                 return;
127         }
128
129         if (!pnors) pnors = MEM_callocN(sizeof(float[3]) * (size_t)numPolys, __func__);
130         /* if (!fnors) fnors = MEM_callocN(sizeof(float[3]) * numFaces, "face nors mesh.c"); */ /* NO NEED TO ALLOC YET */
131
132
133         if (only_face_normals == false) {
134                 /* vertex normals are optional, they require some extra calculations,
135                  * so make them optional */
136                 BKE_mesh_calc_normals_poly(mverts, NULL, numVerts, mloop, mpolys, numLoops, numPolys, pnors, false);
137         }
138         else {
139                 /* only calc poly normals */
140                 mp = mpolys;
141                 for (i = 0; i < numPolys; i++, mp++) {
142                         BKE_mesh_calc_poly_normal(mp, mloop + mp->loopstart, mverts, pnors[i]);
143                 }
144         }
145
146         if (origIndexFace &&
147             /* fnors == r_faceNors */ /* NO NEED TO ALLOC YET */
148             fnors != NULL &&
149             numFaces)
150         {
151                 mf = mfaces;
152                 for (i = 0; i < numFaces; i++, mf++, origIndexFace++) {
153                         if (*origIndexFace < numPolys) {
154                                 copy_v3_v3(fnors[i], pnors[*origIndexFace]);
155                         }
156                         else {
157                                 /* eek, we're not corresponding to polys */
158                                 printf("error in %s: tessellation face indices are incorrect.  normals may look bad.\n", __func__);
159                         }
160                 }
161         }
162
163         if (pnors != r_polyNors) MEM_freeN(pnors);
164         /* if (fnors != r_faceNors) MEM_freeN(fnors); */ /* NO NEED TO ALLOC YET */
165
166         fnors = pnors = NULL;
167         
168 }
169
170 typedef struct MeshCalcNormalsData {
171         const MPoly *mpolys;
172         const MLoop *mloop;
173         MVert *mverts;
174         float (*pnors)[3];
175         float (*vnors)[3];
176 } MeshCalcNormalsData;
177
178 static void mesh_calc_normals_poly_task_cb(void *userdata, const int pidx)
179 {
180         MeshCalcNormalsData *data = userdata;
181         const MPoly *mp = &data->mpolys[pidx];
182
183         BKE_mesh_calc_poly_normal(mp, data->mloop + mp->loopstart, data->mverts, data->pnors[pidx]);
184 }
185
186 static void mesh_calc_normals_poly_accum_task_cb(void *userdata, const int pidx)
187 {
188         MeshCalcNormalsData *data = userdata;
189         const MPoly *mp = &data->mpolys[pidx];
190         const MLoop *ml = &data->mloop[mp->loopstart];
191         const MVert *mverts = data->mverts;
192
193         float pnor_temp[3];
194         float *pnor = data->pnors ? data->pnors[pidx] : pnor_temp;
195         float (*vnors)[3] = data->vnors;
196
197         const int nverts = mp->totloop;
198         float (*edgevecbuf)[3] = BLI_array_alloca(edgevecbuf, (size_t)nverts);
199         int i;
200
201         /* Polygon Normal and edge-vector */
202         /* inline version of #BKE_mesh_calc_poly_normal, also does edge-vectors */
203         {
204                 int i_prev = nverts - 1;
205                 const float *v_prev = mverts[ml[i_prev].v].co;
206                 const float *v_curr;
207
208                 zero_v3(pnor);
209                 /* Newell's Method */
210                 for (i = 0; i < nverts; i++) {
211                         v_curr = mverts[ml[i].v].co;
212                         add_newell_cross_v3_v3v3(pnor, v_prev, v_curr);
213
214                         /* Unrelated to normalize, calculate edge-vector */
215                         sub_v3_v3v3(edgevecbuf[i_prev], v_prev, v_curr);
216                         normalize_v3(edgevecbuf[i_prev]);
217                         i_prev = i;
218
219                         v_prev = v_curr;
220                 }
221                 if (UNLIKELY(normalize_v3(pnor) == 0.0f)) {
222                         pnor[2] = 1.0f; /* other axis set to 0.0 */
223                 }
224         }
225
226         /* accumulate angle weighted face normal */
227         /* inline version of #accumulate_vertex_normals_poly */
228         {
229                 const float *prev_edge = edgevecbuf[nverts - 1];
230
231                 for (i = 0; i < nverts; i++) {
232                         const float *cur_edge = edgevecbuf[i];
233
234                         /* calculate angle between the two poly edges incident on
235                          * this vertex */
236                         const float fac = saacos(-dot_v3v3(cur_edge, prev_edge));
237
238                         /* accumulate */
239                         madd_v3_v3fl(vnors[ml[i].v], pnor, fac);
240                         prev_edge = cur_edge;
241                 }
242         }
243
244 }
245
246 void BKE_mesh_calc_normals_poly(
247         MVert *mverts, float (*r_vertnors)[3], int numVerts,
248         const MLoop *mloop, const MPoly *mpolys,
249         int UNUSED(numLoops), int numPolys, float (*r_polynors)[3],
250         const bool only_face_normals)
251 {
252         float (*pnors)[3] = r_polynors;
253         float (*vnors)[3] = r_vertnors;
254         bool free_vnors = false;
255         int i;
256
257         if (only_face_normals) {
258                 BLI_assert((pnors != NULL) || (numPolys == 0));
259                 BLI_assert(r_vertnors == NULL);
260
261                 MeshCalcNormalsData data = {
262                     .mpolys = mpolys, .mloop = mloop, .mverts = mverts, .pnors = pnors,
263                 };
264
265                 BLI_task_parallel_range(0, numPolys, &data, mesh_calc_normals_poly_task_cb, (numPolys > BKE_MESH_OMP_LIMIT));
266                 return;
267         }
268
269         /* first go through and calculate normals for all the polys */
270         if (vnors == NULL) {
271                 vnors = MEM_callocN(sizeof(*vnors) * (size_t)numVerts, __func__);
272                 free_vnors = true;
273         }
274         else {
275                 memset(vnors, 0, sizeof(*vnors) * (size_t)numVerts);
276         }
277
278         MeshCalcNormalsData data = {
279             .mpolys = mpolys, .mloop = mloop, .mverts = mverts, .pnors = pnors, .vnors = vnors,
280         };
281
282         BLI_task_parallel_range(0, numPolys, &data, mesh_calc_normals_poly_accum_task_cb, (numPolys > BKE_MESH_OMP_LIMIT));
283
284         for (i = 0; i < numVerts; i++) {
285                 MVert *mv = &mverts[i];
286                 float *no = vnors[i];
287
288                 if (UNLIKELY(normalize_v3(no) == 0.0f)) {
289                         /* following Mesh convention; we use vertex coordinate itself for normal in this case */
290                         normalize_v3_v3(no, mv->co);
291                 }
292
293                 normal_float_to_short_v3(mv->no, no);
294         }
295
296         if (free_vnors) {
297                 MEM_freeN(vnors);
298         }
299 }
300
301 void BKE_mesh_calc_normals(Mesh *mesh)
302 {
303 #ifdef DEBUG_TIME
304         TIMEIT_START(BKE_mesh_calc_normals);
305 #endif
306         BKE_mesh_calc_normals_poly(mesh->mvert, NULL, mesh->totvert,
307                                    mesh->mloop, mesh->mpoly, mesh->totloop, mesh->totpoly,
308                                    NULL, false);
309 #ifdef DEBUG_TIME
310         TIMEIT_END(BKE_mesh_calc_normals);
311 #endif
312 }
313
314 void BKE_mesh_calc_normals_tessface(
315         MVert *mverts, int numVerts,
316         const MFace *mfaces, int numFaces,
317         float (*r_faceNors)[3])
318 {
319         float (*tnorms)[3] = MEM_callocN(sizeof(*tnorms) * (size_t)numVerts, "tnorms");
320         float (*fnors)[3] = (r_faceNors) ? r_faceNors : MEM_callocN(sizeof(*fnors) * (size_t)numFaces, "meshnormals");
321         int i;
322
323         for (i = 0; i < numFaces; i++) {
324                 const MFace *mf = &mfaces[i];
325                 float *f_no = fnors[i];
326                 float *n4 = (mf->v4) ? tnorms[mf->v4] : NULL;
327                 const float *c4 = (mf->v4) ? mverts[mf->v4].co : NULL;
328
329                 if (mf->v4)
330                         normal_quad_v3(f_no, mverts[mf->v1].co, mverts[mf->v2].co, mverts[mf->v3].co, mverts[mf->v4].co);
331                 else
332                         normal_tri_v3(f_no, mverts[mf->v1].co, mverts[mf->v2].co, mverts[mf->v3].co);
333
334                 accumulate_vertex_normals(tnorms[mf->v1], tnorms[mf->v2], tnorms[mf->v3], n4,
335                                           f_no, mverts[mf->v1].co, mverts[mf->v2].co, mverts[mf->v3].co, c4);
336         }
337
338         /* following Mesh convention; we use vertex coordinate itself for normal in this case */
339         for (i = 0; i < numVerts; i++) {
340                 MVert *mv = &mverts[i];
341                 float *no = tnorms[i];
342                 
343                 if (UNLIKELY(normalize_v3(no) == 0.0f)) {
344                         normalize_v3_v3(no, mv->co);
345                 }
346
347                 normal_float_to_short_v3(mv->no, no);
348         }
349         
350         MEM_freeN(tnorms);
351
352         if (fnors != r_faceNors)
353                 MEM_freeN(fnors);
354 }
355
356 void BKE_mesh_calc_normals_looptri(
357         MVert *mverts, int numVerts,
358         const MLoop *mloop,
359         const MLoopTri *looptri, int looptri_num,
360         float (*r_tri_nors)[3])
361 {
362         float (*tnorms)[3] = MEM_callocN(sizeof(*tnorms) * (size_t)numVerts, "tnorms");
363         float (*fnors)[3] = (r_tri_nors) ? r_tri_nors : MEM_callocN(sizeof(*fnors) * (size_t)looptri_num, "meshnormals");
364         int i;
365
366         for (i = 0; i < looptri_num; i++) {
367                 const MLoopTri *lt = &looptri[i];
368                 float *f_no = fnors[i];
369                 const unsigned int vtri[3] = {
370                     mloop[lt->tri[0]].v,
371                     mloop[lt->tri[1]].v,
372                     mloop[lt->tri[2]].v,
373                 };
374
375                 normal_tri_v3(
376                         f_no,
377                         mverts[vtri[0]].co, mverts[vtri[1]].co, mverts[vtri[2]].co);
378
379                 accumulate_vertex_normals_tri(
380                         tnorms[vtri[0]], tnorms[vtri[1]], tnorms[vtri[2]],
381                         f_no, mverts[vtri[0]].co, mverts[vtri[1]].co, mverts[vtri[2]].co);
382         }
383
384         /* following Mesh convention; we use vertex coordinate itself for normal in this case */
385         for (i = 0; i < numVerts; i++) {
386                 MVert *mv = &mverts[i];
387                 float *no = tnorms[i];
388
389                 if (UNLIKELY(normalize_v3(no) == 0.0f)) {
390                         normalize_v3_v3(no, mv->co);
391                 }
392
393                 normal_float_to_short_v3(mv->no, no);
394         }
395
396         MEM_freeN(tnorms);
397
398         if (fnors != r_tri_nors)
399                 MEM_freeN(fnors);
400 }
401
402 void BKE_lnor_spacearr_init(MLoopNorSpaceArray *lnors_spacearr, const int numLoops)
403 {
404         if (!(lnors_spacearr->lspacearr && lnors_spacearr->loops_pool)) {
405                 MemArena *mem;
406
407                 if (!lnors_spacearr->mem) {
408                         lnors_spacearr->mem = BLI_memarena_new(BLI_MEMARENA_STD_BUFSIZE, __func__);
409                 }
410                 mem = lnors_spacearr->mem;
411                 lnors_spacearr->lspacearr = BLI_memarena_calloc(mem, sizeof(MLoopNorSpace *) * (size_t)numLoops);
412                 lnors_spacearr->loops_pool = BLI_memarena_alloc(mem, sizeof(LinkNode) * (size_t)numLoops);
413         }
414 }
415
416 void BKE_lnor_spacearr_clear(MLoopNorSpaceArray *lnors_spacearr)
417 {
418         BLI_memarena_clear(lnors_spacearr->mem);
419         lnors_spacearr->lspacearr = NULL;
420         lnors_spacearr->loops_pool = NULL;
421 }
422
423 void BKE_lnor_spacearr_free(MLoopNorSpaceArray *lnors_spacearr)
424 {
425         BLI_memarena_free(lnors_spacearr->mem);
426         lnors_spacearr->lspacearr = NULL;
427         lnors_spacearr->loops_pool = NULL;
428         lnors_spacearr->mem = NULL;
429 }
430
431 MLoopNorSpace *BKE_lnor_space_create(MLoopNorSpaceArray *lnors_spacearr)
432 {
433         return BLI_memarena_calloc(lnors_spacearr->mem, sizeof(MLoopNorSpace));
434 }
435
436 /* This threshold is a bit touchy (usual float precision issue), this value seems OK. */
437 #define LNOR_SPACE_TRIGO_THRESHOLD (1.0f - 1e-6f)
438
439 /* Should only be called once.
440  * Beware, this modifies ref_vec and other_vec in place!
441  * In case no valid space can be generated, ref_alpha and ref_beta are set to zero (which means 'use auto lnors').
442  */
443 void BKE_lnor_space_define(MLoopNorSpace *lnor_space, const float lnor[3],
444                            float vec_ref[3], float vec_other[3], BLI_Stack *edge_vectors)
445 {
446         const float pi2 = (float)M_PI * 2.0f;
447         float tvec[3], dtp;
448         const float dtp_ref = dot_v3v3(vec_ref, lnor);
449         const float dtp_other = dot_v3v3(vec_other, lnor);
450
451         if (UNLIKELY(fabsf(dtp_ref) >= LNOR_SPACE_TRIGO_THRESHOLD || fabsf(dtp_other) >= LNOR_SPACE_TRIGO_THRESHOLD)) {
452                 /* If vec_ref or vec_other are too much aligned with lnor, we can't build lnor space,
453                  * tag it as invalid and abort. */
454                 lnor_space->ref_alpha = lnor_space->ref_beta = 0.0f;
455
456                 if (edge_vectors) {
457                         BLI_stack_clear(edge_vectors);
458                 }
459                 return;
460         }
461
462         copy_v3_v3(lnor_space->vec_lnor, lnor);
463
464         /* Compute ref alpha, average angle of all available edge vectors to lnor. */
465         if (edge_vectors) {
466                 float alpha = 0.0f;
467                 int nbr = 0;
468                 while (!BLI_stack_is_empty(edge_vectors)) {
469                         const float *vec = BLI_stack_peek(edge_vectors);
470                         alpha += saacosf(dot_v3v3(vec, lnor));
471                         BLI_stack_discard(edge_vectors);
472                         nbr++;
473                 }
474                 /* Note: In theory, this could be 'nbr > 2', but there is one case where we only have two edges for
475                  *       two loops: a smooth vertex with only two edges and two faces (our Monkey's nose has that, e.g.). */
476                 BLI_assert(nbr >= 2);  /* This piece of code shall only be called for more than one loop... */
477                 lnor_space->ref_alpha = alpha / (float)nbr;
478         }
479         else {
480                 lnor_space->ref_alpha = (saacosf(dot_v3v3(vec_ref, lnor)) + saacosf(dot_v3v3(vec_other, lnor))) / 2.0f;
481         }
482
483         /* Project vec_ref on lnor's ortho plane. */
484         mul_v3_v3fl(tvec, lnor, dtp_ref);
485         sub_v3_v3(vec_ref, tvec);
486         normalize_v3_v3(lnor_space->vec_ref, vec_ref);
487
488         cross_v3_v3v3(tvec, lnor, lnor_space->vec_ref);
489         normalize_v3_v3(lnor_space->vec_ortho, tvec);
490
491         /* Project vec_other on lnor's ortho plane. */
492         mul_v3_v3fl(tvec, lnor, dtp_other);
493         sub_v3_v3(vec_other, tvec);
494         normalize_v3(vec_other);
495
496         /* Beta is angle between ref_vec and other_vec, around lnor. */
497         dtp = dot_v3v3(lnor_space->vec_ref, vec_other);
498         if (LIKELY(dtp < LNOR_SPACE_TRIGO_THRESHOLD)) {
499                 const float beta = saacos(dtp);
500                 lnor_space->ref_beta = (dot_v3v3(lnor_space->vec_ortho, vec_other) < 0.0f) ? pi2 - beta : beta;
501         }
502         else {
503                 lnor_space->ref_beta = pi2;
504         }
505 }
506
507 void BKE_lnor_space_add_loop(MLoopNorSpaceArray *lnors_spacearr, MLoopNorSpace *lnor_space, const int ml_index,
508                              const bool do_add_loop)
509 {
510         lnors_spacearr->lspacearr[ml_index] = lnor_space;
511         if (do_add_loop) {
512                 BLI_linklist_prepend_nlink(&lnor_space->loops, SET_INT_IN_POINTER(ml_index), &lnors_spacearr->loops_pool[ml_index]);
513         }
514 }
515
516 MINLINE float unit_short_to_float(const short val)
517 {
518         return (float)val / (float)SHRT_MAX;
519 }
520
521 MINLINE short unit_float_to_short(const float val)
522 {
523         /* Rounding... */
524         return (short)floorf(val * (float)SHRT_MAX + 0.5f);
525 }
526
527 void BKE_lnor_space_custom_data_to_normal(MLoopNorSpace *lnor_space, const short clnor_data[2], float r_custom_lnor[3])
528 {
529         /* NOP custom normal data or invalid lnor space, return. */
530         if (clnor_data[0] == 0 || lnor_space->ref_alpha == 0.0f || lnor_space->ref_beta == 0.0f) {
531                 copy_v3_v3(r_custom_lnor, lnor_space->vec_lnor);
532                 return;
533         }
534
535         {
536                 /* TODO Check whether using sincosf() gives any noticeable benefit
537                  *      (could not even get it working under linux though)! */
538                 const float pi2 = (float)(M_PI * 2.0);
539                 const float alphafac = unit_short_to_float(clnor_data[0]);
540                 const float alpha = (alphafac > 0.0f ? lnor_space->ref_alpha : pi2 - lnor_space->ref_alpha) * alphafac;
541                 const float betafac = unit_short_to_float(clnor_data[1]);
542
543                 mul_v3_v3fl(r_custom_lnor, lnor_space->vec_lnor, cosf(alpha));
544
545                 if (betafac == 0.0f) {
546                         madd_v3_v3fl(r_custom_lnor, lnor_space->vec_ref, sinf(alpha));
547                 }
548                 else {
549                         const float sinalpha = sinf(alpha);
550                         const float beta = (betafac > 0.0f ? lnor_space->ref_beta : pi2 - lnor_space->ref_beta) * betafac;
551                         madd_v3_v3fl(r_custom_lnor, lnor_space->vec_ref, sinalpha * cosf(beta));
552                         madd_v3_v3fl(r_custom_lnor, lnor_space->vec_ortho, sinalpha * sinf(beta));
553                 }
554         }
555 }
556
557 void BKE_lnor_space_custom_normal_to_data(MLoopNorSpace *lnor_space, const float custom_lnor[3], short r_clnor_data[2])
558 {
559         /* We use null vector as NOP custom normal (can be simpler than giving autocomputed lnor...). */
560         if (is_zero_v3(custom_lnor) || compare_v3v3(lnor_space->vec_lnor, custom_lnor, 1e-4f)) {
561                 r_clnor_data[0] = r_clnor_data[1] = 0;
562                 return;
563         }
564
565         {
566                 const float pi2 = (float)(M_PI * 2.0);
567                 const float cos_alpha = dot_v3v3(lnor_space->vec_lnor, custom_lnor);
568                 float vec[3], cos_beta;
569                 float alpha;
570
571                 alpha = saacosf(cos_alpha);
572                 if (alpha > lnor_space->ref_alpha) {
573                         /* Note we could stick to [0, pi] range here, but makes decoding more complex, not worth it. */
574                         r_clnor_data[0] = unit_float_to_short(-(pi2 - alpha) / (pi2 - lnor_space->ref_alpha));
575                 }
576                 else {
577                         r_clnor_data[0] = unit_float_to_short(alpha / lnor_space->ref_alpha);
578                 }
579
580                 /* Project custom lnor on (vec_ref, vec_ortho) plane. */
581                 mul_v3_v3fl(vec, lnor_space->vec_lnor, -cos_alpha);
582                 add_v3_v3(vec, custom_lnor);
583                 normalize_v3(vec);
584
585                 cos_beta = dot_v3v3(lnor_space->vec_ref, vec);
586
587                 if (cos_beta < LNOR_SPACE_TRIGO_THRESHOLD) {
588                         float beta = saacosf(cos_beta);
589                         if (dot_v3v3(lnor_space->vec_ortho, vec) < 0.0f) {
590                                 beta = pi2 - beta;
591                         }
592
593                         if (beta > lnor_space->ref_beta) {
594                                 r_clnor_data[1] = unit_float_to_short(-(pi2 - beta) / (pi2 - lnor_space->ref_beta));
595                         }
596                         else {
597                                 r_clnor_data[1] = unit_float_to_short(beta / lnor_space->ref_beta);
598                         }
599                 }
600                 else {
601                         r_clnor_data[1] = 0;
602                 }
603         }
604 }
605
606 #define LOOP_SPLIT_TASK_BLOCK_SIZE 1024
607
608 typedef struct LoopSplitTaskData {
609         /* Specific to each instance (each task). */
610         MLoopNorSpace *lnor_space;  /* We have to create those outside of tasks, since afaik memarena is not threadsafe. */
611         float (*lnor)[3];
612         const MLoop *ml_curr;
613         const MLoop *ml_prev;
614         int ml_curr_index;
615         int ml_prev_index;
616         const int *e2l_prev;  /* Also used a flag to switch between single or fan process! */
617         int mp_index;
618
619         /* This one is special, it's owned and managed by worker tasks, avoid to have to create it for each fan! */
620         BLI_Stack *edge_vectors;
621
622         char pad_c;
623 } LoopSplitTaskData;
624
625 typedef struct LoopSplitTaskDataCommon {
626         /* Read/write.
627          * Note we do not need to protect it, though, since two different tasks will *always* affect different
628          * elements in the arrays. */
629         MLoopNorSpaceArray *lnors_spacearr;
630         BLI_bitmap *sharp_verts;
631         float (*loopnors)[3];
632         short (*clnors_data)[2];
633
634         /* Read-only. */
635         const MVert *mverts;
636         const MEdge *medges;
637         const MLoop *mloops;
638         const MPoly *mpolys;
639         const int (*edge_to_loops)[2];
640         const int *loop_to_poly;
641         const float (*polynors)[3];
642
643         int numPolys;
644
645         /* ***** Workers communication. ***** */
646         ThreadQueue *task_queue;
647
648 } LoopSplitTaskDataCommon;
649
650 #define INDEX_UNSET INT_MIN
651 #define INDEX_INVALID -1
652 /* See comment about edge_to_loops below. */
653 #define IS_EDGE_SHARP(_e2l) (ELEM((_e2l)[1], INDEX_UNSET, INDEX_INVALID))
654
655 static void split_loop_nor_single_do(LoopSplitTaskDataCommon *common_data, LoopSplitTaskData *data)
656 {
657         MLoopNorSpaceArray *lnors_spacearr = common_data->lnors_spacearr;
658         short (*clnors_data)[2] = common_data->clnors_data;
659
660         const MVert *mverts = common_data->mverts;
661         const MEdge *medges = common_data->medges;
662         const float (*polynors)[3] = common_data->polynors;
663
664         MLoopNorSpace *lnor_space = data->lnor_space;
665         float (*lnor)[3] = data->lnor;
666         const MLoop *ml_curr = data->ml_curr;
667         const MLoop *ml_prev = data->ml_prev;
668         const int ml_curr_index = data->ml_curr_index;
669 #if 0  /* Not needed for 'single' loop. */
670         const int ml_prev_index = data->ml_prev_index;
671         const int *e2l_prev = data->e2l_prev;
672 #endif
673         const int mp_index = data->mp_index;
674
675         /* Simple case (both edges around that vertex are sharp in current polygon),
676          * this loop just takes its poly normal.
677          */
678         copy_v3_v3(*lnor, polynors[mp_index]);
679
680         /* printf("BASIC: handling loop %d / edge %d / vert %d\n", ml_curr_index, ml_curr->e, ml_curr->v); */
681
682         /* If needed, generate this (simple!) lnor space. */
683         if (lnors_spacearr) {
684                 float vec_curr[3], vec_prev[3];
685
686                 const unsigned int mv_pivot_index = ml_curr->v;  /* The vertex we are "fanning" around! */
687                 const MVert *mv_pivot = &mverts[mv_pivot_index];
688                 const MEdge *me_curr = &medges[ml_curr->e];
689                 const MVert *mv_2 = (me_curr->v1 == mv_pivot_index) ? &mverts[me_curr->v2] : &mverts[me_curr->v1];
690                 const MEdge *me_prev = &medges[ml_prev->e];
691                 const MVert *mv_3 = (me_prev->v1 == mv_pivot_index) ? &mverts[me_prev->v2] : &mverts[me_prev->v1];
692
693                 sub_v3_v3v3(vec_curr, mv_2->co, mv_pivot->co);
694                 normalize_v3(vec_curr);
695                 sub_v3_v3v3(vec_prev, mv_3->co, mv_pivot->co);
696                 normalize_v3(vec_prev);
697
698                 BKE_lnor_space_define(lnor_space, *lnor, vec_curr, vec_prev, NULL);
699                 /* We know there is only one loop in this space, no need to create a linklist in this case... */
700                 BKE_lnor_space_add_loop(lnors_spacearr, lnor_space, ml_curr_index, false);
701
702                 if (clnors_data) {
703                         BKE_lnor_space_custom_data_to_normal(lnor_space, clnors_data[ml_curr_index], *lnor);
704                 }
705         }
706 }
707
708 static void split_loop_nor_fan_do(LoopSplitTaskDataCommon *common_data, LoopSplitTaskData *data)
709 {
710         MLoopNorSpaceArray *lnors_spacearr = common_data->lnors_spacearr;
711         float (*loopnors)[3] = common_data->loopnors;
712         short (*clnors_data)[2] = common_data->clnors_data;
713
714         const MVert *mverts = common_data->mverts;
715         const MEdge *medges = common_data->medges;
716         const MLoop *mloops = common_data->mloops;
717         const MPoly *mpolys = common_data->mpolys;
718         const int (*edge_to_loops)[2] = common_data->edge_to_loops;
719         const int *loop_to_poly = common_data->loop_to_poly;
720         const float (*polynors)[3] = common_data->polynors;
721
722         MLoopNorSpace *lnor_space = data->lnor_space;
723 #if 0  /* Not needed for 'fan' loops. */
724         float (*lnor)[3] = data->lnor;
725 #endif
726         const MLoop *ml_curr = data->ml_curr;
727         const MLoop *ml_prev = data->ml_prev;
728         const int ml_curr_index = data->ml_curr_index;
729         const int ml_prev_index = data->ml_prev_index;
730         const int mp_index = data->mp_index;
731         const int *e2l_prev = data->e2l_prev;
732
733         BLI_Stack *edge_vectors = data->edge_vectors;
734
735         /* Gah... We have to fan around current vertex, until we find the other non-smooth edge,
736          * and accumulate face normals into the vertex!
737          * Note in case this vertex has only one sharp edges, this is a waste because the normal is the same as
738          * the vertex normal, but I do not see any easy way to detect that (would need to count number
739          * of sharp edges per vertex, I doubt the additional memory usage would be worth it, especially as
740          * it should not be a common case in real-life meshes anyway).
741          */
742         const unsigned int mv_pivot_index = ml_curr->v;  /* The vertex we are "fanning" around! */
743         const MVert *mv_pivot = &mverts[mv_pivot_index];
744         const MEdge *me_org = &medges[ml_curr->e];  /* ml_curr would be mlfan_prev if we needed that one */
745         const int *e2lfan_curr;
746         float vec_curr[3], vec_prev[3], vec_org[3];
747         const MLoop *mlfan_curr, *mlfan_next;
748         const MPoly *mpfan_next;
749         float lnor[3] = {0.0f, 0.0f, 0.0f};
750         /* mlfan_vert_index: the loop of our current edge might not be the loop of our current vertex! */
751         int mlfan_curr_index, mlfan_vert_index, mpfan_curr_index;
752
753         /* We validate clnors data on the fly - cheapest way to do! */
754         int clnors_avg[2] = {0, 0};
755         short (*clnor_ref)[2] = NULL;
756         int clnors_nbr = 0;
757         bool clnors_invalid = false;
758
759         /* Temp loop normal stack. */
760         BLI_SMALLSTACK_DECLARE(normal, float *);
761         /* Temp clnors stack. */
762         BLI_SMALLSTACK_DECLARE(clnors, short *);
763
764         e2lfan_curr = e2l_prev;
765         mlfan_curr = ml_prev;
766         mlfan_curr_index = ml_prev_index;
767         mlfan_vert_index = ml_curr_index;
768         mpfan_curr_index = mp_index;
769
770         BLI_assert(mlfan_curr_index >= 0);
771         BLI_assert(mlfan_vert_index >= 0);
772         BLI_assert(mpfan_curr_index >= 0);
773
774         /* Only need to compute previous edge's vector once, then we can just reuse old current one! */
775         {
776                 const MVert *mv_2 = (me_org->v1 == mv_pivot_index) ? &mverts[me_org->v2] : &mverts[me_org->v1];
777
778                 sub_v3_v3v3(vec_org, mv_2->co, mv_pivot->co);
779                 normalize_v3(vec_org);
780                 copy_v3_v3(vec_prev, vec_org);
781
782                 if (lnors_spacearr) {
783                         BLI_stack_push(edge_vectors, vec_org);
784                 }
785         }
786
787         /* printf("FAN: vert %d, start edge %d\n", mv_pivot_index, ml_curr->e); */
788
789         while (true) {
790                 const MEdge *me_curr = &medges[mlfan_curr->e];
791                 /* Compute edge vectors.
792                  * NOTE: We could pre-compute those into an array, in the first iteration, instead of computing them
793                  *       twice (or more) here. However, time gained is not worth memory and time lost,
794                  *       given the fact that this code should not be called that much in real-life meshes...
795                  */
796                 {
797                         const MVert *mv_2 = (me_curr->v1 == mv_pivot_index) ? &mverts[me_curr->v2] : &mverts[me_curr->v1];
798
799                         sub_v3_v3v3(vec_curr, mv_2->co, mv_pivot->co);
800                         normalize_v3(vec_curr);
801                 }
802
803                 /* printf("\thandling edge %d / loop %d\n", mlfan_curr->e, mlfan_curr_index); */
804
805                 {
806                         /* Code similar to accumulate_vertex_normals_poly. */
807                         /* Calculate angle between the two poly edges incident on this vertex. */
808                         const float fac = saacos(dot_v3v3(vec_curr, vec_prev));
809                         /* Accumulate */
810                         madd_v3_v3fl(lnor, polynors[mpfan_curr_index], fac);
811
812                         if (clnors_data) {
813                                 /* Accumulate all clnors, if they are not all equal we have to fix that! */
814                                 short (*clnor)[2] = &clnors_data[mlfan_vert_index];
815                                 if (clnors_nbr) {
816                                         clnors_invalid |= ((*clnor_ref)[0] != (*clnor)[0] || (*clnor_ref)[1] != (*clnor)[1]);
817                                 }
818                                 else {
819                                         clnor_ref = clnor;
820                                 }
821                                 clnors_avg[0] += (*clnor)[0];
822                                 clnors_avg[1] += (*clnor)[1];
823                                 clnors_nbr++;
824                                 /* We store here a pointer to all custom lnors processed. */
825                                 BLI_SMALLSTACK_PUSH(clnors, (short *)*clnor);
826                         }
827                 }
828
829                 /* We store here a pointer to all loop-normals processed. */
830                 BLI_SMALLSTACK_PUSH(normal, (float *)(loopnors[mlfan_vert_index]));
831
832                 if (lnors_spacearr) {
833                         /* Assign current lnor space to current 'vertex' loop. */
834                         BKE_lnor_space_add_loop(lnors_spacearr, lnor_space, mlfan_vert_index, true);
835                         if (me_curr != me_org) {
836                                 /* We store here all edges-normalized vectors processed. */
837                                 BLI_stack_push(edge_vectors, vec_curr);
838                         }
839                 }
840
841                 if (IS_EDGE_SHARP(e2lfan_curr) || (me_curr == me_org)) {
842                         /* Current edge is sharp and we have finished with this fan of faces around this vert,
843                          * or this vert is smooth, and we have completed a full turn around it.
844                          */
845                         /* printf("FAN: Finished!\n"); */
846                         break;
847                 }
848
849                 copy_v3_v3(vec_prev, vec_curr);
850
851                 /* Warning! This is rather complex!
852                  * We have to find our next edge around the vertex (fan mode).
853                  * First we find the next loop, which is either previous or next to mlfan_curr_index, depending
854                  * whether both loops using current edge are in the same direction or not, and whether
855                  * mlfan_curr_index actually uses the vertex we are fanning around!
856                  * mlfan_curr_index is the index of mlfan_next here, and mlfan_next is not the real next one
857                  * (i.e. not the future mlfan_curr)...
858                  */
859                 mlfan_curr_index = (e2lfan_curr[0] == mlfan_curr_index) ? e2lfan_curr[1] : e2lfan_curr[0];
860                 mpfan_curr_index = loop_to_poly[mlfan_curr_index];
861
862                 BLI_assert(mlfan_curr_index >= 0);
863                 BLI_assert(mpfan_curr_index >= 0);
864
865                 mlfan_next = &mloops[mlfan_curr_index];
866                 mpfan_next = &mpolys[mpfan_curr_index];
867                 if ((mlfan_curr->v == mlfan_next->v && mlfan_curr->v == mv_pivot_index) ||
868                     (mlfan_curr->v != mlfan_next->v && mlfan_curr->v != mv_pivot_index))
869                 {
870                         /* We need the previous loop, but current one is our vertex's loop. */
871                         mlfan_vert_index = mlfan_curr_index;
872                         if (--mlfan_curr_index < mpfan_next->loopstart) {
873                                 mlfan_curr_index = mpfan_next->loopstart + mpfan_next->totloop - 1;
874                         }
875                 }
876                 else {
877                         /* We need the next loop, which is also our vertex's loop. */
878                         if (++mlfan_curr_index >= mpfan_next->loopstart + mpfan_next->totloop) {
879                                 mlfan_curr_index = mpfan_next->loopstart;
880                         }
881                         mlfan_vert_index = mlfan_curr_index;
882                 }
883                 mlfan_curr = &mloops[mlfan_curr_index];
884                 /* And now we are back in sync, mlfan_curr_index is the index of mlfan_curr! Pff! */
885
886                 e2lfan_curr = edge_to_loops[mlfan_curr->e];
887         }
888
889         {
890                 float lnor_len = normalize_v3(lnor);
891
892                 /* If we are generating lnor spacearr, we can now define the one for this fan,
893                  * and optionally compute final lnor from custom data too!
894                  */
895                 if (lnors_spacearr) {
896                         if (UNLIKELY(lnor_len == 0.0f)) {
897                                 /* Use vertex normal as fallback! */
898                                 copy_v3_v3(lnor, loopnors[mlfan_vert_index]);
899                                 lnor_len = 1.0f;
900                         }
901
902                         BKE_lnor_space_define(lnor_space, lnor, vec_org, vec_curr, edge_vectors);
903
904                         if (clnors_data) {
905                                 if (clnors_invalid) {
906                                         short *clnor;
907
908                                         clnors_avg[0] /= clnors_nbr;
909                                         clnors_avg[1] /= clnors_nbr;
910                                         /* Fix/update all clnors of this fan with computed average value. */
911                                         if (G.debug & G_DEBUG) {
912                                                 printf("Invalid clnors in this fan!\n");
913                                         }
914                                         while ((clnor = BLI_SMALLSTACK_POP(clnors))) {
915                                                 //print_v2("org clnor", clnor);
916                                                 clnor[0] = (short)clnors_avg[0];
917                                                 clnor[1] = (short)clnors_avg[1];
918                                         }
919                                         //print_v2("new clnors", clnors_avg);
920                                 }
921                                 /* Extra bonus: since smallstack is local to this func, no more need to empty it at all cost! */
922
923                                 BKE_lnor_space_custom_data_to_normal(lnor_space, *clnor_ref, lnor);
924                         }
925                 }
926
927                 /* In case we get a zero normal here, just use vertex normal already set! */
928                 if (LIKELY(lnor_len != 0.0f)) {
929                         /* Copy back the final computed normal into all related loop-normals. */
930                         float *nor;
931
932                         while ((nor = BLI_SMALLSTACK_POP(normal))) {
933                                 copy_v3_v3(nor, lnor);
934                         }
935                 }
936                 /* Extra bonus: since smallstack is local to this func, no more need to empty it at all cost! */
937         }
938 }
939
940 static void loop_split_worker_do(
941         LoopSplitTaskDataCommon *common_data, LoopSplitTaskData *data, BLI_Stack *edge_vectors)
942 {
943         BLI_assert(data->ml_curr);
944         if (data->e2l_prev) {
945                 BLI_assert((edge_vectors == NULL) || BLI_stack_is_empty(edge_vectors));
946                 data->edge_vectors = edge_vectors;
947                 split_loop_nor_fan_do(common_data, data);
948         }
949         else {
950                 /* No need for edge_vectors for 'single' case! */
951                 split_loop_nor_single_do(common_data, data);
952         }
953 }
954
955 static void loop_split_worker(TaskPool * __restrict UNUSED(pool), void *taskdata, int UNUSED(threadid))
956 {
957         LoopSplitTaskDataCommon *common_data = taskdata;
958         LoopSplitTaskData *data_buff;
959
960         /* Temp edge vectors stack, only used when computing lnor spacearr. */
961         BLI_Stack *edge_vectors = common_data->lnors_spacearr ? BLI_stack_new(sizeof(float[3]), __func__) : NULL;
962
963 #ifdef DEBUG_TIME
964         TIMEIT_START(loop_split_worker);
965 #endif
966
967         while ((data_buff = BLI_thread_queue_pop(common_data->task_queue))) {
968                 LoopSplitTaskData *data = data_buff;
969                 int i;
970
971                 for (i = 0; i < LOOP_SPLIT_TASK_BLOCK_SIZE; i++, data++) {
972                         /* A NULL ml_curr is used to tag ended data! */
973                         if (data->ml_curr == NULL) {
974                                 break;
975                         }
976                         loop_split_worker_do(common_data, data, edge_vectors);
977                 }
978
979                 MEM_freeN(data_buff);
980         }
981
982         if (edge_vectors) {
983                 BLI_stack_free(edge_vectors);
984         }
985
986 #ifdef DEBUG_TIME
987         TIMEIT_END(loop_split_worker);
988 #endif
989 }
990
991 /* Note we use data_buff to detect whether we are in threaded context or not, in later case it is NULL. */
992 static void loop_split_generator_do(LoopSplitTaskDataCommon *common_data, const bool threaded)
993 {
994         MLoopNorSpaceArray *lnors_spacearr = common_data->lnors_spacearr;
995         BLI_bitmap *sharp_verts = common_data->sharp_verts;
996         float (*loopnors)[3] = common_data->loopnors;
997
998         const MLoop *mloops = common_data->mloops;
999         const MPoly *mpolys = common_data->mpolys;
1000         const int (*edge_to_loops)[2] = common_data->edge_to_loops;
1001         const int numPolys = common_data->numPolys;
1002
1003         const MPoly *mp;
1004         int mp_index;
1005
1006         LoopSplitTaskData *data, *data_buff = NULL, data_mem;
1007         int data_idx = 0;
1008
1009         /* Temp edge vectors stack, only used when computing lnor spacearr (and we are not multi-threading). */
1010         BLI_Stack *edge_vectors = (lnors_spacearr && !data_buff) ? BLI_stack_new(sizeof(float[3]), __func__) : NULL;
1011
1012 #ifdef DEBUG_TIME
1013         TIMEIT_START(loop_split_generator);
1014 #endif
1015
1016         if (!threaded) {
1017                 memset(&data_mem, 0, sizeof(data_mem));
1018                 data = &data_mem;
1019         }
1020
1021         /* We now know edges that can be smoothed (with their vector, and their two loops), and edges that will be hard!
1022          * Now, time to generate the normals.
1023          */
1024         for (mp = mpolys, mp_index = 0; mp_index < numPolys; mp++, mp_index++) {
1025                 const MLoop *ml_curr, *ml_prev;
1026                 float (*lnors)[3];
1027                 const int ml_last_index = (mp->loopstart + mp->totloop) - 1;
1028                 int ml_curr_index = mp->loopstart;
1029                 int ml_prev_index = ml_last_index;
1030
1031                 ml_curr = &mloops[ml_curr_index];
1032                 ml_prev = &mloops[ml_prev_index];
1033                 lnors = &loopnors[ml_curr_index];
1034
1035                 for (; ml_curr_index <= ml_last_index; ml_curr++, ml_curr_index++, lnors++) {
1036                         const int *e2l_curr = edge_to_loops[ml_curr->e];
1037                         const int *e2l_prev = edge_to_loops[ml_prev->e];
1038
1039                         if (!IS_EDGE_SHARP(e2l_curr) && (!lnors_spacearr || BLI_BITMAP_TEST_BOOL(sharp_verts, ml_curr->v))) {
1040                                 /* A smooth edge, and we are not generating lnor_spacearr, or the related vertex is sharp.
1041                                  * We skip it because it is either:
1042                                  * - in the middle of a 'smooth fan' already computed (or that will be as soon as we hit
1043                                  *   one of its ends, i.e. one of its two sharp edges), or...
1044                                  * - the related vertex is a "full smooth" one, in which case pre-populated normals from vertex
1045                                  *   are just fine (or it has already be handled in a previous loop in case of needed lnors spacearr)!
1046                                  */
1047                                 /* printf("Skipping loop %d / edge %d / vert %d(%d)\n", ml_curr_index, ml_curr->e, ml_curr->v, sharp_verts[ml_curr->v]); */
1048                         }
1049                         else {
1050                                 if (threaded) {
1051                                         if (data_idx == 0) {
1052                                                 data_buff = MEM_callocN(sizeof(*data_buff) * LOOP_SPLIT_TASK_BLOCK_SIZE, __func__);
1053                                         }
1054                                         data = &data_buff[data_idx];
1055                                 }
1056
1057                                 if (IS_EDGE_SHARP(e2l_curr) && IS_EDGE_SHARP(e2l_prev)) {
1058                                         data->lnor = lnors;
1059                                         data->ml_curr = ml_curr;
1060                                         data->ml_prev = ml_prev;
1061                                         data->ml_curr_index = ml_curr_index;
1062 #if 0  /* Not needed for 'single' loop. */
1063                                         data->ml_prev_index = ml_prev_index;
1064                                         data->e2l_prev = NULL;  /* Tag as 'single' task. */
1065 #endif
1066                                         data->mp_index = mp_index;
1067                                         if (lnors_spacearr) {
1068                                                 data->lnor_space = BKE_lnor_space_create(lnors_spacearr);
1069                                         }
1070                                 }
1071                                 /* We *do not need* to check/tag loops as already computed!
1072                                  * Due to the fact a loop only links to one of its two edges, a same fan *will never be walked
1073                                  * more than once!*
1074                                  * Since we consider edges having neighbor polys with inverted (flipped) normals as sharp, we are sure
1075                                  * that no fan will be skipped, even only considering the case (sharp curr_edge, smooth prev_edge),
1076                                  * and not the alternative (smooth curr_edge, sharp prev_edge).
1077                                  * All this due/thanks to link between normals and loop ordering (i.e. winding).
1078                                  */
1079                                 else {
1080 #if 0  /* Not needed for 'fan' loops. */
1081                                         data->lnor = lnors;
1082 #endif
1083                                         data->ml_curr = ml_curr;
1084                                         data->ml_prev = ml_prev;
1085                                         data->ml_curr_index = ml_curr_index;
1086                                         data->ml_prev_index = ml_prev_index;
1087                                         data->e2l_prev = e2l_prev;  /* Also tag as 'fan' task. */
1088                                         data->mp_index = mp_index;
1089                                         if (lnors_spacearr) {
1090                                                 data->lnor_space = BKE_lnor_space_create(lnors_spacearr);
1091                                                 /* Tag related vertex as sharp, to avoid fanning around it again (in case it was a smooth one).
1092                                                  * This *has* to be done outside of workers tasks! */
1093                                                 BLI_BITMAP_ENABLE(sharp_verts, ml_curr->v);
1094                                         }
1095                                 }
1096
1097                                 if (threaded) {
1098                                         data_idx++;
1099                                         if (data_idx == LOOP_SPLIT_TASK_BLOCK_SIZE) {
1100                                                 BLI_thread_queue_push(common_data->task_queue, data_buff);
1101                                                 data_idx = 0;
1102                                         }
1103                                 }
1104                                 else {
1105                                         loop_split_worker_do(common_data, data, edge_vectors);
1106                                         memset(data, 0, sizeof(data_mem));
1107                                 }
1108                         }
1109
1110                         ml_prev = ml_curr;
1111                         ml_prev_index = ml_curr_index;
1112                 }
1113         }
1114
1115         if (threaded) {
1116                 /* Last block of data... Since it is calloc'ed and we use first NULL item as stopper, everything is fine. */
1117                 if (LIKELY(data_idx)) {
1118                         BLI_thread_queue_push(common_data->task_queue, data_buff);
1119                 }
1120
1121                 /* This will signal all other worker threads to wake up and finish! */
1122                 BLI_thread_queue_nowait(common_data->task_queue);
1123         }
1124
1125         if (edge_vectors) {
1126                 BLI_stack_free(edge_vectors);
1127         }
1128
1129 #ifdef DEBUG_TIME
1130         TIMEIT_END(loop_split_generator);
1131 #endif
1132 }
1133
1134 static void loop_split_generator(TaskPool * __restrict UNUSED(pool), void *taskdata, int UNUSED(threadid))
1135 {
1136         LoopSplitTaskDataCommon *common_data = taskdata;
1137
1138         loop_split_generator_do(common_data, true);
1139 }
1140
1141 /**
1142  * Compute split normals, i.e. vertex normals associated with each poly (hence 'loop normals').
1143  * Useful to materialize sharp edges (or non-smooth faces) without actually modifying the geometry (splitting edges).
1144  */
1145 void BKE_mesh_normals_loop_split(
1146         const MVert *mverts, const int numVerts, MEdge *medges, const int numEdges,
1147         MLoop *mloops, float (*r_loopnors)[3], const int numLoops,
1148         MPoly *mpolys, const float (*polynors)[3], const int numPolys,
1149         const bool use_split_normals, float split_angle,
1150         MLoopNorSpaceArray *r_lnors_spacearr, short (*clnors_data)[2], int *r_loop_to_poly)
1151 {
1152
1153         /* For now this is not supported. If we do not use split normals, we do not generate anything fancy! */
1154         BLI_assert(use_split_normals || !(r_lnors_spacearr));
1155
1156         if (!use_split_normals) {
1157                 /* In this case, we simply fill lnors with vnors (or fnors for flat faces), quite simple!
1158                  * Note this is done here to keep some logic and consistency in this quite complex code,
1159                  * since we may want to use lnors even when mesh's 'autosmooth' is disabled (see e.g. mesh mapping code).
1160                  * As usual, we could handle that on case-by-case basis, but simpler to keep it well confined here.
1161                  */
1162                 int mp_index;
1163
1164                 for (mp_index = 0; mp_index < numPolys; mp_index++) {
1165                         MPoly *mp = &mpolys[mp_index];
1166                         int ml_index = mp->loopstart;
1167                         const int ml_index_end = ml_index + mp->totloop;
1168                         const bool is_poly_flat = ((mp->flag & ME_SMOOTH) == 0);
1169
1170                         for (; ml_index < ml_index_end; ml_index++) {
1171                                 if (r_loop_to_poly) {
1172                                         r_loop_to_poly[ml_index] = mp_index;
1173                                 }
1174                                 if (is_poly_flat) {
1175                                         copy_v3_v3(r_loopnors[ml_index], polynors[mp_index]);
1176                                 }
1177                                 else {
1178                                         normal_short_to_float_v3(r_loopnors[ml_index], mverts[mloops[ml_index].v].no);
1179                                 }
1180                         }
1181                 }
1182                 return;
1183         }
1184
1185         {
1186
1187         /* Mapping edge -> loops.
1188          * If that edge is used by more than two loops (polys), it is always sharp (and tagged as such, see below).
1189          * We also use the second loop index as a kind of flag: smooth edge: > 0,
1190          *                                                      sharp edge: < 0 (INDEX_INVALID || INDEX_UNSET),
1191          *                                                      unset: INDEX_UNSET
1192          * Note that currently we only have two values for second loop of sharp edges. However, if needed, we can
1193          * store the negated value of loop index instead of INDEX_INVALID to retrieve the real value later in code).
1194          * Note also that lose edges always have both values set to 0!
1195          */
1196         int (*edge_to_loops)[2] = MEM_callocN(sizeof(int[2]) * (size_t)numEdges, __func__);
1197
1198         /* Simple mapping from a loop to its polygon index. */
1199         int *loop_to_poly = r_loop_to_poly ? r_loop_to_poly : MEM_mallocN(sizeof(int) * (size_t)numLoops, __func__);
1200
1201         MPoly *mp;
1202         int mp_index, me_index;
1203         bool check_angle = (split_angle < (float)M_PI);
1204         int i;
1205
1206         BLI_bitmap *sharp_verts = NULL;
1207         MLoopNorSpaceArray _lnors_spacearr = {NULL};
1208
1209         LoopSplitTaskDataCommon common_data = {NULL};
1210
1211 #ifdef DEBUG_TIME
1212         TIMEIT_START(BKE_mesh_normals_loop_split);
1213 #endif
1214
1215         if (check_angle) {
1216                 /* When using custom loop normals, disable the angle feature! */
1217                 if (clnors_data) {
1218                         check_angle = false;
1219                 }
1220                 else {
1221                         split_angle = cosf(split_angle);
1222                 }
1223         }
1224
1225         if (!r_lnors_spacearr && clnors_data) {
1226                 /* We need to compute lnor spacearr if some custom lnor data are given to us! */
1227                 r_lnors_spacearr = &_lnors_spacearr;
1228         }
1229         if (r_lnors_spacearr) {
1230                 BKE_lnor_spacearr_init(r_lnors_spacearr, numLoops);
1231                 sharp_verts = BLI_BITMAP_NEW((size_t)numVerts, __func__);
1232         }
1233
1234         /* This first loop check which edges are actually smooth, and compute edge vectors. */
1235         for (mp = mpolys, mp_index = 0; mp_index < numPolys; mp++, mp_index++) {
1236                 MLoop *ml_curr;
1237                 int *e2l;
1238                 int ml_curr_index = mp->loopstart;
1239                 const int ml_last_index = (ml_curr_index + mp->totloop) - 1;
1240
1241                 ml_curr = &mloops[ml_curr_index];
1242
1243                 for (; ml_curr_index <= ml_last_index; ml_curr++, ml_curr_index++) {
1244                         e2l = edge_to_loops[ml_curr->e];
1245
1246                         loop_to_poly[ml_curr_index] = mp_index;
1247
1248                         /* Pre-populate all loop normals as if their verts were all-smooth, this way we don't have to compute
1249                          * those later!
1250                          */
1251                         normal_short_to_float_v3(r_loopnors[ml_curr_index], mverts[ml_curr->v].no);
1252
1253                         /* Check whether current edge might be smooth or sharp */
1254                         if ((e2l[0] | e2l[1]) == 0) {
1255                                 /* 'Empty' edge until now, set e2l[0] (and e2l[1] to INDEX_UNSET to tag it as unset). */
1256                                 e2l[0] = ml_curr_index;
1257                                 /* We have to check this here too, else we might miss some flat faces!!! */
1258                                 e2l[1] = (mp->flag & ME_SMOOTH) ? INDEX_UNSET : INDEX_INVALID;
1259                         }
1260                         else if (e2l[1] == INDEX_UNSET) {
1261                                 /* Second loop using this edge, time to test its sharpness.
1262                                  * An edge is sharp if it is tagged as such, or its face is not smooth,
1263                                  * or both poly have opposed (flipped) normals, i.e. both loops on the same edge share the same vertex,
1264                                  * or angle between both its polys' normals is above split_angle value.
1265                                  */
1266                                 if (!(mp->flag & ME_SMOOTH) || (medges[ml_curr->e].flag & ME_SHARP) ||
1267                                     ml_curr->v == mloops[e2l[0]].v ||
1268                                     (check_angle && dot_v3v3(polynors[loop_to_poly[e2l[0]]], polynors[mp_index]) < split_angle))
1269                                 {
1270                                         /* Note: we are sure that loop != 0 here ;) */
1271                                         e2l[1] = INDEX_INVALID;
1272                                 }
1273                                 else {
1274                                         e2l[1] = ml_curr_index;
1275                                 }
1276                         }
1277                         else if (!IS_EDGE_SHARP(e2l)) {
1278                                 /* More than two loops using this edge, tag as sharp if not yet done. */
1279                                 e2l[1] = INDEX_INVALID;
1280                         }
1281                         /* Else, edge is already 'disqualified' (i.e. sharp)! */
1282                 }
1283         }
1284
1285         if (r_lnors_spacearr) {
1286                 /* Tag vertices that have at least one sharp edge as 'sharp' (used for the lnor spacearr computation).
1287                  * XXX This third loop over edges is a bit disappointing, could not find any other way yet.
1288                  *     Not really performance-critical anyway.
1289                  */
1290                 for (me_index = 0; me_index < numEdges; me_index++) {
1291                         const int *e2l = edge_to_loops[me_index];
1292                         const MEdge *me = &medges[me_index];
1293                         if (IS_EDGE_SHARP(e2l)) {
1294                                 BLI_BITMAP_ENABLE(sharp_verts, me->v1);
1295                                 BLI_BITMAP_ENABLE(sharp_verts, me->v2);
1296                         }
1297                 }
1298         }
1299
1300         /* Init data common to all tasks. */
1301         common_data.lnors_spacearr = r_lnors_spacearr;
1302         common_data.loopnors = r_loopnors;
1303         common_data.clnors_data = clnors_data;
1304
1305         common_data.mverts = mverts;
1306         common_data.medges = medges;
1307         common_data.mloops = mloops;
1308         common_data.mpolys = mpolys;
1309         common_data.sharp_verts = sharp_verts;
1310         common_data.edge_to_loops = (const int(*)[2])edge_to_loops;
1311         common_data.loop_to_poly = loop_to_poly;
1312         common_data.polynors = polynors;
1313         common_data.numPolys = numPolys;
1314
1315         if (numLoops < LOOP_SPLIT_TASK_BLOCK_SIZE * 8) {
1316                 /* Not enough loops to be worth the whole threading overhead... */
1317                 loop_split_generator_do(&common_data, false);
1318         }
1319         else {
1320                 TaskScheduler *task_scheduler;
1321                 TaskPool *task_pool;
1322                 int nbr_workers;
1323
1324                 common_data.task_queue = BLI_thread_queue_init();
1325
1326                 task_scheduler = BLI_task_scheduler_get();
1327                 task_pool = BLI_task_pool_create(task_scheduler, NULL);
1328
1329                 nbr_workers = max_ii(2, BLI_task_scheduler_num_threads(task_scheduler));
1330                 for (i = 1; i < nbr_workers; i++) {
1331                         BLI_task_pool_push(task_pool, loop_split_worker, &common_data, false, TASK_PRIORITY_HIGH);
1332                 }
1333                 BLI_task_pool_push(task_pool, loop_split_generator, &common_data, false, TASK_PRIORITY_HIGH);
1334                 BLI_task_pool_work_and_wait(task_pool);
1335
1336                 BLI_task_pool_free(task_pool);
1337
1338                 BLI_thread_queue_free(common_data.task_queue);
1339         }
1340
1341         MEM_freeN(edge_to_loops);
1342         if (!r_loop_to_poly) {
1343                 MEM_freeN(loop_to_poly);
1344         }
1345
1346         if (r_lnors_spacearr) {
1347                 MEM_freeN(sharp_verts);
1348                 if (r_lnors_spacearr == &_lnors_spacearr) {
1349                         BKE_lnor_spacearr_free(r_lnors_spacearr);
1350                 }
1351         }
1352
1353 #ifdef DEBUG_TIME
1354         TIMEIT_END(BKE_mesh_normals_loop_split);
1355 #endif
1356
1357         }
1358 }
1359
1360 #undef INDEX_UNSET
1361 #undef INDEX_INVALID
1362 #undef IS_EDGE_SHARP
1363
1364 /**
1365  * Compute internal representation of given custom normals (as an array of float[2]).
1366  * It also makes sure the mesh matches those custom normals, by setting sharp edges flag as needed to get a
1367  * same custom lnor for all loops sharing a same smooth fan.
1368  * If use_vertices if true, r_custom_loopnors is assumed to be per-vertex, not per-loop
1369  * (this allows to set whole vert's normals at once, useful in some cases).
1370  * r_custom_loopnors is expected to have normalized normals, or zero ones, in which case they will be replaced
1371  * by default loop/vertex normal.
1372  */
1373 static void mesh_normals_loop_custom_set(
1374         const MVert *mverts, const int numVerts, MEdge *medges, const int numEdges,
1375         MLoop *mloops, float (*r_custom_loopnors)[3], const int numLoops,
1376         MPoly *mpolys, const float (*polynors)[3], const int numPolys,
1377         short (*r_clnors_data)[2], const bool use_vertices)
1378 {
1379         /* We *may* make that poor BKE_mesh_normals_loop_split() even more complex by making it handling that
1380          * feature too, would probably be more efficient in absolute.
1381          * However, this function *is not* performance-critical, since it is mostly expected to be called
1382          * by io addons when importing custom normals, and modifier (and perhaps from some editing tools later?).
1383          * So better to keep some simplicity here, and just call BKE_mesh_normals_loop_split() twice!
1384          */
1385         MLoopNorSpaceArray lnors_spacearr = {NULL};
1386         BLI_bitmap *done_loops = BLI_BITMAP_NEW((size_t)numLoops, __func__);
1387         float (*lnors)[3] = MEM_callocN(sizeof(*lnors) * (size_t)numLoops, __func__);
1388         int *loop_to_poly = MEM_mallocN(sizeof(int) * (size_t)numLoops, __func__);
1389         /* In this case we always consider split nors as ON, and do not want to use angle to define smooth fans! */
1390         const bool use_split_normals = true;
1391         const float split_angle = (float)M_PI;
1392         int i;
1393
1394         BLI_SMALLSTACK_DECLARE(clnors_data, short *);
1395
1396         /* Compute current lnor spacearr. */
1397         BKE_mesh_normals_loop_split(mverts, numVerts, medges, numEdges, mloops, lnors, numLoops,
1398                                     mpolys, polynors, numPolys, use_split_normals, split_angle,
1399                                     &lnors_spacearr, NULL, loop_to_poly);
1400
1401         /* Set all given zero vectors to their default value. */
1402         if (use_vertices) {
1403                 for (i = 0; i < numVerts; i++) {
1404                         if (is_zero_v3(r_custom_loopnors[i])) {
1405                                 normal_short_to_float_v3(r_custom_loopnors[i], mverts[i].no);
1406                         }
1407                 }
1408         }
1409         else {
1410                 for (i = 0; i < numLoops; i++) {
1411                         if (is_zero_v3(r_custom_loopnors[i])) {
1412                                 copy_v3_v3(r_custom_loopnors[i], lnors[i]);
1413                         }
1414                 }
1415         }
1416
1417         /* Now, check each current smooth fan (one lnor space per smooth fan!), and if all its matching custom lnors
1418          * are not (enough) equal, add sharp edges as needed.
1419          * This way, next time we run BKE_mesh_normals_loop_split(), we'll get lnor spacearr/smooth fans matching
1420          * given custom lnors.
1421          * Note this code *will never* unsharp edges!
1422          * And quite obviously, when we set custom normals per vertices, running this is absolutely useless.
1423          */
1424         if (!use_vertices) {
1425                 for (i = 0; i < numLoops; i++) {
1426                         if (!lnors_spacearr.lspacearr[i]) {
1427                                 /* This should not happen in theory, but in some rare case (probably ugly geometry)
1428                                  * we can get some NULL loopspacearr at this point. :/
1429                                  * Maybe we should set those loops' edges as sharp?
1430                                  */
1431                                 BLI_BITMAP_ENABLE(done_loops, i);
1432                                 if (G.debug & G_DEBUG) {
1433                                         printf("WARNING! Getting invalid NULL loop space for loop %d!\n", i);
1434                                 }
1435                                 continue;
1436                         }
1437
1438                         if (!BLI_BITMAP_TEST(done_loops, i)) {
1439                                 /* Notes:
1440                                  *     * In case of mono-loop smooth fan, loops is NULL, so everything is fine (we have nothing to do).
1441                                  *     * Loops in this linklist are ordered (in reversed order compared to how they were discovered by
1442                                  *       BKE_mesh_normals_loop_split(), but this is not a problem). Which means if we find a
1443                                  *       mismatching clnor, we know all remaining loops will have to be in a new, different smooth fan/
1444                                  *       lnor space.
1445                                  *     * In smooth fan case, we compare each clnor against a ref one, to avoid small differences adding
1446                                  *       up into a real big one in the end!
1447                                  */
1448                                 LinkNode *loops = lnors_spacearr.lspacearr[i]->loops;
1449                                 MLoop *prev_ml = NULL;
1450                                 const float *org_nor = NULL;
1451
1452                                 while (loops) {
1453                                         const int lidx = GET_INT_FROM_POINTER(loops->link);
1454                                         MLoop *ml = &mloops[lidx];
1455                                         const int nidx = lidx;
1456                                         float *nor = r_custom_loopnors[nidx];
1457
1458                                         if (!org_nor) {
1459                                                 org_nor = nor;
1460                                         }
1461                                         else if (dot_v3v3(org_nor, nor) < LNOR_SPACE_TRIGO_THRESHOLD) {
1462                                                 /* Current normal differs too much from org one, we have to tag the edge between
1463                                                  * previous loop's face and current's one as sharp.
1464                                                  * We know those two loops do not point to the same edge, since we do not allow reversed winding
1465                                                  * in a same smooth fan.
1466                                                  */
1467                                                 const MPoly *mp = &mpolys[loop_to_poly[lidx]];
1468                                                 const MLoop *mlp = &mloops[(lidx == mp->loopstart) ? mp->loopstart + mp->totloop - 1 : lidx - 1];
1469                                                 medges[(prev_ml->e == mlp->e) ? prev_ml->e : ml->e].flag |= ME_SHARP;
1470
1471                                                 org_nor = nor;
1472                                         }
1473
1474                                         prev_ml = ml;
1475                                         loops = loops->next;
1476                                         BLI_BITMAP_ENABLE(done_loops, lidx);
1477                                 }
1478
1479                                 /* We also have to check between last and first loops, otherwise we may miss some sharp edges here!
1480                                  * This is just a simplified version of above while loop.
1481                                  * See T45984. */
1482                                 loops = lnors_spacearr.lspacearr[i]->loops;
1483                                 if (loops && org_nor) {
1484                                         const int lidx = GET_INT_FROM_POINTER(loops->link);
1485                                         MLoop *ml = &mloops[lidx];
1486                                         const int nidx = lidx;
1487                                         float *nor = r_custom_loopnors[nidx];
1488
1489                                         if (dot_v3v3(org_nor, nor) < LNOR_SPACE_TRIGO_THRESHOLD) {
1490                                                 const MPoly *mp = &mpolys[loop_to_poly[lidx]];
1491                                                 const MLoop *mlp = &mloops[(lidx == mp->loopstart) ? mp->loopstart + mp->totloop - 1 : lidx - 1];
1492                                                 medges[(prev_ml->e == mlp->e) ? prev_ml->e : ml->e].flag |= ME_SHARP;
1493                                         }
1494                                 }
1495
1496                                 /* For single loops, where lnors_spacearr.lspacearr[i]->loops is NULL. */
1497                                 BLI_BITMAP_ENABLE(done_loops, i);
1498                         }
1499                 }
1500
1501                 /* And now, recompute our new auto lnors and lnor spacearr! */
1502                 BKE_lnor_spacearr_clear(&lnors_spacearr);
1503                 BKE_mesh_normals_loop_split(mverts, numVerts, medges, numEdges, mloops, lnors, numLoops,
1504                                             mpolys, polynors, numPolys, use_split_normals, split_angle,
1505                                             &lnors_spacearr, NULL, loop_to_poly);
1506         }
1507         else {
1508                 BLI_BITMAP_SET_ALL(done_loops, true, (size_t)numLoops);
1509         }
1510
1511         /* And we just have to convert plain object-space custom normals to our lnor space-encoded ones. */
1512         for (i = 0; i < numLoops; i++) {
1513                 if (!lnors_spacearr.lspacearr[i]) {
1514                         BLI_BITMAP_DISABLE(done_loops, i);
1515                         if (G.debug & G_DEBUG) {
1516                                 printf("WARNING! Still getting invalid NULL loop space in second loop for loop %d!\n", i);
1517                         }
1518                         continue;
1519                 }
1520
1521                 if (BLI_BITMAP_TEST_BOOL(done_loops, i)) {
1522                         /* Note we accumulate and average all custom normals in current smooth fan, to avoid getting different
1523                          * clnors data (tiny differences in plain custom normals can give rather huge differences in
1524                          * computed 2D factors).
1525                          */
1526                         LinkNode *loops = lnors_spacearr.lspacearr[i]->loops;
1527                         if (loops) {
1528                                 int nbr_nors = 0;
1529                                 float avg_nor[3];
1530                                 short clnor_data_tmp[2], *clnor_data;
1531
1532                                 zero_v3(avg_nor);
1533                                 while (loops) {
1534                                         const int lidx = GET_INT_FROM_POINTER(loops->link);
1535                                         const int nidx = use_vertices ? (int)mloops[lidx].v : lidx;
1536                                         float *nor = r_custom_loopnors[nidx];
1537
1538                                         nbr_nors++;
1539                                         add_v3_v3(avg_nor, nor);
1540                                         BLI_SMALLSTACK_PUSH(clnors_data, (short *)r_clnors_data[lidx]);
1541
1542                                         loops = loops->next;
1543                                         BLI_BITMAP_DISABLE(done_loops, lidx);
1544                                 }
1545
1546                                 mul_v3_fl(avg_nor, 1.0f / (float)nbr_nors);
1547                                 BKE_lnor_space_custom_normal_to_data(lnors_spacearr.lspacearr[i], avg_nor, clnor_data_tmp);
1548
1549                                 while ((clnor_data = BLI_SMALLSTACK_POP(clnors_data))) {
1550                                         clnor_data[0] = clnor_data_tmp[0];
1551                                         clnor_data[1] = clnor_data_tmp[1];
1552                                 }
1553                         }
1554                         else {
1555                                 const int nidx = use_vertices ? (int)mloops[i].v : i;
1556                                 float *nor = r_custom_loopnors[nidx];
1557
1558                                 BKE_lnor_space_custom_normal_to_data(lnors_spacearr.lspacearr[i], nor, r_clnors_data[i]);
1559                                 BLI_BITMAP_DISABLE(done_loops, i);
1560                         }
1561                 }
1562         }
1563
1564         MEM_freeN(lnors);
1565         MEM_freeN(loop_to_poly);
1566         MEM_freeN(done_loops);
1567         BKE_lnor_spacearr_free(&lnors_spacearr);
1568 }
1569
1570 void BKE_mesh_normals_loop_custom_set(
1571         const MVert *mverts, const int numVerts, MEdge *medges, const int numEdges,
1572         MLoop *mloops, float (*r_custom_loopnors)[3], const int numLoops,
1573         MPoly *mpolys, const float (*polynors)[3], const int numPolys,
1574         short (*r_clnors_data)[2])
1575 {
1576         mesh_normals_loop_custom_set(mverts, numVerts, medges, numEdges, mloops, r_custom_loopnors, numLoops,
1577                                      mpolys, polynors, numPolys, r_clnors_data, false);
1578 }
1579
1580 void BKE_mesh_normals_loop_custom_from_vertices_set(
1581         const MVert *mverts, float (*r_custom_vertnors)[3], const int numVerts,
1582         MEdge *medges, const int numEdges, MLoop *mloops, const int numLoops,
1583         MPoly *mpolys, const float (*polynors)[3], const int numPolys,
1584         short (*r_clnors_data)[2])
1585 {
1586         mesh_normals_loop_custom_set(mverts, numVerts, medges, numEdges, mloops, r_custom_vertnors, numLoops,
1587                                      mpolys, polynors, numPolys, r_clnors_data, true);
1588 }
1589
1590 /**
1591  * Computes average per-vertex normals from given custom loop normals.
1592  *
1593  * @param clnors The computed custom loop normals.
1594  * @param r_vert_clnors The (already allocated) array where to store averaged per-vertex normals.
1595  */
1596 void BKE_mesh_normals_loop_to_vertex(
1597         const int numVerts, const MLoop *mloops, const int numLoops,
1598         const float (*clnors)[3], float (*r_vert_clnors)[3])
1599 {
1600         const MLoop *ml;
1601         int i;
1602
1603         int *vert_loops_nbr = MEM_callocN(sizeof(*vert_loops_nbr) * (size_t)numVerts, __func__);
1604
1605         copy_vn_fl((float *)r_vert_clnors, 3 * numVerts, 0.0f);
1606
1607         for (i = 0, ml = mloops; i < numLoops; i++, ml++) {
1608                 const unsigned int v = ml->v;
1609
1610                 add_v3_v3(r_vert_clnors[v], clnors[i]);
1611                 vert_loops_nbr[v]++;
1612         }
1613
1614         for (i = 0; i < numVerts; i++) {
1615                 mul_v3_fl(r_vert_clnors[i], 1.0f / (float)vert_loops_nbr[i]);
1616         }
1617
1618         MEM_freeN(vert_loops_nbr);
1619 }
1620
1621
1622 #undef LNOR_SPACE_TRIGO_THRESHOLD
1623
1624 /** \} */
1625
1626
1627 /* -------------------------------------------------------------------- */
1628
1629 /** \name Mesh Tangent Calculations
1630  * \{ */
1631
1632 /* Tangent space utils. */
1633
1634 /* User data. */
1635 typedef struct {
1636         const MPoly *mpolys;   /* faces */
1637         const MLoop *mloops;   /* faces's vertices */
1638         const MVert *mverts;   /* vertices */
1639         const MLoopUV *luvs;   /* texture coordinates */
1640         float (*lnors)[3];     /* loops' normals */
1641         float (*tangents)[4];  /* output tangents */
1642         int num_polys;         /* number of polygons */
1643 } BKEMeshToTangent;
1644
1645 /* Mikktspace's API */
1646 static int get_num_faces(const SMikkTSpaceContext *pContext)
1647 {
1648         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1649         return p_mesh->num_polys;
1650 }
1651
1652 static int get_num_verts_of_face(const SMikkTSpaceContext *pContext, const int face_idx)
1653 {
1654         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1655         return p_mesh->mpolys[face_idx].totloop;
1656 }
1657
1658 static void get_position(const SMikkTSpaceContext *pContext, float r_co[3], const int face_idx, const int vert_idx)
1659 {
1660         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1661         const int loop_idx = p_mesh->mpolys[face_idx].loopstart + vert_idx;
1662         copy_v3_v3(r_co, p_mesh->mverts[p_mesh->mloops[loop_idx].v].co);
1663 }
1664
1665 static void get_texture_coordinate(const SMikkTSpaceContext *pContext, float r_uv[2], const int face_idx,
1666                                    const int vert_idx)
1667 {
1668         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1669         copy_v2_v2(r_uv, p_mesh->luvs[p_mesh->mpolys[face_idx].loopstart + vert_idx].uv);
1670 }
1671
1672 static void get_normal(const SMikkTSpaceContext *pContext, float r_no[3], const int face_idx, const int vert_idx)
1673 {
1674         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1675         copy_v3_v3(r_no, p_mesh->lnors[p_mesh->mpolys[face_idx].loopstart + vert_idx]);
1676 }
1677
1678 static void set_tspace(const SMikkTSpaceContext *pContext, const float fv_tangent[3], const float face_sign,
1679                        const int face_idx, const int vert_idx)
1680 {
1681         BKEMeshToTangent *p_mesh = (BKEMeshToTangent *)pContext->m_pUserData;
1682         float *p_res = p_mesh->tangents[p_mesh->mpolys[face_idx].loopstart + vert_idx];
1683         copy_v3_v3(p_res, fv_tangent);
1684         p_res[3] = face_sign;
1685 }
1686
1687 /**
1688  * Compute simplified tangent space normals, i.e. tangent vector + sign of bi-tangent one, which combined with
1689  * split normals can be used to recreate the full tangent space.
1690  * Note: * The mesh should be made of only tris and quads!
1691  */
1692 void BKE_mesh_loop_tangents_ex(
1693         const MVert *mverts, const int UNUSED(numVerts), const MLoop *mloops,
1694         float (*r_looptangent)[4], float (*loopnors)[3], const MLoopUV *loopuvs,
1695         const int UNUSED(numLoops), const MPoly *mpolys, const int numPolys,
1696         ReportList *reports)
1697 {
1698         BKEMeshToTangent mesh_to_tangent = {NULL};
1699         SMikkTSpaceContext s_context = {NULL};
1700         SMikkTSpaceInterface s_interface = {NULL};
1701
1702         const MPoly *mp;
1703         int mp_index;
1704
1705         /* First check we do have a tris/quads only mesh. */
1706         for (mp = mpolys, mp_index = 0; mp_index < numPolys; mp++, mp_index++) {
1707                 if (mp->totloop > 4) {
1708                         BKE_report(reports, RPT_ERROR, "Tangent space can only be computed for tris/quads, aborting");
1709                         return;
1710                 }
1711         }
1712
1713         /* Compute Mikktspace's tangent normals. */
1714         mesh_to_tangent.mpolys = mpolys;
1715         mesh_to_tangent.mloops = mloops;
1716         mesh_to_tangent.mverts = mverts;
1717         mesh_to_tangent.luvs = loopuvs;
1718         mesh_to_tangent.lnors = loopnors;
1719         mesh_to_tangent.tangents = r_looptangent;
1720         mesh_to_tangent.num_polys = numPolys;
1721
1722         s_context.m_pUserData = &mesh_to_tangent;
1723         s_context.m_pInterface = &s_interface;
1724         s_interface.m_getNumFaces = get_num_faces;
1725         s_interface.m_getNumVerticesOfFace = get_num_verts_of_face;
1726         s_interface.m_getPosition = get_position;
1727         s_interface.m_getTexCoord = get_texture_coordinate;
1728         s_interface.m_getNormal = get_normal;
1729         s_interface.m_setTSpaceBasic = set_tspace;
1730
1731         /* 0 if failed */
1732         if (genTangSpaceDefault(&s_context) == false) {
1733                 BKE_report(reports, RPT_ERROR, "Mikktspace failed to generate tangents for this mesh!");
1734         }
1735 }
1736
1737 /**
1738  * Wrapper around BKE_mesh_loop_tangents_ex, which takes care of most boiling code.
1739  * \note
1740  * - There must be a valid loop's CD_NORMALS available.
1741  * - The mesh should be made of only tris and quads!
1742  */
1743 void BKE_mesh_loop_tangents(Mesh *mesh, const char *uvmap, float (*r_looptangents)[4], ReportList *reports)
1744 {
1745         MLoopUV *loopuvs;
1746         float (*loopnors)[3];
1747
1748         /* Check we have valid texture coordinates first! */
1749         if (uvmap) {
1750                 loopuvs = CustomData_get_layer_named(&mesh->ldata, CD_MLOOPUV, uvmap);
1751         }
1752         else {
1753                 loopuvs = CustomData_get_layer(&mesh->ldata, CD_MLOOPUV);
1754         }
1755         if (!loopuvs) {
1756                 BKE_reportf(reports, RPT_ERROR, "Tangent space computation needs an UVMap, \"%s\" not found, aborting", uvmap);
1757                 return;
1758         }
1759
1760         loopnors = CustomData_get_layer(&mesh->ldata, CD_NORMAL);
1761         if (!loopnors) {
1762                 BKE_report(reports, RPT_ERROR, "Tangent space computation needs loop normals, none found, aborting");
1763                 return;
1764         }
1765
1766         BKE_mesh_loop_tangents_ex(mesh->mvert, mesh->totvert, mesh->mloop, r_looptangents,
1767                                   loopnors, loopuvs, mesh->totloop, mesh->mpoly, mesh->totpoly, reports);
1768 }
1769
1770 /** \} */
1771
1772
1773 /* -------------------------------------------------------------------- */
1774
1775 /** \name Polygon Calculations
1776  * \{ */
1777
1778 /*
1779  * COMPUTE POLY NORMAL
1780  *
1781  * Computes the normal of a planar
1782  * polygon See Graphics Gems for
1783  * computing newell normal.
1784  *
1785  */
1786 static void mesh_calc_ngon_normal(
1787         const MPoly *mpoly, const MLoop *loopstart,
1788         const MVert *mvert, float normal[3])
1789 {
1790         const int nverts = mpoly->totloop;
1791         const float *v_prev = mvert[loopstart[nverts - 1].v].co;
1792         const float *v_curr;
1793         int i;
1794
1795         zero_v3(normal);
1796
1797         /* Newell's Method */
1798         for (i = 0; i < nverts; i++) {
1799                 v_curr = mvert[loopstart[i].v].co;
1800                 add_newell_cross_v3_v3v3(normal, v_prev, v_curr);
1801                 v_prev = v_curr;
1802         }
1803
1804         if (UNLIKELY(normalize_v3(normal) == 0.0f)) {
1805                 normal[2] = 1.0f; /* other axis set to 0.0 */
1806         }
1807 }
1808
1809 void BKE_mesh_calc_poly_normal(
1810         const MPoly *mpoly, const MLoop *loopstart,
1811         const MVert *mvarray, float r_no[3])
1812 {
1813         if (mpoly->totloop > 4) {
1814                 mesh_calc_ngon_normal(mpoly, loopstart, mvarray, r_no);
1815         }
1816         else if (mpoly->totloop == 3) {
1817                 normal_tri_v3(r_no,
1818                               mvarray[loopstart[0].v].co,
1819                               mvarray[loopstart[1].v].co,
1820                               mvarray[loopstart[2].v].co
1821                               );
1822         }
1823         else if (mpoly->totloop == 4) {
1824                 normal_quad_v3(r_no,
1825                                mvarray[loopstart[0].v].co,
1826                                mvarray[loopstart[1].v].co,
1827                                mvarray[loopstart[2].v].co,
1828                                mvarray[loopstart[3].v].co
1829                                );
1830         }
1831         else { /* horrible, two sided face! */
1832                 r_no[0] = 0.0;
1833                 r_no[1] = 0.0;
1834                 r_no[2] = 1.0;
1835         }
1836 }
1837 /* duplicate of function above _but_ takes coords rather then mverts */
1838 static void mesh_calc_ngon_normal_coords(
1839         const MPoly *mpoly, const MLoop *loopstart,
1840         const float (*vertex_coords)[3], float r_normal[3])
1841 {
1842         const int nverts = mpoly->totloop;
1843         const float *v_prev = vertex_coords[loopstart[nverts - 1].v];
1844         const float *v_curr;
1845         int i;
1846
1847         zero_v3(r_normal);
1848
1849         /* Newell's Method */
1850         for (i = 0; i < nverts; i++) {
1851                 v_curr = vertex_coords[loopstart[i].v];
1852                 add_newell_cross_v3_v3v3(r_normal, v_prev, v_curr);
1853                 v_prev = v_curr;
1854         }
1855
1856         if (UNLIKELY(normalize_v3(r_normal) == 0.0f)) {
1857                 r_normal[2] = 1.0f; /* other axis set to 0.0 */
1858         }
1859 }
1860
1861 void BKE_mesh_calc_poly_normal_coords(
1862         const MPoly *mpoly, const MLoop *loopstart,
1863         const float (*vertex_coords)[3], float r_no[3])
1864 {
1865         if (mpoly->totloop > 4) {
1866                 mesh_calc_ngon_normal_coords(mpoly, loopstart, vertex_coords, r_no);
1867         }
1868         else if (mpoly->totloop == 3) {
1869                 normal_tri_v3(r_no,
1870                               vertex_coords[loopstart[0].v],
1871                               vertex_coords[loopstart[1].v],
1872                               vertex_coords[loopstart[2].v]
1873                               );
1874         }
1875         else if (mpoly->totloop == 4) {
1876                 normal_quad_v3(r_no,
1877                                vertex_coords[loopstart[0].v],
1878                                vertex_coords[loopstart[1].v],
1879                                vertex_coords[loopstart[2].v],
1880                                vertex_coords[loopstart[3].v]
1881                                );
1882         }
1883         else { /* horrible, two sided face! */
1884                 r_no[0] = 0.0;
1885                 r_no[1] = 0.0;
1886                 r_no[2] = 1.0;
1887         }
1888 }
1889
1890 static void mesh_calc_ngon_center(
1891         const MPoly *mpoly, const MLoop *loopstart,
1892         const MVert *mvert, float cent[3])
1893 {
1894         const float w = 1.0f / (float)mpoly->totloop;
1895         int i;
1896
1897         zero_v3(cent);
1898
1899         for (i = 0; i < mpoly->totloop; i++) {
1900                 madd_v3_v3fl(cent, mvert[(loopstart++)->v].co, w);
1901         }
1902 }
1903
1904 void BKE_mesh_calc_poly_center(
1905         const MPoly *mpoly, const MLoop *loopstart,
1906         const MVert *mvarray, float r_cent[3])
1907 {
1908         if (mpoly->totloop == 3) {
1909                 cent_tri_v3(r_cent,
1910                             mvarray[loopstart[0].v].co,
1911                             mvarray[loopstart[1].v].co,
1912                             mvarray[loopstart[2].v].co
1913                             );
1914         }
1915         else if (mpoly->totloop == 4) {
1916                 cent_quad_v3(r_cent,
1917                              mvarray[loopstart[0].v].co,
1918                              mvarray[loopstart[1].v].co,
1919                              mvarray[loopstart[2].v].co,
1920                              mvarray[loopstart[3].v].co
1921                              );
1922         }
1923         else {
1924                 mesh_calc_ngon_center(mpoly, loopstart, mvarray, r_cent);
1925         }
1926 }
1927
1928 /* note, passing polynormal is only a speedup so we can skip calculating it */
1929 float BKE_mesh_calc_poly_area(
1930         const MPoly *mpoly, const MLoop *loopstart,
1931         const MVert *mvarray)
1932 {
1933         if (mpoly->totloop == 3) {
1934                 return area_tri_v3(mvarray[loopstart[0].v].co,
1935                                    mvarray[loopstart[1].v].co,
1936                                    mvarray[loopstart[2].v].co
1937                                    );
1938         }
1939         else {
1940                 int i;
1941                 const MLoop *l_iter = loopstart;
1942                 float area;
1943                 float (*vertexcos)[3] = BLI_array_alloca(vertexcos, (size_t)mpoly->totloop);
1944
1945                 /* pack vertex cos into an array for area_poly_v3 */
1946                 for (i = 0; i < mpoly->totloop; i++, l_iter++) {
1947                         copy_v3_v3(vertexcos[i], mvarray[l_iter->v].co);
1948                 }
1949
1950                 /* finally calculate the area */
1951                 area = area_poly_v3((const float (*)[3])vertexcos, (unsigned int)mpoly->totloop);
1952
1953                 return area;
1954         }
1955 }
1956
1957 /* note, results won't be correct if polygon is non-planar */
1958 static float mesh_calc_poly_planar_area_centroid(
1959         const MPoly *mpoly, const MLoop *loopstart, const MVert *mvarray,
1960         float r_cent[3])
1961 {
1962         int i;
1963         float tri_area;
1964         float total_area = 0.0f;
1965         float v1[3], v2[3], v3[3], normal[3], tri_cent[3];
1966
1967         BKE_mesh_calc_poly_normal(mpoly, loopstart, mvarray, normal);
1968         copy_v3_v3(v1, mvarray[loopstart[0].v].co);
1969         copy_v3_v3(v2, mvarray[loopstart[1].v].co);
1970         zero_v3(r_cent);
1971
1972         for (i = 2; i < mpoly->totloop; i++) {
1973                 copy_v3_v3(v3, mvarray[loopstart[i].v].co);
1974
1975                 tri_area = area_tri_signed_v3(v1, v2, v3, normal);
1976                 total_area += tri_area;
1977
1978                 cent_tri_v3(tri_cent, v1, v2, v3);
1979                 madd_v3_v3fl(r_cent, tri_cent, tri_area);
1980
1981                 copy_v3_v3(v2, v3);
1982         }
1983
1984         mul_v3_fl(r_cent, 1.0f / total_area);
1985
1986         return total_area;
1987 }
1988
1989 #if 0 /* slow version of the function below */
1990 void BKE_mesh_calc_poly_angles(MPoly *mpoly, MLoop *loopstart,
1991                                MVert *mvarray, float angles[])
1992 {
1993         MLoop *ml;
1994         MLoop *mloop = &loopstart[-mpoly->loopstart];
1995
1996         int j;
1997         for (j = 0, ml = loopstart; j < mpoly->totloop; j++, ml++) {
1998                 MLoop *ml_prev = ME_POLY_LOOP_PREV(mloop, mpoly, j);
1999                 MLoop *ml_next = ME_POLY_LOOP_NEXT(mloop, mpoly, j);
2000
2001                 float e1[3], e2[3];
2002
2003                 sub_v3_v3v3(e1, mvarray[ml_next->v].co, mvarray[ml->v].co);
2004                 sub_v3_v3v3(e2, mvarray[ml_prev->v].co, mvarray[ml->v].co);
2005
2006                 angles[j] = (float)M_PI - angle_v3v3(e1, e2);
2007         }
2008 }
2009
2010 #else /* equivalent the function above but avoid multiple subtractions + normalize */
2011
2012 void BKE_mesh_calc_poly_angles(
2013         const MPoly *mpoly, const MLoop *loopstart,
2014         const MVert *mvarray, float angles[])
2015 {
2016         float nor_prev[3];
2017         float nor_next[3];
2018
2019         int i_this = mpoly->totloop - 1;
2020         int i_next = 0;
2021
2022         sub_v3_v3v3(nor_prev, mvarray[loopstart[i_this - 1].v].co, mvarray[loopstart[i_this].v].co);
2023         normalize_v3(nor_prev);
2024
2025         while (i_next < mpoly->totloop) {
2026                 sub_v3_v3v3(nor_next, mvarray[loopstart[i_this].v].co, mvarray[loopstart[i_next].v].co);
2027                 normalize_v3(nor_next);
2028                 angles[i_this] = angle_normalized_v3v3(nor_prev, nor_next);
2029
2030                 /* step */
2031                 copy_v3_v3(nor_prev, nor_next);
2032                 i_this = i_next;
2033                 i_next++;
2034         }
2035 }
2036 #endif
2037
2038 void BKE_mesh_poly_edgehash_insert(EdgeHash *ehash, const MPoly *mp, const MLoop *mloop)
2039 {
2040         const MLoop *ml, *ml_next;
2041         int i = mp->totloop;
2042
2043         ml_next = mloop;       /* first loop */
2044         ml = &ml_next[i - 1];  /* last loop */
2045
2046         while (i-- != 0) {
2047                 BLI_edgehash_reinsert(ehash, ml->v, ml_next->v, NULL);
2048
2049                 ml = ml_next;
2050                 ml_next++;
2051         }
2052 }
2053
2054 void BKE_mesh_poly_edgebitmap_insert(unsigned int *edge_bitmap, const MPoly *mp, const MLoop *mloop)
2055 {
2056         const MLoop *ml;
2057         int i = mp->totloop;
2058
2059         ml = mloop;
2060
2061         while (i-- != 0) {
2062                 BLI_BITMAP_ENABLE(edge_bitmap, ml->e);
2063                 ml++;
2064         }
2065 }
2066
2067 /** \} */
2068
2069
2070 /* -------------------------------------------------------------------- */
2071
2072 /** \name Mesh Center Calculation
2073  * \{ */
2074
2075 bool BKE_mesh_center_median(const Mesh *me, float r_cent[3])
2076 {
2077         int i = me->totvert;
2078         const MVert *mvert;
2079         zero_v3(r_cent);
2080         for (mvert = me->mvert; i--; mvert++) {
2081                 add_v3_v3(r_cent, mvert->co);
2082         }
2083         /* otherwise we get NAN for 0 verts */
2084         if (me->totvert) {
2085                 mul_v3_fl(r_cent, 1.0f / (float)me->totvert);
2086         }
2087
2088         return (me->totvert != 0);
2089 }
2090
2091 bool BKE_mesh_center_bounds(const Mesh *me, float r_cent[3])
2092 {
2093         float min[3], max[3];
2094         INIT_MINMAX(min, max);
2095         if (BKE_mesh_minmax(me, min, max)) {
2096                 mid_v3_v3v3(r_cent, min, max);
2097                 return true;
2098         }
2099
2100         return false;
2101 }
2102
2103 bool BKE_mesh_center_centroid(const Mesh *me, float r_cent[3])
2104 {
2105         int i = me->totpoly;
2106         MPoly *mpoly;
2107         float poly_area;
2108         float total_area = 0.0f;
2109         float poly_cent[3];
2110
2111         zero_v3(r_cent);
2112
2113         /* calculate a weighted average of polygon centroids */
2114         for (mpoly = me->mpoly; i--; mpoly++) {
2115                 poly_area = mesh_calc_poly_planar_area_centroid(mpoly, me->mloop + mpoly->loopstart, me->mvert, poly_cent);
2116
2117                 madd_v3_v3fl(r_cent, poly_cent, poly_area);
2118                 total_area += poly_area;
2119         }
2120         /* otherwise we get NAN for 0 polys */
2121         if (me->totpoly) {
2122                 mul_v3_fl(r_cent, 1.0f / total_area);
2123         }
2124
2125         /* zero area faces cause this, fallback to median */
2126         if (UNLIKELY(!is_finite_v3(r_cent))) {
2127                 return BKE_mesh_center_median(me, r_cent);
2128         }
2129
2130         return (me->totpoly != 0);
2131 }
2132 /** \} */
2133
2134
2135 /* -------------------------------------------------------------------- */
2136
2137 /** \name Mesh Volume Calculation
2138  * \{ */
2139
2140 static bool mesh_calc_center_centroid_ex(
2141         const MVert *mverts, int UNUSED(mverts_num),
2142         const MLoopTri *looptri, int looptri_num,
2143         const MLoop *mloop, float r_center[3])
2144 {
2145         const MLoopTri *lt;
2146         float totweight;
2147         int i;
2148         
2149         zero_v3(r_center);
2150         
2151         if (looptri_num == 0)
2152                 return false;
2153         
2154         totweight = 0.0f;
2155         for (i = 0, lt = looptri; i < looptri_num; i++, lt++) {
2156                 const MVert *v1 = &mverts[mloop[lt->tri[0]].v];
2157                 const MVert *v2 = &mverts[mloop[lt->tri[1]].v];
2158                 const MVert *v3 = &mverts[mloop[lt->tri[2]].v];
2159                 float area;
2160                 
2161                 area = area_tri_v3(v1->co, v2->co, v3->co);
2162                 madd_v3_v3fl(r_center, v1->co, area);
2163                 madd_v3_v3fl(r_center, v2->co, area);
2164                 madd_v3_v3fl(r_center, v3->co, area);
2165                 totweight += area;
2166         }
2167         if (totweight == 0.0f)
2168                 return false;
2169         
2170         mul_v3_fl(r_center, 1.0f / (3.0f * totweight));
2171         
2172         return true;
2173 }
2174
2175 /**
2176  * Calculate the volume and center.
2177  *
2178  * \param r_volume: Volume (unsigned).
2179  * \param r_center: Center of mass.
2180  */
2181 void BKE_mesh_calc_volume(
2182         const MVert *mverts, const int mverts_num,
2183         const MLoopTri *looptri, const int looptri_num,
2184         const MLoop *mloop,
2185         float *r_volume, float r_center[3])
2186 {
2187         const MLoopTri *lt;
2188         float center[3];
2189         float totvol;
2190         int i;
2191         
2192         if (r_volume)
2193                 *r_volume = 0.0f;
2194         if (r_center)
2195                 zero_v3(r_center);
2196         
2197         if (looptri_num == 0)
2198                 return;
2199         
2200         if (!mesh_calc_center_centroid_ex(mverts, mverts_num, looptri, looptri_num, mloop, center))
2201                 return;
2202         
2203         totvol = 0.0f;
2204
2205         for (i = 0, lt = looptri; i < looptri_num; i++, lt++) {
2206                 const MVert *v1 = &mverts[mloop[lt->tri[0]].v];
2207                 const MVert *v2 = &mverts[mloop[lt->tri[1]].v];
2208                 const MVert *v3 = &mverts[mloop[lt->tri[2]].v];
2209                 float vol;
2210                 
2211                 vol = volume_tetrahedron_signed_v3(center, v1->co, v2->co, v3->co);
2212                 if (r_volume) {
2213                         totvol += vol;
2214                 }
2215                 if (r_center) {
2216                         /* averaging factor 1/3 is applied in the end */
2217                         madd_v3_v3fl(r_center, v1->co, vol);
2218                         madd_v3_v3fl(r_center, v2->co, vol);
2219                         madd_v3_v3fl(r_center, v3->co, vol);
2220                 }
2221         }
2222         
2223         /* Note: Depending on arbitrary centroid position,
2224          * totvol can become negative even for a valid mesh.
2225          * The true value is always the positive value.
2226          */
2227         if (r_volume) {
2228                 *r_volume = fabsf(totvol);
2229         }
2230         if (r_center) {
2231                 /* Note: Factor 1/3 is applied once for all vertices here.
2232                  * This also automatically negates the vector if totvol is negative.
2233                  */
2234                 if (totvol != 0.0f)
2235                         mul_v3_fl(r_center, (1.0f / 3.0f) / totvol);
2236         }
2237 }
2238
2239
2240 /* -------------------------------------------------------------------- */
2241
2242 /** \name NGon Tessellation (NGon/Tessface Conversion)
2243  * \{ */
2244
2245 /**
2246  * Convert a triangle or quadrangle of loop/poly data to tessface data
2247  */
2248 void BKE_mesh_loops_to_mface_corners(
2249         CustomData *fdata, CustomData *ldata,
2250         CustomData *pdata, unsigned int lindex[4], int findex,
2251         const int polyindex,
2252         const int mf_len, /* 3 or 4 */
2253
2254         /* cache values to avoid lookups every time */
2255         const int numTex, /* CustomData_number_of_layers(pdata, CD_MTEXPOLY) */
2256         const int numCol, /* CustomData_number_of_layers(ldata, CD_MLOOPCOL) */
2257         const bool hasPCol, /* CustomData_has_layer(ldata, CD_PREVIEW_MLOOPCOL) */
2258         const bool hasOrigSpace, /* CustomData_has_layer(ldata, CD_ORIGSPACE_MLOOP) */
2259         const bool hasLNor /* CustomData_has_layer(ldata, CD_NORMAL) */
2260 )
2261 {
2262         MTFace *texface;
2263         MTexPoly *texpoly;
2264         MCol *mcol;
2265         MLoopCol *mloopcol;
2266         MLoopUV *mloopuv;
2267         int i, j;
2268
2269         for (i = 0; i < numTex; i++) {
2270                 texface = CustomData_get_n(fdata, CD_MTFACE, findex, i);
2271                 texpoly = CustomData_get_n(pdata, CD_MTEXPOLY, polyindex, i);
2272
2273                 ME_MTEXFACE_CPY(texface, texpoly);
2274
2275                 for (j = 0; j < mf_len; j++) {
2276                         mloopuv = CustomData_get_n(ldata, CD_MLOOPUV, (int)lindex[j], i);
2277                         copy_v2_v2(texface->uv[j], mloopuv->uv);
2278                 }
2279         }
2280
2281         for (i = 0; i < numCol; i++) {
2282                 mcol = CustomData_get_n(fdata, CD_MCOL, findex, i);
2283
2284                 for (j = 0; j < mf_len; j++) {
2285                         mloopcol = CustomData_get_n(ldata, CD_MLOOPCOL, (int)lindex[j], i);
2286                         MESH_MLOOPCOL_TO_MCOL(mloopcol, &mcol[j]);
2287                 }
2288         }
2289
2290         if (hasPCol) {
2291                 mcol = CustomData_get(fdata,  findex, CD_PREVIEW_MCOL);
2292
2293                 for (j = 0; j < mf_len; j++) {
2294                         mloopcol = CustomData_get(ldata, (int)lindex[j], CD_PREVIEW_MLOOPCOL);
2295                         MESH_MLOOPCOL_TO_MCOL(mloopcol, &mcol[j]);
2296                 }
2297         }
2298
2299         if (hasOrigSpace) {
2300                 OrigSpaceFace *of = CustomData_get(fdata, findex, CD_ORIGSPACE);
2301                 OrigSpaceLoop *lof;
2302
2303                 for (j = 0; j < mf_len; j++) {
2304                         lof = CustomData_get(ldata, (int)lindex[j], CD_ORIGSPACE_MLOOP);
2305                         copy_v2_v2(of->uv[j], lof->uv);
2306                 }
2307         }
2308
2309         if (hasLNor) {
2310                 short (*tlnors)[3] = CustomData_get(fdata, findex, CD_TESSLOOPNORMAL);
2311
2312                 for (j = 0; j < mf_len; j++) {
2313                         normal_float_to_short_v3(tlnors[j], CustomData_get(ldata, (int)lindex[j], CD_NORMAL));
2314                 }
2315         }
2316 }
2317
2318 /**
2319  * Convert all CD layers from loop/poly to tessface data.
2320  *
2321  * \param loopindices is an array of an int[4] per tessface, mapping tessface's verts to loops indices.
2322  *
2323  * \note when mface is not NULL, mface[face_index].v4 is used to test quads, else, loopindices[face_index][3] is used.
2324  */
2325 void BKE_mesh_loops_to_tessdata(CustomData *fdata, CustomData *ldata, CustomData *pdata, MFace *mface,
2326                                 int *polyindices, unsigned int (*loopindices)[4], const int num_faces)
2327 {
2328         /* Note: performances are sub-optimal when we get a NULL mface, we could be ~25% quicker with dedicated code...
2329          *       Issue is, unless having two different functions with nearly the same code, there's not much ways to solve
2330          *       this. Better imho to live with it for now. :/ --mont29
2331          */
2332         const int numTex = CustomData_number_of_layers(pdata, CD_MTEXPOLY);
2333         const int numCol = CustomData_number_of_layers(ldata, CD_MLOOPCOL);
2334         const bool hasPCol = CustomData_has_layer(ldata, CD_PREVIEW_MLOOPCOL);
2335         const bool hasOrigSpace = CustomData_has_layer(ldata, CD_ORIGSPACE_MLOOP);
2336         const bool hasLoopNormal = CustomData_has_layer(ldata, CD_NORMAL);
2337         const bool hasLoopTangent = CustomData_has_layer(ldata, CD_TANGENT);
2338         int findex, i, j;
2339         const int *pidx;
2340         unsigned int (*lidx)[4];
2341
2342         for (i = 0; i < numTex; i++) {
2343                 MTFace *texface = CustomData_get_layer_n(fdata, CD_MTFACE, i);
2344                 MTexPoly *texpoly = CustomData_get_layer_n(pdata, CD_MTEXPOLY, i);
2345                 MLoopUV *mloopuv = CustomData_get_layer_n(ldata, CD_MLOOPUV, i);
2346
2347                 for (findex = 0, pidx = polyindices, lidx = loopindices;
2348                      findex < num_faces;
2349                      pidx++, lidx++, findex++, texface++)
2350                 {
2351                         ME_MTEXFACE_CPY(texface, &texpoly[*pidx]);
2352
2353                         for (j = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3; j--;) {
2354                                 copy_v2_v2(texface->uv[j], mloopuv[(*lidx)[j]].uv);
2355                         }
2356                 }
2357         }
2358
2359         for (i = 0; i < numCol; i++) {
2360                 MCol (*mcol)[4] = CustomData_get_layer_n(fdata, CD_MCOL, i);
2361                 MLoopCol *mloopcol = CustomData_get_layer_n(ldata, CD_MLOOPCOL, i);
2362
2363                 for (findex = 0, lidx = loopindices; findex < num_faces; lidx++, findex++, mcol++) {
2364                         for (j = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3; j--;) {
2365                                 MESH_MLOOPCOL_TO_MCOL(&mloopcol[(*lidx)[j]], &(*mcol)[j]);
2366                         }
2367                 }
2368         }
2369
2370         if (hasPCol) {
2371                 MCol (*mcol)[4] = CustomData_get_layer(fdata, CD_PREVIEW_MCOL);
2372                 MLoopCol *mloopcol = CustomData_get_layer(ldata, CD_PREVIEW_MLOOPCOL);
2373
2374                 for (findex = 0, lidx = loopindices; findex < num_faces; lidx++, findex++, mcol++) {
2375                         for (j = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3; j--;) {
2376                                 MESH_MLOOPCOL_TO_MCOL(&mloopcol[(*lidx)[j]], &(*mcol)[j]);
2377                         }
2378                 }
2379         }
2380
2381         if (hasOrigSpace) {
2382                 OrigSpaceFace *of = CustomData_get_layer(fdata, CD_ORIGSPACE);
2383                 OrigSpaceLoop *lof = CustomData_get_layer(ldata, CD_ORIGSPACE_MLOOP);
2384
2385                 for (findex = 0, lidx = loopindices; findex < num_faces; lidx++, findex++, of++) {
2386                         for (j = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3; j--;) {
2387                                 copy_v2_v2(of->uv[j], lof[(*lidx)[j]].uv);
2388                         }
2389                 }
2390         }
2391
2392         if (hasLoopNormal) {
2393                 short (*fnors)[4][3] = CustomData_get_layer(fdata, CD_TESSLOOPNORMAL);
2394                 float (*lnors)[3] = CustomData_get_layer(ldata, CD_NORMAL);
2395
2396                 for (findex = 0, lidx = loopindices; findex < num_faces; lidx++, findex++, fnors++) {
2397                         for (j = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3; j--;) {
2398                                 normal_float_to_short_v3((*fnors)[j], lnors[(*lidx)[j]]);
2399                         }
2400                 }
2401         }
2402
2403         if (hasLoopTangent) {
2404                 /* need to do for all uv maps at some point */
2405                 float (*ftangents)[4] = CustomData_get_layer(fdata, CD_TANGENT);
2406                 float (*ltangents)[4] = CustomData_get_layer(ldata, CD_TANGENT);
2407
2408                 for (findex = 0, pidx = polyindices, lidx = loopindices;
2409                      findex < num_faces;
2410                      pidx++, lidx++, findex++)
2411                 {
2412                         int nverts = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3;
2413                         for (j = nverts; j--;) {
2414                                 copy_v4_v4(ftangents[findex * 4 + j], ltangents[(*lidx)[j]]);
2415                         }
2416                 }
2417         }
2418 }
2419
2420 void BKE_mesh_tangent_loops_to_tessdata(CustomData *fdata, CustomData *ldata, MFace *mface,
2421                                         int *polyindices, unsigned int (*loopindices)[4], const int num_faces)
2422 {
2423         /* Note: performances are sub-optimal when we get a NULL mface, we could be ~25% quicker with dedicated code...
2424          *       Issue is, unless having two different functions with nearly the same code, there's not much ways to solve
2425          *       this. Better imho to live with it for now. :/ --mont29
2426          */
2427         const bool hasLoopTangent = CustomData_has_layer(ldata, CD_TANGENT);
2428         int findex, j;
2429         const int *pidx;
2430         unsigned int (*lidx)[4];
2431
2432         if (hasLoopTangent) {
2433                 /* need to do for all uv maps at some point */
2434                 float (*ftangents)[4] = CustomData_get_layer(fdata, CD_TANGENT);
2435                 float (*ltangents)[4] = CustomData_get_layer(ldata, CD_TANGENT);
2436
2437                 for (findex = 0, pidx = polyindices, lidx = loopindices;
2438                      findex < num_faces;
2439                      pidx++, lidx++, findex++)
2440                 {
2441                         int nverts = (mface ? mface[findex].v4 : (*lidx)[3]) ? 4 : 3;
2442                         for (j = nverts; j--;) {
2443                                 copy_v4_v4(ftangents[findex * 4 + j], ltangents[(*lidx)[j]]);
2444                         }
2445                 }
2446         }
2447 }
2448
2449 /**
2450  * Recreate tessellation.
2451  *
2452  * \param do_face_nor_copy: Controls whether the normals from the poly are copied to the tessellated faces.
2453  *
2454  * \return number of tessellation faces.
2455  */
2456 int BKE_mesh_recalc_tessellation(
2457         CustomData *fdata, CustomData *ldata, CustomData *pdata,
2458         MVert *mvert,
2459         int totface, int totloop, int totpoly,
2460         const bool do_face_nor_copy)
2461 {
2462         /* use this to avoid locking pthread for _every_ polygon
2463          * and calling the fill function */
2464
2465 #define USE_TESSFACE_SPEEDUP
2466 #define USE_TESSFACE_QUADS  /* NEEDS FURTHER TESTING */
2467
2468 /* We abuse MFace->edcode to tag quad faces. See below for details. */
2469 #define TESSFACE_IS_QUAD 1
2470
2471         const int looptri_num = poly_to_tri_count(totpoly, totloop);
2472
2473         MPoly *mp, *mpoly;
2474         MLoop *ml, *mloop;
2475         MFace *mface, *mf;
2476         MemArena *arena = NULL;
2477         int *mface_to_poly_map;
2478         unsigned int (*lindices)[4];
2479         int poly_index, mface_index;
2480         unsigned int j;
2481
2482         mpoly = CustomData_get_layer(pdata, CD_MPOLY);
2483         mloop = CustomData_get_layer(ldata, CD_MLOOP);
2484
2485         /* allocate the length of totfaces, avoid many small reallocs,
2486          * if all faces are tri's it will be correct, quads == 2x allocs */
2487         /* take care. we are _not_ calloc'ing so be sure to initialize each field */
2488         mface_to_poly_map = MEM_mallocN(sizeof(*mface_to_poly_map) * (size_t)looptri_num, __func__);
2489         mface             = MEM_mallocN(sizeof(*mface) *             (size_t)looptri_num, __func__);
2490         lindices          = MEM_mallocN(sizeof(*lindices) *          (size_t)looptri_num, __func__);
2491
2492         mface_index = 0;
2493         mp = mpoly;
2494         for (poly_index = 0; poly_index < totpoly; poly_index++, mp++) {
2495                 const unsigned int mp_loopstart = (unsigned int)mp->loopstart;
2496                 const unsigned int mp_totloop = (unsigned int)mp->totloop;
2497                 unsigned int l1, l2, l3, l4;
2498                 unsigned int *lidx;
2499                 if (mp_totloop < 3) {
2500                         /* do nothing */
2501                 }
2502
2503 #ifdef USE_TESSFACE_SPEEDUP
2504
2505 #define ML_TO_MF(i1, i2, i3)                                                  \
2506                 mface_to_poly_map[mface_index] = poly_index;                          \
2507                 mf = &mface[mface_index];                                             \
2508                 lidx = lindices[mface_index];                                         \
2509                 /* set loop indices, transformed to vert indices later */             \
2510                 l1 = mp_loopstart + i1;                                               \
2511                 l2 = mp_loopstart + i2;                                               \
2512                 l3 = mp_loopstart + i3;                                               \
2513                 mf->v1 = mloop[l1].v;                                                 \
2514                 mf->v2 = mloop[l2].v;                                                 \
2515                 mf->v3 = mloop[l3].v;                                                 \
2516                 mf->v4 = 0;                                                           \
2517                 lidx[0] = l1;                                                         \
2518                 lidx[1] = l2;                                                         \
2519                 lidx[2] = l3;                                                         \
2520                 lidx[3] = 0;                                                          \
2521                 mf->mat_nr = mp->mat_nr;                                              \
2522                 mf->flag = mp->flag;                                                  \
2523                 mf->edcode = 0;                                                       \
2524                 (void)0
2525
2526 /* ALMOST IDENTICAL TO DEFINE ABOVE (see EXCEPTION) */
2527 #define ML_TO_MF_QUAD()                                                       \
2528                 mface_to_poly_map[mface_index] = poly_index;                          \
2529                 mf = &mface[mface_index];                                             \
2530                 lidx = lindices[mface_index];                                         \
2531                 /* set loop indices, transformed to vert indices later */             \
2532                 l1 = mp_loopstart + 0; /* EXCEPTION */                                \
2533                 l2 = mp_loopstart + 1; /* EXCEPTION */                                \
2534                 l3 = mp_loopstart + 2; /* EXCEPTION */                                \
2535                 l4 = mp_loopstart + 3; /* EXCEPTION */                                \
2536                 mf->v1 = mloop[l1].v;                                                 \
2537                 mf->v2 = mloop[l2].v;                                                 \
2538                 mf->v3 = mloop[l3].v;                                                 \
2539                 mf->v4 = mloop[l4].v;                                                 \
2540                 lidx[0] = l1;                                                         \
2541                 lidx[1] = l2;                                                         \
2542                 lidx[2] = l3;                                                         \
2543                 lidx[3] = l4;                                                         \
2544                 mf->mat_nr = mp->mat_nr;                                              \
2545                 mf->flag = mp->flag;                                                  \
2546                 mf->edcode = TESSFACE_IS_QUAD;                                        \
2547                 (void)0
2548
2549
2550                 else if (mp_totloop == 3) {
2551                         ML_TO_MF(0, 1, 2);
2552                         mface_index++;
2553                 }
2554                 else if (mp_totloop == 4) {
2555 #ifdef USE_TESSFACE_QUADS
2556                         ML_TO_MF_QUAD();
2557                         mface_index++;
2558 #else
2559                         ML_TO_MF(0, 1, 2);
2560                         mface_index++;
2561                         ML_TO_MF(0, 2, 3);
2562                         mface_index++;
2563 #endif
2564                 }
2565 #endif /* USE_TESSFACE_SPEEDUP */
2566                 else {
2567                         const float *co_curr, *co_prev;
2568
2569                         float normal[3];
2570
2571                         float axis_mat[3][3];
2572                         float (*projverts)[2];
2573                         unsigned int (*tris)[3];
2574
2575                         const unsigned int totfilltri = mp_totloop - 2;
2576
2577                         if (UNLIKELY(arena == NULL)) {
2578                                 arena = BLI_memarena_new(BLI_MEMARENA_STD_BUFSIZE, __func__);
2579                         }
2580
2581                         tris = BLI_memarena_alloc(arena, sizeof(*tris) * (size_t)totfilltri);
2582                         projverts = BLI_memarena_alloc(arena, sizeof(*projverts) * (size_t)mp_totloop);
2583
2584                         zero_v3(normal);
2585
2586                         /* calc normal, flipped: to get a positive 2d cross product */
2587                         ml = mloop + mp_loopstart;
2588                         co_prev = mvert[ml[mp_totloop - 1].v].co;
2589                         for (j = 0; j < mp_totloop; j++, ml++) {
2590                                 co_curr = mvert[ml->v].co;
2591                                 add_newell_cross_v3_v3v3(normal, co_prev, co_curr);
2592                                 co_prev = co_curr;
2593                         }
2594                         if (UNLIKELY(normalize_v3(normal) == 0.0f)) {
2595                                 normal[2] = 1.0f;
2596                         }
2597
2598                         /* project verts to 2d */
2599                         axis_dominant_v3_to_m3_negate(axis_mat, normal);
2600
2601                         ml = mloop + mp_loopstart;
2602                         for (j = 0; j < mp_totloop; j++, ml++) {
2603                                 mul_v2_m3v3(projverts[j], axis_mat, mvert[ml->v].co);
2604                         }
2605
2606                         BLI_polyfill_calc_arena((const float (*)[2])projverts, mp_totloop, 1, tris, arena);
2607
2608                         /* apply fill */
2609                         for (j = 0; j < totfilltri; j++) {
2610                                 unsigned int *tri = tris[j];
2611                                 lidx = lindices[mface_index];
2612
2613                                 mface_to_poly_map[mface_index] = poly_index;
2614                                 mf = &mface[mface_index];
2615
2616                                 /* set loop indices, transformed to vert indices later */
2617                                 l1 = mp_loopstart + tri[0];
2618                                 l2 = mp_loopstart + tri[1];
2619                                 l3 = mp_loopstart + tri[2];
2620
2621                                 mf->v1 = mloop[l1].v;
2622                                 mf->v2 = mloop[l2].v;
2623                                 mf->v3 = mloop[l3].v;
2624                                 mf->v4 = 0;
2625
2626                                 lidx[0] = l1;
2627                                 lidx[1] = l2;
2628                                 lidx[2] = l3;
2629                                 lidx[3] = 0;
2630
2631                                 mf->mat_nr = mp->mat_nr;
2632                                 mf->flag = mp->flag;
2633                                 mf->edcode = 0;
2634
2635                                 mface_index++;
2636                         }
2637
2638                         BLI_memarena_clear(arena);
2639                 }
2640         }
2641
2642         if (arena) {
2643                 BLI_memarena_free(arena);
2644                 arena = NULL;
2645         }
2646
2647         CustomData_free(fdata, totface);
2648         totface = mface_index;
2649
2650         BLI_assert(totface <= looptri_num);
2651
2652         /* not essential but without this we store over-alloc'd memory in the CustomData layers */
2653         if (LIKELY(looptri_num != totface)) {
2654                 mface = MEM_reallocN(mface, sizeof(*mface) * (size_t)totface);
2655                 mface_to_poly_map = MEM_reallocN(mface_to_poly_map, sizeof(*mface_to_poly_map) * (size_t)totface);
2656         }
2657
2658         CustomData_add_layer(fdata, CD_MFACE, CD_ASSIGN, mface, totface);
2659
2660         /* CD_ORIGINDEX will contain an array of indices from tessfaces to the polygons
2661          * they are directly tessellated from */
2662         CustomData_add_layer(fdata, CD_ORIGINDEX, CD_ASSIGN, mface_to_poly_map, totface);
2663         CustomData_from_bmeshpoly(fdata, pdata, ldata, totface);
2664
2665         if (do_face_nor_copy) {
2666                 /* If polys have a normals layer, copying that to faces can help
2667                  * avoid the need to recalculate normals later */
2668                 if (CustomData_has_layer(pdata, CD_NORMAL)) {
2669                         float (*pnors)[3] = CustomData_get_layer(pdata, CD_NORMAL);
2670                         float (*fnors)[3] = CustomData_add_layer(fdata, CD_NORMAL, CD_CALLOC, NULL, totface);
2671                         for (mface_index = 0; mface_index < totface; mface_index++) {
2672                                 copy_v3_v3(fnors[mface_index], pnors[mface_to_poly_map[mface_index]]);
2673                         }
2674                 }
2675         }
2676
2677         /* NOTE: quad detection issue - fourth vertidx vs fourth loopidx:
2678          * Polygons take care of their loops ordering, hence not of their vertices ordering.
2679          * Currently, our tfaces' fourth vertex index might be 0 even for a quad. However, we know our fourth loop index is
2680          * never 0 for quads (because they are sorted for polygons, and our quads are still mere copies of their polygons).
2681          * So we pass NULL as MFace pointer, and BKE_mesh_loops_to_tessdata will use the fourth loop index as quad test.
2682          * ...
2683          */
2684         BKE_mesh_loops_to_tessdata(fdata, ldata, pdata, NULL, mface_to_poly_map, lindices, totface);
2685
2686         /* NOTE: quad detection issue - fourth vertidx vs fourth loopidx:
2687          * ...However, most TFace code uses 'MFace->v4 == 0' test to check whether it is a tri or quad.
2688          * test_index_face() will check this and rotate the tessellated face if needed.
2689          */
2690 #ifdef USE_TESSFACE_QUADS
2691         mf = mface;
2692         for (mface_index = 0; mface_index < totface; mface_index++, mf++) {
2693                 if (mf->edcode == TESSFACE_IS_QUAD) {
2694                         test_index_face(mf, fdata, mface_index, 4);
2695                         mf->edcode = 0;
2696                 }
2697         }
2698 #endif
2699
2700         MEM_freeN(lindices);
2701
2702         return totface;
2703
2704 #undef USE_TESSFACE_SPEEDUP
2705 #undef USE_TESSFACE_QUADS
2706
2707 #undef ML_TO_MF
2708 #undef ML_TO_MF_QUAD
2709
2710 }
2711
2712 /**
2713  * Calculate tessellation into #MLoopTri which exist only for this purpose.
2714  */
2715 void BKE_mesh_recalc_looptri(
2716         const MLoop *mloop, const MPoly *mpoly,
2717         const MVert *mvert,
2718         int totloop, int totpoly,
2719         MLoopTri *mlooptri)
2720 {
2721         /* use this to avoid locking pthread for _every_ polygon
2722          * and calling the fill function */
2723
2724 #define USE_TESSFACE_SPEEDUP
2725
2726         const MPoly *mp;
2727         const MLoop *ml;
2728         MLoopTri *mlt;
2729         MemArena *arena = NULL;
2730         int poly_index, mlooptri_index;
2731         unsigned int j;
2732
2733         mlooptri_index = 0;
2734         mp = mpoly;
2735         for (poly_index = 0; poly_index < totpoly; poly_index++, mp++) {
2736                 const unsigned int mp_loopstart = (unsigned int)mp->loopstart;
2737                 const unsigned int mp_totloop = (unsigned int)mp->totloop;
2738                 unsigned int l1, l2, l3;
2739                 if (mp_totloop < 3) {
2740                         /* do nothing */
2741                 }
2742
2743 #ifdef USE_TESSFACE_SPEEDUP
2744
2745 #define ML_TO_MLT(i1, i2, i3)  { \
2746                         mlt = &mlooptri[mlooptri_index]; \
2747                         l1 = mp_loopstart + i1; \
2748                         l2 = mp_loopstart + i2; \
2749                         l3 = mp_loopstart + i3; \
2750                         ARRAY_SET_ITEMS(mlt->tri, l1, l2, l3); \
2751                         mlt->poly = (unsigned int)poly_index; \
2752                 } ((void)0)
2753
2754                 else if (mp_totloop == 3) {
2755                         ML_TO_MLT(0, 1, 2);
2756                         mlooptri_index++;
2757                 }
2758                 else if (mp_totloop == 4) {
2759                         ML_TO_MLT(0, 1, 2);
2760                         mlooptri_index++;
2761                         ML_TO_MLT(0, 2, 3);
2762                         mlooptri_index++;
2763                 }
2764 #endif /* USE_TESSFACE_SPEEDUP */
2765                 else {
2766                         const float *co_curr, *co_prev;
2767
2768                         float normal[3];
2769
2770                         float axis_mat[3][3];
2771                         float (*projverts)[2];
2772                         unsigned int (*tris)[3];
2773
2774                         const unsigned int totfilltri = mp_totloop - 2;
2775
2776                         if (UNLIKELY(arena == NULL)) {
2777                                 arena = BLI_memarena_new(BLI_MEMARENA_STD_BUFSIZE, __func__);
2778                         }
2779
2780                         tris = BLI_memarena_alloc(arena, sizeof(*tris) * (size_t)totfilltri);
2781                         projverts = BLI_memarena_alloc(arena, sizeof(*projverts) * (size_t)mp_totloop);
2782
2783                         zero_v3(normal);
2784
2785                         /* calc normal, flipped: to get a positive 2d cross product */
2786                         ml = mloop + mp_loopstart;
2787                         co_prev = mvert[ml[mp_totloop - 1].v].co;
2788                         for (j = 0; j < mp_totloop; j++, ml++) {
2789                                 co_curr = mvert[ml->v].co;
2790                                 add_newell_cross_v3_v3v3(normal, co_prev, co_curr);
2791                                 co_prev = co_curr;
2792                         }
2793                         if (UNLIKELY(normalize_v3(normal) == 0.0f)) {
2794                                 normal[2] = 1.0f;
2795                         }
2796
2797                         /* project verts to 2d */
2798                         axis_dominant_v3_to_m3_negate(axis_mat, normal);
2799
2800                         ml = mloop + mp_loopstart;
2801                         for (j = 0; j < mp_totloop; j++, ml++) {
2802                                 mul_v2_m3v3(projverts[j], axis_mat, mvert[ml->v].co);
2803                         }
2804
2805                         BLI_polyfill_calc_arena((const float (*)[2])projverts, mp_totloop, 1, tris, arena);
2806
2807                         /* apply fill */
2808                         for (j = 0; j < totfilltri; j++) {
2809                                 unsigned int *tri = tris[j];
2810
2811                                 mlt = &mlooptri[mlooptri_index];
2812
2813                                 /* set loop indices, transformed to vert indices later */
2814                                 l1 = mp_loopstart + tri[0];
2815                                 l2 = mp_loopstart + tri[1];
2816                                 l3 = mp_loopstart + tri[2];
2817
2818                                 ARRAY_SET_ITEMS(mlt->tri, l1, l2, l3);
2819                                 mlt->poly = (unsigned int)poly_index;
2820
2821                                 mlooptri_index++;
2822                         }
2823
2824                         BLI_memarena_clear(arena);
2825                 }
2826         }
2827
2828         if (arena) {
2829                 BLI_memarena_free(arena);
2830                 arena = NULL;
2831         }
2832
2833         BLI_assert(mlooptri_index == poly_to_tri_count(totpoly, totloop));
2834         UNUSED_VARS_NDEBUG(totloop);
2835
2836 #undef USE_TESSFACE_SPEEDUP
2837 #undef ML_TO_MLT
2838 }
2839
2840 /* -------------------------------------------------------------------- */
2841
2842
2843 #ifdef USE_BMESH_SAVE_AS_COMPAT
2844
2845 /**
2846  * This function recreates a tessellation.
2847  * returns number of tessellation faces.
2848  *
2849  * for forwards compat only quad->tri polys to mface, skip ngons.
2850  */
2851 int BKE_mesh_mpoly_to_mface(struct CustomData *fdata, struct CustomData *ldata,
2852                             struct CustomData *pdata, int totface, int UNUSED(totloop), int totpoly)
2853 {
2854         MLoop *mloop;
2855
2856         unsigned int lindex[4];
2857         int i;
2858         int k;
2859
2860         MPoly *mp, *mpoly;
2861         MFace *mface, *mf;
2862
2863         const int numTex = CustomData_number_of_layers(pdata, CD_MTEXPOLY);
2864         const int numCol = CustomData_number_of_layers(ldata, CD_MLOOPCOL);
2865         const bool hasPCol = CustomData_has_layer(ldata, CD_PREVIEW_MLOOPCOL);
2866         const bool hasOrigSpace = CustomData_has_layer(ldata, CD_ORIGSPACE_MLOOP);
2867         const bool hasLNor = CustomData_has_layer(ldata, CD_NORMAL);
2868
2869         /* over-alloc, ngons will be skipped */
2870         mface = MEM_mallocN(sizeof(*mface) * (size_t)totpoly, __func__);
2871
2872         mpoly = CustomData_get_layer(pdata, CD_MPOLY);
2873         mloop = CustomData_get_layer(ldata, CD_MLOOP);
2874
2875         mp = mpoly;
2876         k = 0;
2877         for (i = 0; i < totpoly; i++, mp++) {
2878                 if (ELEM(mp->totloop, 3, 4)) {
2879                         const unsigned int mp_loopstart = (unsigned int)mp->loopstart;
2880                         mf = &mface[k];
2881
2882                         mf->mat_nr = mp->mat_nr;
2883                         mf->flag = mp->flag;
2884
2885                         mf->v1 = mp_loopstart + 0;
2886                         mf->v2 = mp_loopstart + 1;
2887                         mf->v3 = mp_loopstart + 2;
2888                         mf->v4 = (mp->totloop == 4) ? (mp_loopstart + 3) : 0;
2889
2890                         /* abuse edcode for temp storage and clear next loop */
2891                         mf->edcode = (char)mp->totloop; /* only ever 3 or 4 */
2892
2893                         k++;
2894                 }
2895         }
2896
2897         CustomData_free(fdata, totface);
2898
2899         totface = k;
2900
2901         CustomData_add_layer(fdata, CD_MFACE, CD_ASSIGN, mface, totface);
2902
2903         CustomData_from_bmeshpoly(fdata, pdata, ldata, totface);
2904
2905         mp = mpoly;
2906         k = 0;
2907         for (i = 0; i < totpoly; i++, mp++) {
2908                 if (ELEM(mp->totloop, 3, 4)) {
2909                         mf = &mface[k];
2910
2911                         if (mf->edcode == 3) {
2912                                 /* sort loop indices to ensure winding is correct */
2913                                 /* NO SORT - looks like we can skip this */
2914
2915                                 lindex[0] = mf->v1;
2916                                 lindex[1] = mf->v2;
2917                                 lindex[2] = mf->v3;
2918                                 lindex[3] = 0; /* unused */
2919
2920                                 /* transform loop indices to vert indices */
2921                                 mf->v1 = mloop[mf->v1].v;
2922                                 mf->v2 = mloop[mf->v2].v;
2923                                 mf->v3 = mloop[mf->v3].v;
2924
2925                                 BKE_mesh_loops_to_mface_corners(fdata, ldata, pdata,
2926                                                                 lindex, k, i, 3,
2927                                                                 numTex, numCol, hasPCol, hasOrigSpace, hasLNor);
2928                                 test_index_face(mf, fdata, k, 3);
2929                         }
2930                         else {
2931                                 /* sort loop indices to ensure winding is correct */
2932                                 /* NO SORT - looks like we can skip this */
2933
2934                                 lindex[0] = mf->v1;
2935                                 lindex[1] = mf->v2;
2936                                 lindex[2] = mf->v3;
2937                                 lindex[3] = mf->v4;
2938
2939                                 /* transform loop indices to vert indices */
2940                                 mf->v1 = mloop[mf->v1].v;
2941                                 mf->v2 = mloop[mf->v2].v;
2942                                 mf->v3 = mloop[mf->v3].v;
2943                                 mf->v4 = mloop[mf->v4].v;
2944
2945                                 BKE_mesh_loops_to_mface_corners(fdata, ldata, pdata,
2946                                                                 lindex, k, i, 4,
2947                                                                 numTex, numCol, hasPCol, hasOrigSpace, hasLNor);
2948                                 test_index_face(mf, fdata, k, 4);
2949                         }
2950
2951                         mf->edcode = 0;
2952
2953                         k++;
2954                 }
2955         }
2956
2957         return k;
2958 }
2959 #endif /* USE_BMESH_SAVE_AS_COMPAT */
2960
2961
2962 static void bm_corners_to_loops_ex(ID *id, CustomData *fdata, CustomData *ldata, CustomData *pdata,
2963                                    MFace *mface, int totloop, int findex, int loopstart, int numTex, int numCol)
2964 {
2965         MTFace *texface;
2966         MTexPoly *texpoly;
2967         MCol *mcol;
2968         MLoopCol *mloopcol;
2969         MLoopUV *mloopuv;
2970         MFace *mf;
2971         int i;
2972
2973         mf = mface + findex;
2974
2975         for (i = 0; i < numTex; i++) {
2976                 texface = CustomData_get_n(fdata, CD_MTFACE, findex, i);
2977                 texpoly = CustomData_get_n(pdata, CD_MTEXPOLY, findex, i);
2978
2979                 ME_MTEXFACE_CPY(texpoly, texface);
2980
2981                 mloopuv = CustomData_get_n(ldata, CD_MLOOPUV, loopstart, i);
2982                 copy_v2_v2(mloopuv->uv, texface->uv[0]); mloopuv++;
2983                 copy_v2_v2(mloopuv->uv, texface->uv[1]); mloopuv++;
2984                 copy_v2_v2(mloopuv->uv, texface->uv[2]); mloopuv++;
2985
2986                 if (mf->v4) {
2987                         copy_v2_v2(mloopuv->uv, texface->uv[3]); mloopuv++;
2988                 }
2989         }
2990
2991         for (i = 0; i < numCol; i++) {
2992                 mloopcol = CustomData_get_n(ldata, CD_MLOOPCOL, loopstart, i);
2993                 mcol = CustomData_get_n(fdata, CD_MCOL, findex, i);
2994
2995                 MESH_MLOOPCOL_FROM_MCOL(mloopcol, &mcol[0]); mloopcol++;
2996                 MESH_MLOOPCOL_FROM_MCOL(mloopcol, &mcol[1]); mloopcol++;
2997                 MESH_MLOOPCOL_FROM_MCOL(mloopcol, &mcol[2]); mloopcol++;
2998                 if (mf->v4) {
2999                         MESH_MLOOPCOL_FROM_MCOL(mloopcol, &mcol[3]); mloopcol++;
3000                 }
3001         }
3002
3003         if (CustomData_has_layer(fdata, CD_TESSLOOPNORMAL)) {
3004                 float (*lnors)[3] = CustomData_get(ldata, loopstart, CD_NORMAL);
3005                 short (*tlnors)[3] = CustomData_get(fdata, findex, CD_TESSLOOPNORMAL);
3006                 const int max = mf->v4 ? 4 : 3;
3007
3008                 for (i = 0; i < max; i++, lnors++, tlnors++) {
3009                         normal_short_to_float_v3(*lnors, *tlnors);
3010                 }
3011         }
3012
3013         if (CustomData_has_layer(fdata, CD_MDISPS)) {
3014                 MDisps *ld = CustomData_get(ldata, loopstart, CD_MDISPS);
3015                 MDisps *fd = CustomData_get(fdata, findex, CD_MDISPS);
3016                 float (*disps)[3] = fd->disps;
3017                 int tot = mf->v4 ? 4 : 3;
3018                 int corners;
3019
3020                 if (CustomData_external_test(fdata, CD_MDISPS)) {
3021                         if (id && fdata->external) {
3022                                 CustomData_external_add(ldata, id, CD_MDISPS,
3023                                                         totloop, fdata->external->filename);
3024                         }
3025                 }
3026
3027                 corners = multires_mdisp_corners(fd);
3028
3029                 if (corners == 0) {
3030                         /* Empty MDisp layers appear in at least one of the sintel.blend files.
3031                          * Not sure w