293d2c34b07009b4db9e336255ed8269ca17a438
[blender.git] / source / blender / blenkernel / intern / curve.c
1 /*
2  * This program is free software; you can redistribute it and/or
3  * modify it under the terms of the GNU General Public License
4  * as published by the Free Software Foundation; either version 2
5  * of the License, or (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program; if not, write to the Free Software Foundation,
14  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
15  *
16  * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
17  * All rights reserved.
18  */
19
20 /** \file
21  * \ingroup bke
22  */
23
24 #include <math.h>  // floor
25 #include <string.h>
26 #include <stdlib.h>
27
28 #include "MEM_guardedalloc.h"
29
30 #include "BLI_blenlib.h"
31 #include "BLI_math.h"
32 #include "BLI_utildefines.h"
33 #include "BLI_ghash.h"
34 #include "BLI_linklist.h"
35
36 #include "DNA_anim_types.h"
37 #include "DNA_curve_types.h"
38 #include "DNA_material_types.h"
39
40 /* for dereferencing pointers */
41 #include "DNA_key_types.h"
42 #include "DNA_scene_types.h"
43 #include "DNA_vfont_types.h"
44 #include "DNA_object_types.h"
45
46 #include "BKE_animsys.h"
47 #include "BKE_curve.h"
48 #include "BKE_displist.h"
49 #include "BKE_font.h"
50 #include "BKE_key.h"
51 #include "BKE_library.h"
52 #include "BKE_main.h"
53 #include "BKE_object.h"
54 #include "BKE_material.h"
55
56 #include "DEG_depsgraph.h"
57
58 #include "CLG_log.h"
59
60 /* globals */
61
62 /* local */
63 static CLG_LogRef LOG = {"bke.curve"};
64
65 static int cu_isectLL(const float v1[3],
66                       const float v2[3],
67                       const float v3[3],
68                       const float v4[3],
69                       short cox,
70                       short coy,
71                       float *lambda,
72                       float *mu,
73                       float vec[3]);
74
75 /* frees editcurve entirely */
76 void BKE_curve_editfont_free(Curve *cu)
77 {
78   if (cu->editfont) {
79     EditFont *ef = cu->editfont;
80
81     if (ef->textbuf) {
82       MEM_freeN(ef->textbuf);
83     }
84     if (ef->textbufinfo) {
85       MEM_freeN(ef->textbufinfo);
86     }
87     if (ef->selboxes) {
88       MEM_freeN(ef->selboxes);
89     }
90
91     MEM_freeN(ef);
92     cu->editfont = NULL;
93   }
94 }
95
96 static void curve_editNurb_keyIndex_cv_free_cb(void *val)
97 {
98   CVKeyIndex *index = val;
99   MEM_freeN(index->orig_cv);
100   MEM_freeN(val);
101 }
102
103 void BKE_curve_editNurb_keyIndex_delCV(GHash *keyindex, const void *cv)
104 {
105   BLI_assert(keyindex != NULL);
106   BLI_ghash_remove(keyindex, cv, NULL, curve_editNurb_keyIndex_cv_free_cb);
107 }
108
109 void BKE_curve_editNurb_keyIndex_free(GHash **keyindex)
110 {
111   if (!(*keyindex)) {
112     return;
113   }
114   BLI_ghash_free(*keyindex, NULL, curve_editNurb_keyIndex_cv_free_cb);
115   *keyindex = NULL;
116 }
117
118 void BKE_curve_editNurb_free(Curve *cu)
119 {
120   if (cu->editnurb) {
121     BKE_nurbList_free(&cu->editnurb->nurbs);
122     BKE_curve_editNurb_keyIndex_free(&cu->editnurb->keyindex);
123     MEM_freeN(cu->editnurb);
124     cu->editnurb = NULL;
125   }
126 }
127
128 /** Free (or release) any data used by this curve (does not free the curve itself). */
129 void BKE_curve_free(Curve *cu)
130 {
131   BKE_animdata_free((ID *)cu, false);
132
133   BKE_curve_batch_cache_free(cu);
134
135   BKE_nurbList_free(&cu->nurb);
136   BKE_curve_editfont_free(cu);
137
138   BKE_curve_editNurb_free(cu);
139
140   MEM_SAFE_FREE(cu->mat);
141   MEM_SAFE_FREE(cu->str);
142   MEM_SAFE_FREE(cu->strinfo);
143   MEM_SAFE_FREE(cu->bb);
144   MEM_SAFE_FREE(cu->tb);
145 }
146
147 void BKE_curve_init(Curve *cu)
148 {
149   /* BLI_assert(MEMCMP_STRUCT_AFTER_IS_ZERO(cu, id)); */ /* cu->type is already initialized... */
150
151   copy_v3_fl(cu->size, 1.0f);
152   cu->flag = CU_DEFORM_BOUNDS_OFF | CU_PATH_RADIUS;
153   cu->pathlen = 100;
154   cu->resolu = cu->resolv = (cu->type == OB_SURF) ? 4 : 12;
155   cu->width = 1.0;
156   cu->wordspace = 1.0;
157   cu->spacing = cu->linedist = 1.0;
158   cu->fsize = 1.0;
159   cu->ulheight = 0.05;
160   cu->texflag = CU_AUTOSPACE;
161   cu->smallcaps_scale = 0.75f;
162   /* XXX: this one seems to be the best one in most cases, at least for curve deform... */
163   cu->twist_mode = CU_TWIST_MINIMUM;
164   cu->bevfac1 = 0.0f;
165   cu->bevfac2 = 1.0f;
166   cu->bevfac1_mapping = CU_BEVFAC_MAP_RESOLU;
167   cu->bevfac2_mapping = CU_BEVFAC_MAP_RESOLU;
168   cu->bevresol = 4;
169
170   cu->bb = BKE_boundbox_alloc_unit();
171
172   if (cu->type == OB_FONT) {
173     cu->flag |= CU_FRONT | CU_BACK;
174     cu->vfont = cu->vfontb = cu->vfonti = cu->vfontbi = BKE_vfont_builtin_get();
175     cu->vfont->id.us += 4;
176     cu->str = MEM_malloc_arrayN(12, sizeof(unsigned char), "str");
177     BLI_strncpy(cu->str, "Text", 12);
178     cu->len = cu->len_wchar = cu->pos = 4;
179     cu->strinfo = MEM_calloc_arrayN(12, sizeof(CharInfo), "strinfo new");
180     cu->totbox = cu->actbox = 1;
181     cu->tb = MEM_calloc_arrayN(MAXTEXTBOX, sizeof(TextBox), "textbox");
182     cu->tb[0].w = cu->tb[0].h = 0.0;
183   }
184 }
185
186 Curve *BKE_curve_add(Main *bmain, const char *name, int type)
187 {
188   Curve *cu;
189
190   cu = BKE_libblock_alloc(bmain, ID_CU, name, 0);
191   cu->type = type;
192
193   BKE_curve_init(cu);
194
195   return cu;
196 }
197
198 /**
199  * Only copy internal data of Curve ID from source to already allocated/initialized destination.
200  * You probably never want to use that directly, use BKE_id_copy or BKE_id_copy_ex for typical needs.
201  *
202  * WARNING! This function will not handle ID user count!
203  *
204  * \param flag: Copying options (see BKE_library.h's LIB_ID_COPY_... flags for more).
205  */
206 void BKE_curve_copy_data(Main *bmain, Curve *cu_dst, const Curve *cu_src, const int flag)
207 {
208   BLI_listbase_clear(&cu_dst->nurb);
209   BKE_nurbList_duplicate(&(cu_dst->nurb), &(cu_src->nurb));
210
211   cu_dst->mat = MEM_dupallocN(cu_src->mat);
212
213   cu_dst->str = MEM_dupallocN(cu_src->str);
214   cu_dst->strinfo = MEM_dupallocN(cu_src->strinfo);
215   cu_dst->tb = MEM_dupallocN(cu_src->tb);
216   cu_dst->bb = MEM_dupallocN(cu_src->bb);
217   cu_dst->batch_cache = NULL;
218
219   if (cu_src->key && (flag & LIB_ID_COPY_SHAPEKEY)) {
220     BKE_id_copy_ex(bmain, &cu_src->key->id, (ID **)&cu_dst->key, flag);
221   }
222
223   cu_dst->editnurb = NULL;
224   cu_dst->editfont = NULL;
225 }
226
227 Curve *BKE_curve_copy(Main *bmain, const Curve *cu)
228 {
229   Curve *cu_copy;
230   BKE_id_copy(bmain, &cu->id, (ID **)&cu_copy);
231   return cu_copy;
232 }
233
234 void BKE_curve_make_local(Main *bmain, Curve *cu, const bool lib_local)
235 {
236   BKE_id_make_local_generic(bmain, &cu->id, true, lib_local);
237 }
238
239 /* Get list of nurbs from editnurbs structure */
240 ListBase *BKE_curve_editNurbs_get(Curve *cu)
241 {
242   if (cu->editnurb) {
243     return &cu->editnurb->nurbs;
244   }
245
246   return NULL;
247 }
248
249 short BKE_curve_type_get(Curve *cu)
250 {
251   Nurb *nu;
252   int type = cu->type;
253
254   if (cu->vfont) {
255     return OB_FONT;
256   }
257
258   if (!cu->type) {
259     type = OB_CURVE;
260
261     for (nu = cu->nurb.first; nu; nu = nu->next) {
262       if (nu->pntsv > 1) {
263         type = OB_SURF;
264       }
265     }
266   }
267
268   return type;
269 }
270
271 void BKE_curve_curve_dimension_update(Curve *cu)
272 {
273   ListBase *nurbs = BKE_curve_nurbs_get(cu);
274   Nurb *nu = nurbs->first;
275
276   if (cu->flag & CU_3D) {
277     for (; nu; nu = nu->next) {
278       nu->flag &= ~CU_2D;
279     }
280   }
281   else {
282     for (; nu; nu = nu->next) {
283       nu->flag |= CU_2D;
284       BKE_nurb_test_2d(nu);
285
286       /* since the handles are moved they need to be auto-located again */
287       if (nu->type == CU_BEZIER) {
288         BKE_nurb_handles_calc(nu);
289       }
290     }
291   }
292 }
293
294 void BKE_curve_type_test(Object *ob)
295 {
296   ob->type = BKE_curve_type_get(ob->data);
297
298   if (ob->type == OB_CURVE) {
299     BKE_curve_curve_dimension_update((Curve *)ob->data);
300   }
301 }
302
303 void BKE_curve_boundbox_calc(Curve *cu, float r_loc[3], float r_size[3])
304 {
305   BoundBox *bb;
306   float min[3], max[3];
307   float mloc[3], msize[3];
308
309   if (cu->bb == NULL) {
310     cu->bb = MEM_callocN(sizeof(BoundBox), "boundbox");
311   }
312   bb = cu->bb;
313
314   if (!r_loc) {
315     r_loc = mloc;
316   }
317   if (!r_size) {
318     r_size = msize;
319   }
320
321   INIT_MINMAX(min, max);
322   if (!BKE_curve_minmax(cu, true, min, max)) {
323     min[0] = min[1] = min[2] = -1.0f;
324     max[0] = max[1] = max[2] = 1.0f;
325   }
326
327   mid_v3_v3v3(r_loc, min, max);
328
329   r_size[0] = (max[0] - min[0]) / 2.0f;
330   r_size[1] = (max[1] - min[1]) / 2.0f;
331   r_size[2] = (max[2] - min[2]) / 2.0f;
332
333   BKE_boundbox_init_from_minmax(bb, min, max);
334
335   bb->flag &= ~BOUNDBOX_DIRTY;
336 }
337
338 BoundBox *BKE_curve_boundbox_get(Object *ob)
339 {
340   /* This is Object-level data access, DO NOT touch to Mesh's bb, would be totally thread-unsafe. */
341   if (ob->runtime.bb == NULL || ob->runtime.bb->flag & BOUNDBOX_DIRTY) {
342     Curve *cu = ob->data;
343     float min[3], max[3];
344
345     INIT_MINMAX(min, max);
346     BKE_curve_minmax(cu, true, min, max);
347
348     if (ob->runtime.bb == NULL) {
349       ob->runtime.bb = MEM_mallocN(sizeof(*ob->runtime.bb), __func__);
350     }
351     BKE_boundbox_init_from_minmax(ob->runtime.bb, min, max);
352     ob->runtime.bb->flag &= ~BOUNDBOX_DIRTY;
353   }
354
355   return ob->runtime.bb;
356 }
357
358 void BKE_curve_texspace_calc(Curve *cu)
359 {
360   float loc[3], size[3];
361   int a;
362
363   BKE_curve_boundbox_calc(cu, loc, size);
364
365   if (cu->texflag & CU_AUTOSPACE) {
366     for (a = 0; a < 3; a++) {
367       if (size[a] == 0.0f) {
368         size[a] = 1.0f;
369       }
370       else if (size[a] > 0.0f && size[a] < 0.00001f) {
371         size[a] = 0.00001f;
372       }
373       else if (size[a] < 0.0f && size[a] > -0.00001f) {
374         size[a] = -0.00001f;
375       }
376     }
377
378     copy_v3_v3(cu->loc, loc);
379     copy_v3_v3(cu->size, size);
380     zero_v3(cu->rot);
381   }
382 }
383
384 BoundBox *BKE_curve_texspace_get(Curve *cu, float r_loc[3], float r_rot[3], float r_size[3])
385 {
386   if (cu->bb == NULL || (cu->bb->flag & BOUNDBOX_DIRTY)) {
387     BKE_curve_texspace_calc(cu);
388   }
389
390   if (r_loc) {
391     copy_v3_v3(r_loc, cu->loc);
392   }
393   if (r_rot) {
394     copy_v3_v3(r_rot, cu->rot);
395   }
396   if (r_size) {
397     copy_v3_v3(r_size, cu->size);
398   }
399
400   return cu->bb;
401 }
402
403 bool BKE_nurbList_index_get_co(ListBase *nurb, const int index, float r_co[3])
404 {
405   Nurb *nu;
406   int tot = 0;
407
408   for (nu = nurb->first; nu; nu = nu->next) {
409     int tot_nu;
410     if (nu->type == CU_BEZIER) {
411       tot_nu = nu->pntsu;
412       if (index - tot < tot_nu) {
413         copy_v3_v3(r_co, nu->bezt[index - tot].vec[1]);
414         return true;
415       }
416     }
417     else {
418       tot_nu = nu->pntsu * nu->pntsv;
419       if (index - tot < tot_nu) {
420         copy_v3_v3(r_co, nu->bp[index - tot].vec);
421         return true;
422       }
423     }
424     tot += tot_nu;
425   }
426
427   return false;
428 }
429
430 int BKE_nurbList_verts_count(ListBase *nurb)
431 {
432   Nurb *nu;
433   int tot = 0;
434
435   nu = nurb->first;
436   while (nu) {
437     if (nu->bezt) {
438       tot += 3 * nu->pntsu;
439     }
440     else if (nu->bp) {
441       tot += nu->pntsu * nu->pntsv;
442     }
443
444     nu = nu->next;
445   }
446   return tot;
447 }
448
449 int BKE_nurbList_verts_count_without_handles(ListBase *nurb)
450 {
451   Nurb *nu;
452   int tot = 0;
453
454   nu = nurb->first;
455   while (nu) {
456     if (nu->bezt) {
457       tot += nu->pntsu;
458     }
459     else if (nu->bp) {
460       tot += nu->pntsu * nu->pntsv;
461     }
462
463     nu = nu->next;
464   }
465   return tot;
466 }
467
468 /* **************** NURBS ROUTINES ******************** */
469
470 void BKE_nurb_free(Nurb *nu)
471 {
472
473   if (nu == NULL) {
474     return;
475   }
476
477   if (nu->bezt) {
478     MEM_freeN(nu->bezt);
479   }
480   nu->bezt = NULL;
481   if (nu->bp) {
482     MEM_freeN(nu->bp);
483   }
484   nu->bp = NULL;
485   if (nu->knotsu) {
486     MEM_freeN(nu->knotsu);
487   }
488   nu->knotsu = NULL;
489   if (nu->knotsv) {
490     MEM_freeN(nu->knotsv);
491   }
492   nu->knotsv = NULL;
493   /* if (nu->trim.first) freeNurblist(&(nu->trim)); */
494
495   MEM_freeN(nu);
496 }
497
498 void BKE_nurbList_free(ListBase *lb)
499 {
500   Nurb *nu, *next;
501
502   if (lb == NULL) {
503     return;
504   }
505
506   nu = lb->first;
507   while (nu) {
508     next = nu->next;
509     BKE_nurb_free(nu);
510     nu = next;
511   }
512   BLI_listbase_clear(lb);
513 }
514
515 Nurb *BKE_nurb_duplicate(const Nurb *nu)
516 {
517   Nurb *newnu;
518   int len;
519
520   newnu = (Nurb *)MEM_mallocN(sizeof(Nurb), "duplicateNurb");
521   if (newnu == NULL) {
522     return NULL;
523   }
524   memcpy(newnu, nu, sizeof(Nurb));
525
526   if (nu->bezt) {
527     newnu->bezt = (BezTriple *)MEM_malloc_arrayN(nu->pntsu, sizeof(BezTriple), "duplicateNurb2");
528     memcpy(newnu->bezt, nu->bezt, nu->pntsu * sizeof(BezTriple));
529   }
530   else {
531     len = nu->pntsu * nu->pntsv;
532     newnu->bp = (BPoint *)MEM_malloc_arrayN(len, sizeof(BPoint), "duplicateNurb3");
533     memcpy(newnu->bp, nu->bp, len * sizeof(BPoint));
534
535     newnu->knotsu = newnu->knotsv = NULL;
536
537     if (nu->knotsu) {
538       len = KNOTSU(nu);
539       if (len) {
540         newnu->knotsu = MEM_malloc_arrayN(len, sizeof(float), "duplicateNurb4");
541         memcpy(newnu->knotsu, nu->knotsu, sizeof(float) * len);
542       }
543     }
544     if (nu->pntsv > 1 && nu->knotsv) {
545       len = KNOTSV(nu);
546       if (len) {
547         newnu->knotsv = MEM_malloc_arrayN(len, sizeof(float), "duplicateNurb5");
548         memcpy(newnu->knotsv, nu->knotsv, sizeof(float) * len);
549       }
550     }
551   }
552   return newnu;
553 }
554
555 /* copy the nurb but allow for different number of points (to be copied after this) */
556 Nurb *BKE_nurb_copy(Nurb *src, int pntsu, int pntsv)
557 {
558   Nurb *newnu = (Nurb *)MEM_mallocN(sizeof(Nurb), "copyNurb");
559   memcpy(newnu, src, sizeof(Nurb));
560
561   if (pntsu == 1) {
562     SWAP(int, pntsu, pntsv);
563   }
564   newnu->pntsu = pntsu;
565   newnu->pntsv = pntsv;
566
567   /* caller can manually handle these arrays */
568   newnu->knotsu = NULL;
569   newnu->knotsv = NULL;
570
571   if (src->bezt) {
572     newnu->bezt = (BezTriple *)MEM_malloc_arrayN(pntsu * pntsv, sizeof(BezTriple), "copyNurb2");
573   }
574   else {
575     newnu->bp = (BPoint *)MEM_malloc_arrayN(pntsu * pntsv, sizeof(BPoint), "copyNurb3");
576   }
577
578   return newnu;
579 }
580
581 void BKE_nurbList_duplicate(ListBase *lb1, const ListBase *lb2)
582 {
583   Nurb *nu, *nun;
584
585   BKE_nurbList_free(lb1);
586
587   nu = lb2->first;
588   while (nu) {
589     nun = BKE_nurb_duplicate(nu);
590     BLI_addtail(lb1, nun);
591
592     nu = nu->next;
593   }
594 }
595
596 void BKE_nurb_test_2d(Nurb *nu)
597 {
598   BezTriple *bezt;
599   BPoint *bp;
600   int a;
601
602   if ((nu->flag & CU_2D) == 0) {
603     return;
604   }
605
606   if (nu->type == CU_BEZIER) {
607     a = nu->pntsu;
608     bezt = nu->bezt;
609     while (a--) {
610       bezt->vec[0][2] = 0.0;
611       bezt->vec[1][2] = 0.0;
612       bezt->vec[2][2] = 0.0;
613       bezt++;
614     }
615   }
616   else {
617     a = nu->pntsu * nu->pntsv;
618     bp = nu->bp;
619     while (a--) {
620       bp->vec[2] = 0.0;
621       bp++;
622     }
623   }
624 }
625
626 /**
627  * if use_radius is truth, minmax will take points' radius into account,
628  * which will make boundbox closer to beveled curve.
629  */
630 void BKE_nurb_minmax(Nurb *nu, bool use_radius, float min[3], float max[3])
631 {
632   BezTriple *bezt;
633   BPoint *bp;
634   int a;
635   float point[3];
636
637   if (nu->type == CU_BEZIER) {
638     a = nu->pntsu;
639     bezt = nu->bezt;
640     while (a--) {
641       if (use_radius) {
642         float radius_vector[3];
643         radius_vector[0] = radius_vector[1] = radius_vector[2] = bezt->radius;
644
645         add_v3_v3v3(point, bezt->vec[1], radius_vector);
646         minmax_v3v3_v3(min, max, point);
647
648         sub_v3_v3v3(point, bezt->vec[1], radius_vector);
649         minmax_v3v3_v3(min, max, point);
650       }
651       else {
652         minmax_v3v3_v3(min, max, bezt->vec[1]);
653       }
654       minmax_v3v3_v3(min, max, bezt->vec[0]);
655       minmax_v3v3_v3(min, max, bezt->vec[2]);
656       bezt++;
657     }
658   }
659   else {
660     a = nu->pntsu * nu->pntsv;
661     bp = nu->bp;
662     while (a--) {
663       if (nu->pntsv == 1 && use_radius) {
664         float radius_vector[3];
665         radius_vector[0] = radius_vector[1] = radius_vector[2] = bp->radius;
666
667         add_v3_v3v3(point, bp->vec, radius_vector);
668         minmax_v3v3_v3(min, max, point);
669
670         sub_v3_v3v3(point, bp->vec, radius_vector);
671         minmax_v3v3_v3(min, max, point);
672       }
673       else {
674         /* Surfaces doesn't use bevel, so no need to take radius into account. */
675         minmax_v3v3_v3(min, max, bp->vec);
676       }
677       bp++;
678     }
679   }
680 }
681
682 float BKE_nurb_calc_length(const Nurb *nu, int resolution)
683 {
684   BezTriple *bezt, *prevbezt;
685   BPoint *bp, *prevbp;
686   int a, b;
687   float length = 0.0f;
688   int resolu = resolution ? resolution : nu->resolu;
689   int pntsu = nu->pntsu;
690   float *points, *pntsit, *prevpntsit;
691
692   if (nu->type == CU_POLY) {
693     a = nu->pntsu - 1;
694     bp = nu->bp;
695     if (nu->flagu & CU_NURB_CYCLIC) {
696       ++a;
697       prevbp = nu->bp + (nu->pntsu - 1);
698     }
699     else {
700       prevbp = bp;
701       bp++;
702     }
703
704     while (a--) {
705       length += len_v3v3(prevbp->vec, bp->vec);
706       prevbp = bp;
707       ++bp;
708     }
709   }
710   else if (nu->type == CU_BEZIER) {
711     points = MEM_mallocN(sizeof(float[3]) * (resolu + 1), "getLength_bezier");
712     a = nu->pntsu - 1;
713     bezt = nu->bezt;
714     if (nu->flagu & CU_NURB_CYCLIC) {
715       ++a;
716       prevbezt = nu->bezt + (nu->pntsu - 1);
717     }
718     else {
719       prevbezt = bezt;
720       ++bezt;
721     }
722
723     while (a--) {
724       if (prevbezt->h2 == HD_VECT && bezt->h1 == HD_VECT) {
725         length += len_v3v3(prevbezt->vec[1], bezt->vec[1]);
726       }
727       else {
728         for (int j = 0; j < 3; j++) {
729           BKE_curve_forward_diff_bezier(prevbezt->vec[1][j],
730                                         prevbezt->vec[2][j],
731                                         bezt->vec[0][j],
732                                         bezt->vec[1][j],
733                                         points + j,
734                                         resolu,
735                                         3 * sizeof(float));
736         }
737
738         prevpntsit = pntsit = points;
739         b = resolu;
740         while (b--) {
741           pntsit += 3;
742           length += len_v3v3(prevpntsit, pntsit);
743           prevpntsit = pntsit;
744         }
745       }
746       prevbezt = bezt;
747       ++bezt;
748     }
749
750     MEM_freeN(points);
751   }
752   else if (nu->type == CU_NURBS) {
753     if (nu->pntsv == 1) {
754       /* important to zero for BKE_nurb_makeCurve. */
755       points = MEM_callocN(sizeof(float[3]) * pntsu * resolu, "getLength_nurbs");
756
757       BKE_nurb_makeCurve(nu, points, NULL, NULL, NULL, resolu, sizeof(float[3]));
758
759       if (nu->flagu & CU_NURB_CYCLIC) {
760         b = pntsu * resolu + 1;
761         prevpntsit = points + 3 * (pntsu * resolu - 1);
762         pntsit = points;
763       }
764       else {
765         b = (pntsu - 1) * resolu;
766         prevpntsit = points;
767         pntsit = points + 3;
768       }
769
770       while (--b) {
771         length += len_v3v3(prevpntsit, pntsit);
772         prevpntsit = pntsit;
773         pntsit += 3;
774       }
775
776       MEM_freeN(points);
777     }
778   }
779
780   return length;
781 }
782
783 /* be sure to call makeknots after this */
784 void BKE_nurb_points_add(Nurb *nu, int number)
785 {
786   BPoint *bp;
787   int i;
788
789   nu->bp = MEM_recallocN(nu->bp, (nu->pntsu + number) * sizeof(BPoint));
790
791   for (i = 0, bp = &nu->bp[nu->pntsu]; i < number; i++, bp++) {
792     bp->radius = 1.0f;
793   }
794
795   nu->pntsu += number;
796 }
797
798 void BKE_nurb_bezierPoints_add(Nurb *nu, int number)
799 {
800   BezTriple *bezt;
801   int i;
802
803   nu->bezt = MEM_recallocN(nu->bezt, (nu->pntsu + number) * sizeof(BezTriple));
804
805   for (i = 0, bezt = &nu->bezt[nu->pntsu]; i < number; i++, bezt++) {
806     bezt->radius = 1.0f;
807   }
808
809   nu->pntsu += number;
810 }
811
812 int BKE_nurb_index_from_uv(Nurb *nu, int u, int v)
813 {
814   const int totu = nu->pntsu;
815   const int totv = nu->pntsv;
816
817   if (nu->flagu & CU_NURB_CYCLIC) {
818     u = mod_i(u, totu);
819   }
820   else if (u < 0 || u >= totu) {
821     return -1;
822   }
823
824   if (nu->flagv & CU_NURB_CYCLIC) {
825     v = mod_i(v, totv);
826   }
827   else if (v < 0 || v >= totv) {
828     return -1;
829   }
830
831   return (v * totu) + u;
832 }
833
834 void BKE_nurb_index_to_uv(Nurb *nu, int index, int *r_u, int *r_v)
835 {
836   const int totu = nu->pntsu;
837   const int totv = nu->pntsv;
838   BLI_assert(index >= 0 && index < (nu->pntsu * nu->pntsv));
839   *r_u = (index % totu);
840   *r_v = (index / totu) % totv;
841 }
842
843 BezTriple *BKE_nurb_bezt_get_next(Nurb *nu, BezTriple *bezt)
844 {
845   BezTriple *bezt_next;
846
847   BLI_assert(ARRAY_HAS_ITEM(bezt, nu->bezt, nu->pntsu));
848
849   if (bezt == &nu->bezt[nu->pntsu - 1]) {
850     if (nu->flagu & CU_NURB_CYCLIC) {
851       bezt_next = nu->bezt;
852     }
853     else {
854       bezt_next = NULL;
855     }
856   }
857   else {
858     bezt_next = bezt + 1;
859   }
860
861   return bezt_next;
862 }
863
864 BPoint *BKE_nurb_bpoint_get_next(Nurb *nu, BPoint *bp)
865 {
866   BPoint *bp_next;
867
868   BLI_assert(ARRAY_HAS_ITEM(bp, nu->bp, nu->pntsu));
869
870   if (bp == &nu->bp[nu->pntsu - 1]) {
871     if (nu->flagu & CU_NURB_CYCLIC) {
872       bp_next = nu->bp;
873     }
874     else {
875       bp_next = NULL;
876     }
877   }
878   else {
879     bp_next = bp + 1;
880   }
881
882   return bp_next;
883 }
884
885 BezTriple *BKE_nurb_bezt_get_prev(Nurb *nu, BezTriple *bezt)
886 {
887   BezTriple *bezt_prev;
888
889   BLI_assert(ARRAY_HAS_ITEM(bezt, nu->bezt, nu->pntsu));
890   BLI_assert(nu->pntsv <= 1);
891
892   if (bezt == nu->bezt) {
893     if (nu->flagu & CU_NURB_CYCLIC) {
894       bezt_prev = &nu->bezt[nu->pntsu - 1];
895     }
896     else {
897       bezt_prev = NULL;
898     }
899   }
900   else {
901     bezt_prev = bezt - 1;
902   }
903
904   return bezt_prev;
905 }
906
907 BPoint *BKE_nurb_bpoint_get_prev(Nurb *nu, BPoint *bp)
908 {
909   BPoint *bp_prev;
910
911   BLI_assert(ARRAY_HAS_ITEM(bp, nu->bp, nu->pntsu));
912   BLI_assert(nu->pntsv == 1);
913
914   if (bp == nu->bp) {
915     if (nu->flagu & CU_NURB_CYCLIC) {
916       bp_prev = &nu->bp[nu->pntsu - 1];
917     }
918     else {
919       bp_prev = NULL;
920     }
921   }
922   else {
923     bp_prev = bp - 1;
924   }
925
926   return bp_prev;
927 }
928
929 void BKE_nurb_bezt_calc_normal(struct Nurb *UNUSED(nu), BezTriple *bezt, float r_normal[3])
930 {
931   /* calculate the axis matrix from the spline */
932   float dir_prev[3], dir_next[3];
933
934   sub_v3_v3v3(dir_prev, bezt->vec[0], bezt->vec[1]);
935   sub_v3_v3v3(dir_next, bezt->vec[1], bezt->vec[2]);
936
937   normalize_v3(dir_prev);
938   normalize_v3(dir_next);
939
940   add_v3_v3v3(r_normal, dir_prev, dir_next);
941   normalize_v3(r_normal);
942 }
943
944 void BKE_nurb_bezt_calc_plane(struct Nurb *nu, BezTriple *bezt, float r_plane[3])
945 {
946   float dir_prev[3], dir_next[3];
947
948   sub_v3_v3v3(dir_prev, bezt->vec[0], bezt->vec[1]);
949   sub_v3_v3v3(dir_next, bezt->vec[1], bezt->vec[2]);
950
951   normalize_v3(dir_prev);
952   normalize_v3(dir_next);
953
954   cross_v3_v3v3(r_plane, dir_prev, dir_next);
955   if (normalize_v3(r_plane) < FLT_EPSILON) {
956     BezTriple *bezt_prev = BKE_nurb_bezt_get_prev(nu, bezt);
957     BezTriple *bezt_next = BKE_nurb_bezt_get_next(nu, bezt);
958
959     if (bezt_prev) {
960       sub_v3_v3v3(dir_prev, bezt_prev->vec[1], bezt->vec[1]);
961       normalize_v3(dir_prev);
962     }
963     if (bezt_next) {
964       sub_v3_v3v3(dir_next, bezt->vec[1], bezt_next->vec[1]);
965       normalize_v3(dir_next);
966     }
967     cross_v3_v3v3(r_plane, dir_prev, dir_next);
968   }
969
970   /* matches with bones more closely */
971   {
972     float dir_mid[3], tvec[3];
973     add_v3_v3v3(dir_mid, dir_prev, dir_next);
974     cross_v3_v3v3(tvec, r_plane, dir_mid);
975     copy_v3_v3(r_plane, tvec);
976   }
977
978   normalize_v3(r_plane);
979 }
980
981 void BKE_nurb_bpoint_calc_normal(struct Nurb *nu, BPoint *bp, float r_normal[3])
982 {
983   BPoint *bp_prev = BKE_nurb_bpoint_get_prev(nu, bp);
984   BPoint *bp_next = BKE_nurb_bpoint_get_next(nu, bp);
985
986   zero_v3(r_normal);
987
988   if (bp_prev) {
989     float dir_prev[3];
990     sub_v3_v3v3(dir_prev, bp_prev->vec, bp->vec);
991     normalize_v3(dir_prev);
992     add_v3_v3(r_normal, dir_prev);
993   }
994   if (bp_next) {
995     float dir_next[3];
996     sub_v3_v3v3(dir_next, bp->vec, bp_next->vec);
997     normalize_v3(dir_next);
998     add_v3_v3(r_normal, dir_next);
999   }
1000
1001   normalize_v3(r_normal);
1002 }
1003
1004 void BKE_nurb_bpoint_calc_plane(struct Nurb *nu, BPoint *bp, float r_plane[3])
1005 {
1006   BPoint *bp_prev = BKE_nurb_bpoint_get_prev(nu, bp);
1007   BPoint *bp_next = BKE_nurb_bpoint_get_next(nu, bp);
1008
1009   float dir_prev[3] = {0.0f}, dir_next[3] = {0.0f};
1010
1011   if (bp_prev) {
1012     sub_v3_v3v3(dir_prev, bp_prev->vec, bp->vec);
1013     normalize_v3(dir_prev);
1014   }
1015   if (bp_next) {
1016     sub_v3_v3v3(dir_next, bp->vec, bp_next->vec);
1017     normalize_v3(dir_next);
1018   }
1019   cross_v3_v3v3(r_plane, dir_prev, dir_next);
1020
1021   /* matches with bones more closely */
1022   {
1023     float dir_mid[3], tvec[3];
1024     add_v3_v3v3(dir_mid, dir_prev, dir_next);
1025     cross_v3_v3v3(tvec, r_plane, dir_mid);
1026     copy_v3_v3(r_plane, tvec);
1027   }
1028
1029   normalize_v3(r_plane);
1030 }
1031
1032 /* ~~~~~~~~~~~~~~~~~~~~Non Uniform Rational B Spline calculations ~~~~~~~~~~~ */
1033
1034 static void calcknots(float *knots, const int pnts, const short order, const short flag)
1035 {
1036   /* knots: number of pnts NOT corrected for cyclic */
1037   const int pnts_order = pnts + order;
1038   float k;
1039   int a;
1040
1041   switch (flag & (CU_NURB_ENDPOINT | CU_NURB_BEZIER)) {
1042     case CU_NURB_ENDPOINT:
1043       k = 0.0;
1044       for (a = 1; a <= pnts_order; a++) {
1045         knots[a - 1] = k;
1046         if (a >= order && a <= pnts) {
1047           k += 1.0f;
1048         }
1049       }
1050       break;
1051     case CU_NURB_BEZIER:
1052       /* Warning, the order MUST be 2 or 4,
1053        * if this is not enforced, the displist will be corrupt */
1054       if (order == 4) {
1055         k = 0.34;
1056         for (a = 0; a < pnts_order; a++) {
1057           knots[a] = floorf(k);
1058           k += (1.0f / 3.0f);
1059         }
1060       }
1061       else if (order == 3) {
1062         k = 0.6f;
1063         for (a = 0; a < pnts_order; a++) {
1064           if (a >= order && a <= pnts) {
1065             k += 0.5f;
1066           }
1067           knots[a] = floorf(k);
1068         }
1069       }
1070       else {
1071         CLOG_ERROR(&LOG, "bez nurb curve order is not 3 or 4, should never happen");
1072       }
1073       break;
1074     default:
1075       for (a = 0; a < pnts_order; a++) {
1076         knots[a] = (float)a;
1077       }
1078       break;
1079   }
1080 }
1081
1082 static void makecyclicknots(float *knots, int pnts, short order)
1083 /* pnts, order: number of pnts NOT corrected for cyclic */
1084 {
1085   int a, b, order2, c;
1086
1087   if (knots == NULL) {
1088     return;
1089   }
1090
1091   order2 = order - 1;
1092
1093   /* do first long rows (order -1), remove identical knots at endpoints */
1094   if (order > 2) {
1095     b = pnts + order2;
1096     for (a = 1; a < order2; a++) {
1097       if (knots[b] != knots[b - a]) {
1098         break;
1099       }
1100     }
1101     if (a == order2) {
1102       knots[pnts + order - 2] += 1.0f;
1103     }
1104   }
1105
1106   b = order;
1107   c = pnts + order + order2;
1108   for (a = pnts + order2; a < c; a++) {
1109     knots[a] = knots[a - 1] + (knots[b] - knots[b - 1]);
1110     b--;
1111   }
1112 }
1113
1114 static void makeknots(Nurb *nu, short uv)
1115 {
1116   if (nu->type == CU_NURBS) {
1117     if (uv == 1) {
1118       if (nu->knotsu) {
1119         MEM_freeN(nu->knotsu);
1120       }
1121       if (BKE_nurb_check_valid_u(nu)) {
1122         nu->knotsu = MEM_calloc_arrayN(KNOTSU(nu) + 1, sizeof(float), "makeknots");
1123         if (nu->flagu & CU_NURB_CYCLIC) {
1124           calcknots(nu->knotsu, nu->pntsu, nu->orderu, 0); /* cyclic should be uniform */
1125           makecyclicknots(nu->knotsu, nu->pntsu, nu->orderu);
1126         }
1127         else {
1128           calcknots(nu->knotsu, nu->pntsu, nu->orderu, nu->flagu);
1129         }
1130       }
1131       else {
1132         nu->knotsu = NULL;
1133       }
1134     }
1135     else if (uv == 2) {
1136       if (nu->knotsv) {
1137         MEM_freeN(nu->knotsv);
1138       }
1139       if (BKE_nurb_check_valid_v(nu)) {
1140         nu->knotsv = MEM_calloc_arrayN(KNOTSV(nu) + 1, sizeof(float), "makeknots");
1141         if (nu->flagv & CU_NURB_CYCLIC) {
1142           calcknots(nu->knotsv, nu->pntsv, nu->orderv, 0); /* cyclic should be uniform */
1143           makecyclicknots(nu->knotsv, nu->pntsv, nu->orderv);
1144         }
1145         else {
1146           calcknots(nu->knotsv, nu->pntsv, nu->orderv, nu->flagv);
1147         }
1148       }
1149       else {
1150         nu->knotsv = NULL;
1151       }
1152     }
1153   }
1154 }
1155
1156 void BKE_nurb_knot_calc_u(Nurb *nu)
1157 {
1158   makeknots(nu, 1);
1159 }
1160
1161 void BKE_nurb_knot_calc_v(Nurb *nu)
1162 {
1163   makeknots(nu, 2);
1164 }
1165
1166 static void basisNurb(
1167     float t, short order, int pnts, float *knots, float *basis, int *start, int *end)
1168 {
1169   float d, e;
1170   int i, i1 = 0, i2 = 0, j, orderpluspnts, opp2, o2;
1171
1172   orderpluspnts = order + pnts;
1173   opp2 = orderpluspnts - 1;
1174
1175   /* this is for float inaccuracy */
1176   if (t < knots[0]) {
1177     t = knots[0];
1178   }
1179   else if (t > knots[opp2]) {
1180     t = knots[opp2];
1181   }
1182
1183   /* this part is order '1' */
1184   o2 = order + 1;
1185   for (i = 0; i < opp2; i++) {
1186     if (knots[i] != knots[i + 1] && t >= knots[i] && t <= knots[i + 1]) {
1187       basis[i] = 1.0;
1188       i1 = i - o2;
1189       if (i1 < 0) {
1190         i1 = 0;
1191       }
1192       i2 = i;
1193       i++;
1194       while (i < opp2) {
1195         basis[i] = 0.0;
1196         i++;
1197       }
1198       break;
1199     }
1200     else {
1201       basis[i] = 0.0;
1202     }
1203   }
1204   basis[i] = 0.0;
1205
1206   /* this is order 2, 3, ... */
1207   for (j = 2; j <= order; j++) {
1208
1209     if (i2 + j >= orderpluspnts) {
1210       i2 = opp2 - j;
1211     }
1212
1213     for (i = i1; i <= i2; i++) {
1214       if (basis[i] != 0.0f) {
1215         d = ((t - knots[i]) * basis[i]) / (knots[i + j - 1] - knots[i]);
1216       }
1217       else {
1218         d = 0.0f;
1219       }
1220
1221       if (basis[i + 1] != 0.0f) {
1222         e = ((knots[i + j] - t) * basis[i + 1]) / (knots[i + j] - knots[i + 1]);
1223       }
1224       else {
1225         e = 0.0;
1226       }
1227
1228       basis[i] = d + e;
1229     }
1230   }
1231
1232   *start = 1000;
1233   *end = 0;
1234
1235   for (i = i1; i <= i2; i++) {
1236     if (basis[i] > 0.0f) {
1237       *end = i;
1238       if (*start == 1000) {
1239         *start = i;
1240       }
1241     }
1242   }
1243 }
1244
1245 /**
1246  * \param coord_array: has to be (3 * 4 * resolu * resolv) in size, and zero-ed.
1247  */
1248 void BKE_nurb_makeFaces(const Nurb *nu, float *coord_array, int rowstride, int resolu, int resolv)
1249 {
1250   BPoint *bp;
1251   float *basisu, *basis, *basisv, *sum, *fp, *in;
1252   float u, v, ustart, uend, ustep, vstart, vend, vstep, sumdiv;
1253   int i, j, iofs, jofs, cycl, len, curu, curv;
1254   int istart, iend, jsta, jen, *jstart, *jend, ratcomp;
1255
1256   int totu = nu->pntsu * resolu, totv = nu->pntsv * resolv;
1257
1258   if (nu->knotsu == NULL || nu->knotsv == NULL) {
1259     return;
1260   }
1261   if (nu->orderu > nu->pntsu) {
1262     return;
1263   }
1264   if (nu->orderv > nu->pntsv) {
1265     return;
1266   }
1267   if (coord_array == NULL) {
1268     return;
1269   }
1270
1271   /* allocate and initialize */
1272   len = totu * totv;
1273   if (len == 0) {
1274     return;
1275   }
1276
1277   sum = (float *)MEM_calloc_arrayN(len, sizeof(float), "makeNurbfaces1");
1278
1279   bp = nu->bp;
1280   i = nu->pntsu * nu->pntsv;
1281   ratcomp = 0;
1282   while (i--) {
1283     if (bp->vec[3] != 1.0f) {
1284       ratcomp = 1;
1285       break;
1286     }
1287     bp++;
1288   }
1289
1290   fp = nu->knotsu;
1291   ustart = fp[nu->orderu - 1];
1292   if (nu->flagu & CU_NURB_CYCLIC) {
1293     uend = fp[nu->pntsu + nu->orderu - 1];
1294   }
1295   else {
1296     uend = fp[nu->pntsu];
1297   }
1298   ustep = (uend - ustart) / ((nu->flagu & CU_NURB_CYCLIC) ? totu : totu - 1);
1299
1300   basisu = (float *)MEM_malloc_arrayN(KNOTSU(nu), sizeof(float), "makeNurbfaces3");
1301
1302   fp = nu->knotsv;
1303   vstart = fp[nu->orderv - 1];
1304
1305   if (nu->flagv & CU_NURB_CYCLIC) {
1306     vend = fp[nu->pntsv + nu->orderv - 1];
1307   }
1308   else {
1309     vend = fp[nu->pntsv];
1310   }
1311   vstep = (vend - vstart) / ((nu->flagv & CU_NURB_CYCLIC) ? totv : totv - 1);
1312
1313   len = KNOTSV(nu);
1314   basisv = (float *)MEM_malloc_arrayN(len * totv, sizeof(float), "makeNurbfaces3");
1315   jstart = (int *)MEM_malloc_arrayN(totv, sizeof(float), "makeNurbfaces4");
1316   jend = (int *)MEM_malloc_arrayN(totv, sizeof(float), "makeNurbfaces5");
1317
1318   /* precalculation of basisv and jstart, jend */
1319   if (nu->flagv & CU_NURB_CYCLIC) {
1320     cycl = nu->orderv - 1;
1321   }
1322   else {
1323     cycl = 0;
1324   }
1325   v = vstart;
1326   basis = basisv;
1327   curv = totv;
1328   while (curv--) {
1329     basisNurb(v, nu->orderv, nu->pntsv + cycl, nu->knotsv, basis, jstart + curv, jend + curv);
1330     basis += KNOTSV(nu);
1331     v += vstep;
1332   }
1333
1334   if (nu->flagu & CU_NURB_CYCLIC) {
1335     cycl = nu->orderu - 1;
1336   }
1337   else {
1338     cycl = 0;
1339   }
1340   in = coord_array;
1341   u = ustart;
1342   curu = totu;
1343   while (curu--) {
1344     basisNurb(u, nu->orderu, nu->pntsu + cycl, nu->knotsu, basisu, &istart, &iend);
1345
1346     basis = basisv;
1347     curv = totv;
1348     while (curv--) {
1349       jsta = jstart[curv];
1350       jen = jend[curv];
1351
1352       /* calculate sum */
1353       sumdiv = 0.0;
1354       fp = sum;
1355
1356       for (j = jsta; j <= jen; j++) {
1357
1358         if (j >= nu->pntsv) {
1359           jofs = (j - nu->pntsv);
1360         }
1361         else {
1362           jofs = j;
1363         }
1364         bp = nu->bp + nu->pntsu * jofs + istart - 1;
1365
1366         for (i = istart; i <= iend; i++, fp++) {
1367           if (i >= nu->pntsu) {
1368             iofs = i - nu->pntsu;
1369             bp = nu->bp + nu->pntsu * jofs + iofs;
1370           }
1371           else {
1372             bp++;
1373           }
1374
1375           if (ratcomp) {
1376             *fp = basisu[i] * basis[j] * bp->vec[3];
1377             sumdiv += *fp;
1378           }
1379           else {
1380             *fp = basisu[i] * basis[j];
1381           }
1382         }
1383       }
1384
1385       if (ratcomp) {
1386         fp = sum;
1387         for (j = jsta; j <= jen; j++) {
1388           for (i = istart; i <= iend; i++, fp++) {
1389             *fp /= sumdiv;
1390           }
1391         }
1392       }
1393
1394       zero_v3(in);
1395
1396       /* one! (1.0) real point now */
1397       fp = sum;
1398       for (j = jsta; j <= jen; j++) {
1399
1400         if (j >= nu->pntsv) {
1401           jofs = (j - nu->pntsv);
1402         }
1403         else {
1404           jofs = j;
1405         }
1406         bp = nu->bp + nu->pntsu * jofs + istart - 1;
1407
1408         for (i = istart; i <= iend; i++, fp++) {
1409           if (i >= nu->pntsu) {
1410             iofs = i - nu->pntsu;
1411             bp = nu->bp + nu->pntsu * jofs + iofs;
1412           }
1413           else {
1414             bp++;
1415           }
1416
1417           if (*fp != 0.0f) {
1418             madd_v3_v3fl(in, bp->vec, *fp);
1419           }
1420         }
1421       }
1422
1423       in += 3;
1424       basis += KNOTSV(nu);
1425     }
1426     u += ustep;
1427     if (rowstride != 0) {
1428       in = (float *)(((unsigned char *)in) + (rowstride - 3 * totv * sizeof(*in)));
1429     }
1430   }
1431
1432   /* free */
1433   MEM_freeN(sum);
1434   MEM_freeN(basisu);
1435   MEM_freeN(basisv);
1436   MEM_freeN(jstart);
1437   MEM_freeN(jend);
1438 }
1439
1440 /**
1441  * \param coord_array: Has to be 3 * 4 * pntsu * resolu in size and zero-ed
1442  * \param tilt_array: set when non-NULL
1443  * \param radius_array: set when non-NULL
1444  */
1445 void BKE_nurb_makeCurve(const Nurb *nu,
1446                         float *coord_array,
1447                         float *tilt_array,
1448                         float *radius_array,
1449                         float *weight_array,
1450                         int resolu,
1451                         int stride)
1452 {
1453   const float eps = 1e-6f;
1454   BPoint *bp;
1455   float u, ustart, uend, ustep, sumdiv;
1456   float *basisu, *sum, *fp;
1457   float *coord_fp = coord_array, *tilt_fp = tilt_array, *radius_fp = radius_array,
1458         *weight_fp = weight_array;
1459   int i, len, istart, iend, cycl;
1460
1461   if (nu->knotsu == NULL) {
1462     return;
1463   }
1464   if (nu->orderu > nu->pntsu) {
1465     return;
1466   }
1467   if (coord_array == NULL) {
1468     return;
1469   }
1470
1471   /* allocate and initialize */
1472   len = nu->pntsu;
1473   if (len == 0) {
1474     return;
1475   }
1476   sum = (float *)MEM_calloc_arrayN(len, sizeof(float), "makeNurbcurve1");
1477
1478   resolu = (resolu * SEGMENTSU(nu));
1479
1480   if (resolu == 0) {
1481     MEM_freeN(sum);
1482     return;
1483   }
1484
1485   fp = nu->knotsu;
1486   ustart = fp[nu->orderu - 1];
1487   if (nu->flagu & CU_NURB_CYCLIC) {
1488     uend = fp[nu->pntsu + nu->orderu - 1];
1489   }
1490   else {
1491     uend = fp[nu->pntsu];
1492   }
1493   ustep = (uend - ustart) / (resolu - ((nu->flagu & CU_NURB_CYCLIC) ? 0 : 1));
1494
1495   basisu = (float *)MEM_malloc_arrayN(KNOTSU(nu), sizeof(float), "makeNurbcurve3");
1496
1497   if (nu->flagu & CU_NURB_CYCLIC) {
1498     cycl = nu->orderu - 1;
1499   }
1500   else {
1501     cycl = 0;
1502   }
1503
1504   u = ustart;
1505   while (resolu--) {
1506     basisNurb(u, nu->orderu, nu->pntsu + cycl, nu->knotsu, basisu, &istart, &iend);
1507
1508     /* calc sum */
1509     sumdiv = 0.0;
1510     fp = sum;
1511     bp = nu->bp + istart - 1;
1512     for (i = istart; i <= iend; i++, fp++) {
1513       if (i >= nu->pntsu) {
1514         bp = nu->bp + (i - nu->pntsu);
1515       }
1516       else {
1517         bp++;
1518       }
1519
1520       *fp = basisu[i] * bp->vec[3];
1521       sumdiv += *fp;
1522     }
1523     if ((sumdiv != 0.0f) && (sumdiv < 1.0f - eps || sumdiv > 1.0f + eps)) {
1524       /* is normalizing needed? */
1525       fp = sum;
1526       for (i = istart; i <= iend; i++, fp++) {
1527         *fp /= sumdiv;
1528       }
1529     }
1530
1531     zero_v3(coord_fp);
1532
1533     /* one! (1.0) real point */
1534     fp = sum;
1535     bp = nu->bp + istart - 1;
1536     for (i = istart; i <= iend; i++, fp++) {
1537       if (i >= nu->pntsu) {
1538         bp = nu->bp + (i - nu->pntsu);
1539       }
1540       else {
1541         bp++;
1542       }
1543
1544       if (*fp != 0.0f) {
1545         madd_v3_v3fl(coord_fp, bp->vec, *fp);
1546
1547         if (tilt_fp) {
1548           (*tilt_fp) += (*fp) * bp->tilt;
1549         }
1550
1551         if (radius_fp) {
1552           (*radius_fp) += (*fp) * bp->radius;
1553         }
1554
1555         if (weight_fp) {
1556           (*weight_fp) += (*fp) * bp->weight;
1557         }
1558       }
1559     }
1560
1561     coord_fp = POINTER_OFFSET(coord_fp, stride);
1562
1563     if (tilt_fp) {
1564       tilt_fp = POINTER_OFFSET(tilt_fp, stride);
1565     }
1566     if (radius_fp) {
1567       radius_fp = POINTER_OFFSET(radius_fp, stride);
1568     }
1569     if (weight_fp) {
1570       weight_fp = POINTER_OFFSET(weight_fp, stride);
1571     }
1572
1573     u += ustep;
1574   }
1575
1576   /* free */
1577   MEM_freeN(sum);
1578   MEM_freeN(basisu);
1579 }
1580
1581 /**
1582  * Calculate the length for arrays filled in by #BKE_curve_calc_coords_axis.
1583  */
1584 unsigned int BKE_curve_calc_coords_axis_len(const unsigned int bezt_array_len,
1585                                             const unsigned int resolu,
1586                                             const bool is_cyclic,
1587                                             const bool use_cyclic_duplicate_endpoint)
1588 {
1589   const unsigned int segments = bezt_array_len - (is_cyclic ? 0 : 1);
1590   const unsigned int points_len = (segments * resolu) +
1591                                   (is_cyclic ? (use_cyclic_duplicate_endpoint) : 1);
1592   return points_len;
1593 }
1594
1595 /**
1596  * Calculate an array for the entire curve (cyclic or non-cyclic).
1597  * \note Call for each axis.
1598  *
1599  * \param use_cyclic_duplicate_endpoint: Duplicate values at the beginning & end of the array.
1600  */
1601 void BKE_curve_calc_coords_axis(const BezTriple *bezt_array,
1602                                 const unsigned int bezt_array_len,
1603                                 const unsigned int resolu,
1604                                 const bool is_cyclic,
1605                                 const bool use_cyclic_duplicate_endpoint,
1606                                 /* array params */
1607                                 const unsigned int axis,
1608                                 const unsigned int stride,
1609                                 float *r_points)
1610 {
1611   const unsigned int points_len = BKE_curve_calc_coords_axis_len(
1612       bezt_array_len, resolu, is_cyclic, use_cyclic_duplicate_endpoint);
1613   float *r_points_offset = r_points;
1614
1615   const unsigned int resolu_stride = resolu * stride;
1616   const unsigned int bezt_array_last = bezt_array_len - 1;
1617
1618   for (unsigned int i = 0; i < bezt_array_last; i++) {
1619     const BezTriple *bezt_curr = &bezt_array[i];
1620     const BezTriple *bezt_next = &bezt_array[i + 1];
1621     BKE_curve_forward_diff_bezier(bezt_curr->vec[1][axis],
1622                                   bezt_curr->vec[2][axis],
1623                                   bezt_next->vec[0][axis],
1624                                   bezt_next->vec[1][axis],
1625                                   r_points_offset,
1626                                   (int)resolu,
1627                                   stride);
1628     r_points_offset = POINTER_OFFSET(r_points_offset, resolu_stride);
1629   }
1630
1631   if (is_cyclic) {
1632     const BezTriple *bezt_curr = &bezt_array[bezt_array_last];
1633     const BezTriple *bezt_next = &bezt_array[0];
1634     BKE_curve_forward_diff_bezier(bezt_curr->vec[1][axis],
1635                                   bezt_curr->vec[2][axis],
1636                                   bezt_next->vec[0][axis],
1637                                   bezt_next->vec[1][axis],
1638                                   r_points_offset,
1639                                   (int)resolu,
1640                                   stride);
1641     r_points_offset = POINTER_OFFSET(r_points_offset, resolu_stride);
1642     if (use_cyclic_duplicate_endpoint) {
1643       *r_points_offset = *r_points;
1644       r_points_offset = POINTER_OFFSET(r_points_offset, stride);
1645     }
1646   }
1647   else {
1648     float *r_points_last = POINTER_OFFSET(r_points, bezt_array_last * resolu_stride);
1649     *r_points_last = bezt_array[bezt_array_last].vec[1][axis];
1650     r_points_offset = POINTER_OFFSET(r_points_offset, stride);
1651   }
1652
1653   BLI_assert(POINTER_OFFSET(r_points, points_len * stride) == r_points_offset);
1654   UNUSED_VARS_NDEBUG(points_len);
1655 }
1656
1657 /* forward differencing method for bezier curve */
1658 void BKE_curve_forward_diff_bezier(
1659     float q0, float q1, float q2, float q3, float *p, int it, int stride)
1660 {
1661   float rt0, rt1, rt2, rt3, f;
1662   int a;
1663
1664   f = (float)it;
1665   rt0 = q0;
1666   rt1 = 3.0f * (q1 - q0) / f;
1667   f *= f;
1668   rt2 = 3.0f * (q0 - 2.0f * q1 + q2) / f;
1669   f *= it;
1670   rt3 = (q3 - q0 + 3.0f * (q1 - q2)) / f;
1671
1672   q0 = rt0;
1673   q1 = rt1 + rt2 + rt3;
1674   q2 = 2 * rt2 + 6 * rt3;
1675   q3 = 6 * rt3;
1676
1677   for (a = 0; a <= it; a++) {
1678     *p = q0;
1679     p = POINTER_OFFSET(p, stride);
1680     q0 += q1;
1681     q1 += q2;
1682     q2 += q3;
1683   }
1684 }
1685
1686 /* forward differencing method for first derivative of cubic bezier curve */
1687 void BKE_curve_forward_diff_tangent_bezier(
1688     float q0, float q1, float q2, float q3, float *p, int it, int stride)
1689 {
1690   float rt0, rt1, rt2, f;
1691   int a;
1692
1693   f = 1.0f / (float)it;
1694
1695   rt0 = 3.0f * (q1 - q0);
1696   rt1 = f * (3.0f * (q3 - q0) + 9.0f * (q1 - q2));
1697   rt2 = 6.0f * (q0 + q2) - 12.0f * q1;
1698
1699   q0 = rt0;
1700   q1 = f * (rt1 + rt2);
1701   q2 = 2.0f * f * rt1;
1702
1703   for (a = 0; a <= it; a++) {
1704     *p = q0;
1705     p = POINTER_OFFSET(p, stride);
1706     q0 += q1;
1707     q1 += q2;
1708   }
1709 }
1710
1711 static void forward_diff_bezier_cotangent(const float p0[3],
1712                                           const float p1[3],
1713                                           const float p2[3],
1714                                           const float p3[3],
1715                                           float p[3],
1716                                           int it,
1717                                           int stride)
1718 {
1719   /* note that these are not perpendicular to the curve
1720    * they need to be rotated for this,
1721    *
1722    * This could also be optimized like BKE_curve_forward_diff_bezier */
1723   int a;
1724   for (a = 0; a <= it; a++) {
1725     float t = (float)a / (float)it;
1726
1727     int i;
1728     for (i = 0; i < 3; i++) {
1729       p[i] = (-6.0f * t + 6.0f) * p0[i] + (18.0f * t - 12.0f) * p1[i] +
1730              (-18.0f * t + 6.0f) * p2[i] + (6.0f * t) * p3[i];
1731     }
1732     normalize_v3(p);
1733     p = POINTER_OFFSET(p, stride);
1734   }
1735 }
1736
1737 /* ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ */
1738
1739 float *BKE_curve_surf_make_orco(Object *ob)
1740 {
1741   /* Note: this function is used in convertblender only atm, so
1742    * suppose nonzero curve's render resolution should always be used */
1743   Curve *cu = ob->data;
1744   Nurb *nu;
1745   int a, b, tot = 0;
1746   int sizeu, sizev;
1747   int resolu, resolv;
1748   float *fp, *coord_array;
1749
1750   /* first calculate the size of the datablock */
1751   nu = cu->nurb.first;
1752   while (nu) {
1753     /* as we want to avoid the seam in a cyclic nurbs
1754      * texture wrapping, reserve extra orco data space to save these extra needed
1755      * vertex based UV coordinates for the meridian vertices.
1756      * Vertices on the 0/2pi boundary are not duplicated inside the displist but later in
1757      * the renderface/vert construction.
1758      *
1759      * See also convertblender.c: init_render_surf()
1760      */
1761
1762     resolu = cu->resolu_ren ? cu->resolu_ren : nu->resolu;
1763     resolv = cu->resolv_ren ? cu->resolv_ren : nu->resolv;
1764
1765     sizeu = nu->pntsu * resolu;
1766     sizev = nu->pntsv * resolv;
1767     if (nu->flagu & CU_NURB_CYCLIC) {
1768       sizeu++;
1769     }
1770     if (nu->flagv & CU_NURB_CYCLIC) {
1771       sizev++;
1772     }
1773     if (nu->pntsv > 1) {
1774       tot += sizeu * sizev;
1775     }
1776
1777     nu = nu->next;
1778   }
1779   /* makeNurbfaces wants zeros */
1780   fp = coord_array = MEM_calloc_arrayN(tot, 3 * sizeof(float), "make_orco");
1781
1782   nu = cu->nurb.first;
1783   while (nu) {
1784     resolu = cu->resolu_ren ? cu->resolu_ren : nu->resolu;
1785     resolv = cu->resolv_ren ? cu->resolv_ren : nu->resolv;
1786
1787     if (nu->pntsv > 1) {
1788       sizeu = nu->pntsu * resolu;
1789       sizev = nu->pntsv * resolv;
1790
1791       if (nu->flagu & CU_NURB_CYCLIC) {
1792         sizeu++;
1793       }
1794       if (nu->flagv & CU_NURB_CYCLIC) {
1795         sizev++;
1796       }
1797
1798       if (cu->flag & CU_UV_ORCO) {
1799         for (b = 0; b < sizeu; b++) {
1800           for (a = 0; a < sizev; a++) {
1801
1802             if (sizev < 2) {
1803               fp[0] = 0.0f;
1804             }
1805             else {
1806               fp[0] = -1.0f + 2.0f * ((float)a) / (sizev - 1);
1807             }
1808
1809             if (sizeu < 2) {
1810               fp[1] = 0.0f;
1811             }
1812             else {
1813               fp[1] = -1.0f + 2.0f * ((float)b) / (sizeu - 1);
1814             }
1815
1816             fp[2] = 0.0;
1817
1818             fp += 3;
1819           }
1820         }
1821       }
1822       else {
1823         int size = (nu->pntsu * resolu) * (nu->pntsv * resolv) * 3 * sizeof(float);
1824         float *_tdata = MEM_mallocN(size, "temp data");
1825         float *tdata = _tdata;
1826
1827         BKE_nurb_makeFaces(nu, tdata, 0, resolu, resolv);
1828
1829         for (b = 0; b < sizeu; b++) {
1830           int use_b = b;
1831           if (b == sizeu - 1 && (nu->flagu & CU_NURB_CYCLIC)) {
1832             use_b = false;
1833           }
1834
1835           for (a = 0; a < sizev; a++) {
1836             int use_a = a;
1837             if (a == sizev - 1 && (nu->flagv & CU_NURB_CYCLIC)) {
1838               use_a = false;
1839             }
1840
1841             tdata = _tdata + 3 * (use_b * (nu->pntsv * resolv) + use_a);
1842
1843             fp[0] = (tdata[0] - cu->loc[0]) / cu->size[0];
1844             fp[1] = (tdata[1] - cu->loc[1]) / cu->size[1];
1845             fp[2] = (tdata[2] - cu->loc[2]) / cu->size[2];
1846             fp += 3;
1847           }
1848         }
1849
1850         MEM_freeN(_tdata);
1851       }
1852     }
1853     nu = nu->next;
1854   }
1855
1856   return coord_array;
1857 }
1858
1859 /* NOTE: This routine is tied to the order of vertex
1860  * built by displist and as passed to the renderer.
1861  */
1862 float *BKE_curve_make_orco(Depsgraph *depsgraph, Scene *scene, Object *ob, int *r_numVerts)
1863 {
1864   Curve *cu = ob->data;
1865   DispList *dl;
1866   int u, v, numVerts;
1867   float *fp, *coord_array;
1868   ListBase disp = {NULL, NULL};
1869
1870   BKE_displist_make_curveTypes_forOrco(depsgraph, scene, ob, &disp, NULL);
1871
1872   numVerts = 0;
1873   for (dl = disp.first; dl; dl = dl->next) {
1874     if (dl->type == DL_INDEX3) {
1875       numVerts += dl->nr;
1876     }
1877     else if (dl->type == DL_SURF) {
1878       /* convertblender.c uses the Surface code for creating renderfaces when cyclic U only
1879        * (closed circle beveling)
1880        */
1881       if (dl->flag & DL_CYCL_U) {
1882         if (dl->flag & DL_CYCL_V) {
1883           numVerts += (dl->parts + 1) * (dl->nr + 1);
1884         }
1885         else {
1886           numVerts += dl->parts * (dl->nr + 1);
1887         }
1888       }
1889       else if (dl->flag & DL_CYCL_V) {
1890         numVerts += (dl->parts + 1) * dl->nr;
1891       }
1892       else {
1893         numVerts += dl->parts * dl->nr;
1894       }
1895     }
1896   }
1897
1898   if (r_numVerts) {
1899     *r_numVerts = numVerts;
1900   }
1901
1902   fp = coord_array = MEM_malloc_arrayN(numVerts, 3 * sizeof(float), "cu_orco");
1903   for (dl = disp.first; dl; dl = dl->next) {
1904     if (dl->type == DL_INDEX3) {
1905       for (u = 0; u < dl->nr; u++, fp += 3) {
1906         if (cu->flag & CU_UV_ORCO) {
1907           fp[0] = 2.0f * u / (dl->nr - 1) - 1.0f;
1908           fp[1] = 0.0;
1909           fp[2] = 0.0;
1910         }
1911         else {
1912           copy_v3_v3(fp, &dl->verts[u * 3]);
1913
1914           fp[0] = (fp[0] - cu->loc[0]) / cu->size[0];
1915           fp[1] = (fp[1] - cu->loc[1]) / cu->size[1];
1916           fp[2] = (fp[2] - cu->loc[2]) / cu->size[2];
1917         }
1918       }
1919     }
1920     else if (dl->type == DL_SURF) {
1921       int sizeu = dl->nr, sizev = dl->parts;
1922
1923       /* exception as handled in convertblender.c too */
1924       if (dl->flag & DL_CYCL_U) {
1925         sizeu++;
1926         if (dl->flag & DL_CYCL_V) {
1927           sizev++;
1928         }
1929       }
1930       else if (dl->flag & DL_CYCL_V) {
1931         sizev++;
1932       }
1933
1934       for (u = 0; u < sizev; u++) {
1935         for (v = 0; v < sizeu; v++, fp += 3) {
1936           if (cu->flag & CU_UV_ORCO) {
1937             fp[0] = 2.0f * u / (sizev - 1) - 1.0f;
1938             fp[1] = 2.0f * v / (sizeu - 1) - 1.0f;
1939             fp[2] = 0.0;
1940           }
1941           else {
1942             const float *vert;
1943             int realv = v % dl->nr;
1944             int realu = u % dl->parts;
1945
1946             vert = dl->verts + 3 * (dl->nr * realu + realv);
1947             copy_v3_v3(fp, vert);
1948
1949             fp[0] = (fp[0] - cu->loc[0]) / cu->size[0];
1950             fp[1] = (fp[1] - cu->loc[1]) / cu->size[1];
1951             fp[2] = (fp[2] - cu->loc[2]) / cu->size[2];
1952           }
1953         }
1954       }
1955     }
1956   }
1957
1958   BKE_displist_free(&disp);
1959
1960   return coord_array;
1961 }
1962
1963 /* ***************** BEVEL ****************** */
1964
1965 void BKE_curve_bevel_make(Depsgraph *depsgraph,
1966                           Scene *scene,
1967                           Object *ob,
1968                           ListBase *disp,
1969                           const bool for_render,
1970                           const bool use_render_resolution,
1971                           LinkNode *ob_cyclic_list)
1972 {
1973   DispList *dl, *dlnew;
1974   Curve *bevcu, *cu;
1975   float *fp, facx, facy, angle, dangle;
1976   int nr, a;
1977
1978   cu = ob->data;
1979   BLI_listbase_clear(disp);
1980
1981   /* if a font object is being edited, then do nothing */
1982   // XXX  if ( ob == obedit && ob->type == OB_FONT ) return;
1983
1984   if (cu->bevobj) {
1985     if (cu->bevobj->type != OB_CURVE) {
1986       return;
1987     }
1988
1989     bevcu = cu->bevobj->data;
1990     if (bevcu->ext1 == 0.0f && bevcu->ext2 == 0.0f) {
1991       ListBase bevdisp = {NULL, NULL};
1992       facx = cu->bevobj->scale[0];
1993       facy = cu->bevobj->scale[1];
1994
1995       if (for_render) {
1996         if (BLI_linklist_index(ob_cyclic_list, cu->bevobj) == -1) {
1997           BKE_displist_make_curveTypes_forRender(depsgraph,
1998                                                  scene,
1999                                                  cu->bevobj,
2000                                                  &bevdisp,
2001                                                  NULL,
2002                                                  false,
2003                                                  use_render_resolution,
2004                                                  &(LinkNode){
2005                                                      .link = ob,
2006                                                      .next = ob_cyclic_list,
2007                                                  });
2008           dl = bevdisp.first;
2009         }
2010         else {
2011           dl = NULL;
2012         }
2013       }
2014       else if (cu->bevobj->runtime.curve_cache) {
2015         dl = cu->bevobj->runtime.curve_cache->disp.first;
2016       }
2017       else {
2018         BLI_assert(cu->bevobj->runtime.curve_cache != NULL);
2019         dl = NULL;
2020       }
2021
2022       while (dl) {
2023         if (ELEM(dl->type, DL_POLY, DL_SEGM)) {
2024           dlnew = MEM_mallocN(sizeof(DispList), "makebevelcurve1");
2025           *dlnew = *dl;
2026           dlnew->verts = MEM_malloc_arrayN(
2027               dl->parts * dl->nr, 3 * sizeof(float), "makebevelcurve1");
2028           memcpy(dlnew->verts, dl->verts, 3 * sizeof(float) * dl->parts * dl->nr);
2029
2030           if (dlnew->type == DL_SEGM) {
2031             dlnew->flag |= (DL_FRONT_CURVE | DL_BACK_CURVE);
2032           }
2033
2034           BLI_addtail(disp, dlnew);
2035           fp = dlnew->verts;
2036           nr = dlnew->parts * dlnew->nr;
2037           while (nr--) {
2038             fp[2] = fp[1] * facy;
2039             fp[1] = -fp[0] * facx;
2040             fp[0] = 0.0;
2041             fp += 3;
2042           }
2043         }
2044         dl = dl->next;
2045       }
2046
2047       BKE_displist_free(&bevdisp);
2048     }
2049   }
2050   else if (cu->ext1 == 0.0f && cu->ext2 == 0.0f) {
2051     /* pass */
2052   }
2053   else if (cu->ext2 == 0.0f) {
2054     dl = MEM_callocN(sizeof(DispList), "makebevelcurve2");
2055     dl->verts = MEM_malloc_arrayN(2, sizeof(float[3]), "makebevelcurve2");
2056     BLI_addtail(disp, dl);
2057     dl->type = DL_SEGM;
2058     dl->parts = 1;
2059     dl->flag = DL_FRONT_CURVE | DL_BACK_CURVE;
2060     dl->nr = 2;
2061
2062     fp = dl->verts;
2063     fp[0] = fp[1] = 0.0;
2064     fp[2] = -cu->ext1;
2065     fp[3] = fp[4] = 0.0;
2066     fp[5] = cu->ext1;
2067   }
2068   else if ((cu->flag & (CU_FRONT | CU_BACK)) == 0 &&
2069            cu->ext1 == 0.0f) { /* we make a full round bevel in that case */
2070     nr = 4 + 2 * cu->bevresol;
2071
2072     dl = MEM_callocN(sizeof(DispList), "makebevelcurve p1");
2073     dl->verts = MEM_malloc_arrayN(nr, sizeof(float[3]), "makebevelcurve p1");
2074     BLI_addtail(disp, dl);
2075     dl->type = DL_POLY;
2076     dl->parts = 1;
2077     dl->flag = DL_BACK_CURVE;
2078     dl->nr = nr;
2079
2080     /* a circle */
2081     fp = dl->verts;
2082     dangle = (2.0f * (float)M_PI / (nr));
2083     angle = -(nr - 1) * dangle;
2084
2085     for (a = 0; a < nr; a++) {
2086       fp[0] = 0.0;
2087       fp[1] = (cosf(angle) * (cu->ext2));
2088       fp[2] = (sinf(angle) * (cu->ext2)) - cu->ext1;
2089       angle += dangle;
2090       fp += 3;
2091     }
2092   }
2093   else {
2094     short dnr;
2095
2096     /* bevel now in three parts, for proper vertex normals */
2097     /* part 1, back */
2098
2099     if ((cu->flag & CU_BACK) || !(cu->flag & CU_FRONT)) {
2100       dnr = nr = 2 + cu->bevresol;
2101       if ((cu->flag & (CU_FRONT | CU_BACK)) == 0) {
2102         nr = 3 + 2 * cu->bevresol;
2103       }
2104       dl = MEM_callocN(sizeof(DispList), "makebevelcurve p1");
2105       dl->verts = MEM_malloc_arrayN(nr, sizeof(float[3]), "makebevelcurve p1");
2106       BLI_addtail(disp, dl);
2107       dl->type = DL_SEGM;
2108       dl->parts = 1;
2109       dl->flag = DL_BACK_CURVE;
2110       dl->nr = nr;
2111
2112       /* half a circle */
2113       fp = dl->verts;
2114       dangle = ((float)M_PI_2 / (dnr - 1));
2115       angle = -(nr - 1) * dangle;
2116
2117       for (a = 0; a < nr; a++) {
2118         fp[0] = 0.0;
2119         fp[1] = (float)(cosf(angle) * (cu->ext2));
2120         fp[2] = (float)(sinf(angle) * (cu->ext2)) - cu->ext1;
2121         angle += dangle;
2122         fp += 3;
2123       }
2124     }
2125
2126     /* part 2, sidefaces */
2127     if (cu->ext1 != 0.0f) {
2128       nr = 2;
2129
2130       dl = MEM_callocN(sizeof(DispList), "makebevelcurve p2");
2131       dl->verts = MEM_malloc_arrayN(nr, sizeof(float[3]), "makebevelcurve p2");
2132       BLI_addtail(disp, dl);
2133       dl->type = DL_SEGM;
2134       dl->parts = 1;
2135       dl->nr = nr;
2136
2137       fp = dl->verts;
2138       fp[1] = cu->ext2;
2139       fp[2] = -cu->ext1;
2140       fp[4] = cu->ext2;
2141       fp[5] = cu->ext1;
2142
2143       if ((cu->flag & (CU_FRONT | CU_BACK)) == 0) {
2144         dl = MEM_dupallocN(dl);
2145         dl->verts = MEM_dupallocN(dl->verts);
2146         BLI_addtail(disp, dl);
2147
2148         fp = dl->verts;
2149         fp[1] = -fp[1];
2150         fp[2] = -fp[2];
2151         fp[4] = -fp[4];
2152         fp[5] = -fp[5];
2153       }
2154     }
2155
2156     /* part 3, front */
2157     if ((cu->flag & CU_FRONT) || !(cu->flag & CU_BACK)) {
2158       dnr = nr = 2 + cu->bevresol;
2159       if ((cu->flag & (CU_FRONT | CU_BACK)) == 0) {
2160         nr = 3 + 2 * cu->bevresol;
2161       }
2162       dl = MEM_callocN(sizeof(DispList), "makebevelcurve p3");
2163       dl->verts = MEM_malloc_arrayN(nr, sizeof(float[3]), "makebevelcurve p3");
2164       BLI_addtail(disp, dl);
2165       dl->type = DL_SEGM;
2166       dl->flag = DL_FRONT_CURVE;
2167       dl->parts = 1;
2168       dl->nr = nr;
2169
2170       /* half a circle */
2171       fp = dl->verts;
2172       angle = 0.0;
2173       dangle = ((float)M_PI_2 / (dnr - 1));
2174
2175       for (a = 0; a < nr; a++) {
2176         fp[0] = 0.0;
2177         fp[1] = (float)(cosf(angle) * (cu->ext2));
2178         fp[2] = (float)(sinf(angle) * (cu->ext2)) + cu->ext1;
2179         angle += dangle;
2180         fp += 3;
2181       }
2182     }
2183   }
2184 }
2185
2186 static int cu_isectLL(const float v1[3],
2187                       const float v2[3],
2188                       const float v3[3],
2189                       const float v4[3],
2190                       short cox,
2191                       short coy,
2192                       float *lambda,
2193                       float *mu,
2194                       float vec[3])
2195 {
2196   /* return:
2197    * -1: collinear
2198    *  0: no intersection of segments
2199    *  1: exact intersection of segments
2200    *  2: cross-intersection of segments
2201    */
2202   float deler;
2203
2204   deler = (v1[cox] - v2[cox]) * (v3[coy] - v4[coy]) - (v3[cox] - v4[cox]) * (v1[coy] - v2[coy]);
2205   if (deler == 0.0f) {
2206     return -1;
2207   }
2208
2209   *lambda = (v1[coy] - v3[coy]) * (v3[cox] - v4[cox]) - (v1[cox] - v3[cox]) * (v3[coy] - v4[coy]);
2210   *lambda = -(*lambda / deler);
2211
2212   deler = v3[coy] - v4[coy];
2213   if (deler == 0) {
2214     deler = v3[cox] - v4[cox];
2215     *mu = -(*lambda * (v2[cox] - v1[cox]) + v1[cox] - v3[cox]) / deler;
2216   }
2217   else {
2218     *mu = -(*lambda * (v2[coy] - v1[coy]) + v1[coy] - v3[coy]) / deler;
2219   }
2220   vec[cox] = *lambda * (v2[cox] - v1[cox]) + v1[cox];
2221   vec[coy] = *lambda * (v2[coy] - v1[coy]) + v1[coy];
2222
2223   if (*lambda >= 0.0f && *lambda <= 1.0f && *mu >= 0.0f && *mu <= 1.0f) {
2224     if (*lambda == 0.0f || *lambda == 1.0f || *mu == 0.0f || *mu == 1.0f) {
2225       return 1;
2226     }
2227     return 2;
2228   }
2229   return 0;
2230 }
2231
2232 static bool bevelinside(BevList *bl1, BevList *bl2)
2233 {
2234   /* is bl2 INSIDE bl1 ? with left-right method and "lambda's" */
2235   /* returns '1' if correct hole  */
2236   BevPoint *bevp, *prevbevp;
2237   float min, max, vec[3], hvec1[3], hvec2[3], lab, mu;
2238   int nr, links = 0, rechts = 0, mode;
2239
2240   /* take first vertex of possible hole */
2241
2242   bevp = bl2->bevpoints;
2243   hvec1[0] = bevp->vec[0];
2244   hvec1[1] = bevp->vec[1];
2245   hvec1[2] = 0.0;
2246   copy_v3_v3(hvec2, hvec1);
2247   hvec2[0] += 1000;
2248
2249   /* test it with all edges of potential surrounding poly */
2250   /* count number of transitions left-right  */
2251
2252   bevp = bl1->bevpoints;
2253   nr = bl1->nr;
2254   prevbevp = bevp + (nr - 1);
2255
2256   while (nr--) {
2257     min = prevbevp->vec[1];
2258     max = bevp->vec[1];
2259     if (max < min) {
2260       min = max;
2261       max = prevbevp->vec[1];
2262     }
2263     if (min != max) {
2264       if (min <= hvec1[1] && max >= hvec1[1]) {
2265         /* there's a transition, calc intersection point */
2266         mode = cu_isectLL(prevbevp->vec, bevp->vec, hvec1, hvec2, 0, 1, &lab, &mu, vec);
2267         /* if lab==0.0 or lab==1.0 then the edge intersects exactly a transition
2268          * only allow for one situation: we choose lab= 1.0
2269          */
2270         if (mode >= 0 && lab != 0.0f) {
2271           if (vec[0] < hvec1[0]) {
2272             links++;
2273           }
2274           else {
2275             rechts++;
2276           }
2277         }
2278       }
2279     }
2280     prevbevp = bevp;
2281     bevp++;
2282   }
2283
2284   return (links & 1) && (rechts & 1);
2285 }
2286
2287 struct BevelSort {
2288   BevList *bl;
2289   float left;
2290   int dir;
2291 };
2292
2293 static int vergxcobev(const void *a1, const void *a2)
2294 {
2295   const struct BevelSort *x1 = a1, *x2 = a2;
2296
2297   if (x1->left > x2->left) {
2298     return 1;
2299   }
2300   else if (x1->left < x2->left) {
2301     return -1;
2302   }
2303   return 0;
2304 }
2305
2306 /* this function cannot be replaced with atan2, but why? */
2307
2308 static void calc_bevel_sin_cos(
2309     float x1, float y1, float x2, float y2, float *r_sina, float *r_cosa)
2310 {
2311   float t01, t02, x3, y3;
2312
2313   t01 = sqrtf(x1 * x1 + y1 * y1);
2314   t02 = sqrtf(x2 * x2 + y2 * y2);
2315   if (t01 == 0.0f) {
2316     t01 = 1.0f;
2317   }
2318   if (t02 == 0.0f) {
2319     t02 = 1.0f;
2320   }
2321
2322   x1 /= t01;
2323   y1 /= t01;
2324   x2 /= t02;
2325   y2 /= t02;
2326
2327   t02 = x1 * x2 + y1 * y2;
2328   if (fabsf(t02) >= 1.0f) {
2329     t02 = M_PI_2;
2330   }
2331   else {
2332     t02 = (saacos(t02)) / 2.0f;
2333   }
2334
2335   t02 = sinf(t02);
2336   if (t02 == 0.0f) {
2337     t02 = 1.0f;
2338   }
2339
2340   x3 = x1 - x2;
2341   y3 = y1 - y2;
2342   if (x3 == 0 && y3 == 0) {
2343     x3 = y1;
2344     y3 = -x1;
2345   }
2346   else {
2347     t01 = sqrtf(x3 * x3 + y3 * y3);
2348     x3 /= t01;
2349     y3 /= t01;
2350   }
2351
2352   *r_sina = -y3 / t02;
2353   *r_cosa = x3 / t02;
2354 }
2355
2356 static void tilt_bezpart(BezTriple *prevbezt,
2357                          BezTriple *bezt,
2358                          Nurb *nu,
2359                          float *tilt_array,
2360                          float *radius_array,
2361                          float *weight_array,
2362                          int resolu,
2363                          int stride)
2364 {
2365   BezTriple *pprev, *next, *last;
2366   float fac, dfac, t[4];
2367   int a;
2368
2369   if (tilt_array == NULL && radius_array == NULL) {
2370     return;
2371   }
2372
2373   last = nu->bezt + (nu->pntsu - 1);
2374
2375   /* returns a point */
2376   if (prevbezt == nu->bezt) {
2377     if (nu->flagu & CU_NURB_CYCLIC) {
2378       pprev = last;
2379     }
2380     else {
2381       pprev = prevbezt;
2382     }
2383   }
2384   else {
2385     pprev = prevbezt - 1;
2386   }
2387
2388   /* next point */
2389   if (bezt == last) {
2390     if (nu->flagu & CU_NURB_CYCLIC) {
2391       next = nu->bezt;
2392     }
2393     else {
2394       next = bezt;
2395     }
2396   }
2397   else {
2398     next = bezt + 1;
2399   }
2400
2401   fac = 0.0;
2402   dfac = 1.0f / (float)resolu;
2403
2404   for (a = 0; a < resolu; a++, fac += dfac) {
2405     if (tilt_array) {
2406       if (nu->tilt_interp ==
2407           KEY_CU_EASE) { /* May as well support for tilt also 2.47 ease interp */
2408         *tilt_array = prevbezt->tilt +
2409                       (bezt->tilt - prevbezt->tilt) * (3.0f * fac * fac - 2.0f * fac * fac * fac);
2410       }
2411       else {
2412         key_curve_position_weights(fac, t, nu->tilt_interp);
2413         *tilt_array = t[0] * pprev->tilt + t[1] * prevbezt->tilt + t[2] * bezt->tilt +
2414                       t[3] * next->tilt;
2415       }
2416
2417       tilt_array = POINTER_OFFSET(tilt_array, stride);
2418     }
2419
2420     if (radius_array) {
2421       if (nu->radius_interp == KEY_CU_EASE) {
2422         /* Support 2.47 ease interp
2423          * Note! - this only takes the 2 points into account,
2424          * giving much more localized results to changes in radius, sometimes you want that */
2425         *radius_array = prevbezt->radius + (bezt->radius - prevbezt->radius) *
2426                                                (3.0f * fac * fac - 2.0f * fac * fac * fac);
2427       }
2428       else {
2429
2430         /* reuse interpolation from tilt if we can */
2431         if (tilt_array == NULL || nu->tilt_interp != nu->radius_interp) {
2432           key_curve_position_weights(fac, t, nu->radius_interp);
2433         }
2434         *radius_array = t[0] * pprev->radius + t[1] * prevbezt->radius + t[2] * bezt->radius +
2435                         t[3] * next->radius;
2436       }
2437
2438       radius_array = POINTER_OFFSET(radius_array, stride);
2439     }
2440
2441     if (weight_array) {
2442       /* basic interpolation for now, could copy tilt interp too  */
2443       *weight_array = prevbezt->weight + (bezt->weight - prevbezt->weight) *
2444                                              (3.0f * fac * fac - 2.0f * fac * fac * fac);
2445
2446       weight_array = POINTER_OFFSET(weight_array, stride);
2447     }
2448   }
2449 }
2450
2451 /* make_bevel_list_3D_* funcs, at a minimum these must
2452  * fill in the bezp->quat and bezp->dir values */
2453
2454 /* utility for make_bevel_list_3D_* funcs */
2455 static void bevel_list_calc_bisect(BevList *bl)
2456 {
2457   BevPoint *bevp2, *bevp1, *bevp0;
2458   int nr;
2459   bool is_cyclic = bl->poly != -1;
2460
2461   if (is_cyclic) {
2462     bevp2 = bl->bevpoints;
2463     bevp1 = bevp2 + (bl->nr - 1);
2464     bevp0 = bevp1 - 1;
2465     nr = bl->nr;
2466   }
2467   else {
2468     /* If spline is not cyclic, direction of first and
2469      * last bevel points matches direction of CV handle.
2470      *
2471      * This is getting calculated earlier when we know
2472      * CV's handles and here we might simply skip evaluation
2473      * of direction for this guys.
2474      */
2475
2476     bevp0 = bl->bevpoints;
2477     bevp1 = bevp0 + 1;
2478     bevp2 = bevp1 + 1;
2479
2480     nr = bl->nr - 2;
2481   }
2482
2483   while (nr--) {
2484     /* totally simple */
2485     bisect_v3_v3v3v3(bevp1->dir, bevp0->vec, bevp1->vec, bevp2->vec);
2486
2487     bevp0 = bevp1;
2488     bevp1 = bevp2;
2489     bevp2++;
2490   }
2491 }
2492 static void bevel_list_flip_tangents(BevList *bl)
2493 {
2494   BevPoint *bevp2, *bevp1, *bevp0;
2495   int nr;
2496
2497   bevp2 = bl->bevpoints;
2498   bevp1 = bevp2 + (bl->nr - 1);
2499   bevp0 = bevp1 - 1;
2500
2501   nr = bl->nr;
2502   while (nr--) {
2503     if (angle_normalized_v3v3(bevp0->tan, bevp1->tan) > DEG2RADF(90.0f)) {
2504       negate_v3(bevp1->tan);
2505     }
2506
2507     bevp0 = bevp1;
2508     bevp1 = bevp2;
2509     bevp2++;
2510   }
2511 }
2512 /* apply user tilt */
2513 static void bevel_list_apply_tilt(BevList *bl)
2514 {
2515   BevPoint *bevp2, *bevp1;
2516   int nr;
2517   float q[4];
2518
2519   bevp2 = bl->bevpoints;
2520   bevp1 = bevp2 + (bl->nr - 1);
2521
2522   nr = bl->nr;
2523   while (nr--) {
2524     axis_angle_to_quat(q, bevp1->dir, bevp1->tilt);
2525     mul_qt_qtqt(bevp1->quat, q, bevp1->quat);
2526     normalize_qt(bevp1->quat);
2527
2528     bevp1 = bevp2;
2529     bevp2++;
2530   }
2531 }
2532 /* smooth quats, this function should be optimized, it can get slow with many iterations. */
2533 static void bevel_list_smooth(BevList *bl, int smooth_iter)
2534 {
2535   BevPoint *bevp2, *bevp1, *bevp0;
2536   int nr;
2537
2538   float q[4];
2539   float bevp0_quat[4];
2540   int a;
2541
2542   for (a = 0; a < smooth_iter; a++) {
2543     bevp2 = bl->bevpoints;
2544     bevp1 = bevp2 + (bl->nr - 1);
2545     bevp0 = bevp1 - 1;
2546
2547     nr = bl->nr;
2548
2549     if (bl->poly == -1) { /* check its not cyclic */
2550       /* skip the first point */
2551       /* bevp0 = bevp1; */
2552       bevp1 = bevp2;
2553       bevp2++;
2554       nr--;
2555
2556       bevp0 = bevp1;
2557       bevp1 = bevp2;
2558       bevp2++;
2559       nr--;
2560     }
2561
2562     copy_qt_qt(bevp0_quat, bevp0->quat);
2563
2564     while (nr--) {
2565       /* interpolate quats */
2566       float zaxis[3] = {0, 0, 1}, cross[3], q2[4];
2567       interp_qt_qtqt(q, bevp0_quat, bevp2->quat, 0.5);
2568       normalize_qt(q);
2569
2570       mul_qt_v3(q, zaxis);
2571       cross_v3_v3v3(cross, zaxis, bevp1->dir);
2572       axis_angle_to_quat(q2, cross, angle_normalized_v3v3(zaxis, bevp1->dir));
2573       normalize_qt(q2);
2574
2575       copy_qt_qt(bevp0_quat, bevp1->quat);
2576       mul_qt_qtqt(q, q2, q);
2577       interp_qt_qtqt(bevp1->quat, bevp1->quat, q, 0.5);
2578       normalize_qt(bevp1->quat);
2579
2580       /* bevp0 = bevp1; */ /* UNUSED */
2581       bevp1 = bevp2;
2582       bevp2++;
2583     }
2584   }
2585 }
2586
2587 static void make_bevel_list_3D_zup(BevList *bl)
2588 {
2589   BevPoint *bevp = bl->bevpoints;
2590   int nr = bl->nr;
2591
2592   bevel_list_calc_bisect(bl);
2593
2594   while (nr--) {
2595     vec_to_quat(bevp->quat, bevp->dir, 5, 1);
2596     bevp++;
2597   }
2598 }
2599
2600 static void minimum_twist_between_two_points(BevPoint *current_point, BevPoint *previous_point)
2601 {
2602   float angle = angle_normalized_v3v3(previous_point->dir, current_point->dir);
2603   float q[4];
2604
2605   if (angle > 0.0f) { /* otherwise we can keep as is */
2606     float cross_tmp[3];
2607     cross_v3_v3v3(cross_tmp, previous_point->dir, current_point->dir);
2608     axis_angle_to_quat(q, cross_tmp, angle);
2609     mul_qt_qtqt(current_point->quat, q, previous_point->quat);
2610   }
2611   else {
2612     copy_qt_qt(current_point->quat, previous_point->quat);
2613   }
2614 }
2615
2616 static void make_bevel_list_3D_minimum_twist(BevList *bl)
2617 {
2618   BevPoint *bevp2, *bevp1, *bevp0; /* standard for all make_bevel_list_3D_* funcs */
2619   int nr;
2620   float q[4];
2621
2622   bevel_list_calc_bisect(bl);
2623
2624   bevp2 = bl->bevpoints;
2625   bevp1 = bevp2 + (bl->nr - 1);
2626   bevp0 = bevp1 - 1;
2627
2628   nr = bl->nr;
2629   while (nr--) {
2630
2631     if (nr + 4 > bl->nr) { /* first time and second time, otherwise first point adjusts last */
2632       vec_to_quat(bevp1->quat, bevp1->dir, 5, 1);
2633     }
2634     else {
2635       minimum_twist_between_two_points(bevp1, bevp0);
2636     }
2637
2638     bevp0 = bevp1;
2639     bevp1 = bevp2;
2640     bevp2++;
2641   }
2642
2643   if (bl->poly != -1) { /* check for cyclic */
2644
2645     /* Need to correct for the start/end points not matching
2646      * do this by calculating the tilt angle difference, then apply
2647      * the rotation gradually over the entire curve
2648      *
2649      * note that the split is between last and second last, rather than first/last as youd expect.
2650      *
2651      * real order is like this
2652      * 0,1,2,3,4 --> 1,2,3,4,0
2653      *
2654      * this is why we compare last with second last
2655      * */
2656     float vec_1[3] = {0, 1, 0}, vec_2[3] = {0, 1, 0}, angle, ang_fac, cross_tmp[3];
2657
2658     BevPoint *bevp_first;
2659     BevPoint *bevp_last;
2660
2661     bevp_first = bl->bevpoints;
2662     bevp_first += bl->nr - 1;
2663     bevp_last = bevp_first;
2664     bevp_last--;
2665
2666     /* quats and vec's are normalized, should not need to re-normalize */
2667     mul_qt_v3(bevp_first->quat, vec_1);
2668     mul_qt_v3(bevp_last->quat, vec_2);
2669     normalize_v3(vec_1);
2670     normalize_v3(vec_2);
2671
2672     /* align the vector, can avoid this and it looks 98% OK but
2673      * better to align the angle quat roll's before comparing */
2674     {
2675       cross_v3_v3v3(cross_tmp, bevp_last->dir, bevp_first->dir);
2676       angle = angle_normalized_v3v3(bevp_first->dir, bevp_last->dir);
2677       axis_angle_to_quat(q, cross_tmp, angle);
2678       mul_qt_v3(q, vec_2);
2679     }
2680
2681     angle = angle_normalized_v3v3(vec_1, vec_2);
2682
2683     /* flip rotation if needs be */
2684     cross_v3_v3v3(cross_tmp, vec_1, vec_2);
2685     normalize_v3(cross_tmp);
2686     if (angle_normalized_v3v3(bevp_first->dir, cross_tmp) < DEG2RADF(90.0f)) {
2687       angle = -angle;
2688     }
2689
2690     bevp2 = bl->bevpoints;
2691     bevp1 = bevp2 + (bl->nr - 1);
2692     bevp0 = bevp1 - 1;
2693
2694     nr = bl->nr;
2695     while (nr--) {
2696       ang_fac = angle * (1.0f - ((float)nr / bl->nr)); /* also works */
2697
2698       axis_angle_to_quat(q, bevp1->dir, ang_fac);
2699       mul_qt_qtqt(bevp1->quat, q, bevp1->quat);
2700
2701       bevp0 = bevp1;
2702       bevp1 = bevp2;
2703       bevp2++;
2704     }
2705   }
2706   else {
2707     /* Need to correct quat for the first/last point,
2708      * this is so because previously it was only calculated
2709      * using it's own direction, which might not correspond
2710      * the twist of neighbor point.
2711      */
2712     bevp1 = bl->bevpoints;
2713     bevp0 = bevp1 + 1;
2714     minimum_twist_between_two_points(bevp1, bevp0);
2715
2716     bevp2 = bl->bevpoints;
2717     bevp1 = bevp2 + (bl->nr - 1);
2718     bevp0 = bevp1 - 1;
2719     minimum_twist_between_two_points(bevp1, bevp0);
2720   }
2721 }
2722
2723 static void make_bevel_list_3D_tangent(BevList *bl)
2724 {
2725   BevPoint *bevp2, *bevp1, *bevp0; /* standard for all make_bevel_list_3D_* funcs */
2726   int nr;
2727
2728   float bevp0_tan[3];
2729
2730   bevel_list_calc_bisect(bl);
2731   bevel_list_flip_tangents(bl);
2732
2733   /* correct the tangents */
2734   bevp2 = bl->bevpoints;
2735   bevp1 = bevp2 + (bl->nr - 1);
2736   bevp0 = bevp1 - 1;
2737
2738   nr = bl->nr;
2739   while (nr--) {
2740     float cross_tmp[3];
2741     cross_v3_v3v3(cross_tmp, bevp1->tan, bevp1->dir);
2742     cross_v3_v3v3(bevp1->tan, cross_tmp, bevp1->dir);
2743     normalize_v3(bevp1->tan);
2744
2745     bevp0 = bevp1;
2746     bevp1 = bevp2;
2747     bevp2++;
2748   }
2749
2750   /* now for the real twist calc */
2751   bevp2 = bl->bevpoints;
2752   bevp1 = bevp2 + (bl->nr - 1);
2753   bevp0 = bevp1 - 1;
2754
2755   copy_v3_v3(bevp0_tan, bevp0->tan);
2756
2757   nr = bl->nr;
2758   while (nr--) {
2759     /* make perpendicular, modify tan in place, is ok */
2760     float cross_tmp[3];
2761     float zero[3] = {0, 0, 0};
2762
2763     cross_v3_v3v3(cross_tmp, bevp1->tan, bevp1->dir);
2764     normalize_v3(cross_tmp);
2765     tri_to_quat(bevp1->quat, zero, cross_tmp, bevp1->tan); /* XXX - could be faster */
2766
2767     /* bevp0 = bevp1; */ /* UNUSED */
2768     bevp1 = bevp2;
2769     bevp2++;
2770   }
2771 }
2772
2773 static void make_bevel_list_3D(BevList *bl, int smooth_iter, int twist_mode)
2774 {
2775   switch (twist_mode) {
2776     case CU_TWIST_TANGENT:
2777       make_bevel_list_3D_tangent(bl);
2778       break;
2779     case CU_TWIST_MINIMUM:
2780       make_bevel_list_3D_minimum_twist(bl);
2781       break;
2782     default: /* CU_TWIST_Z_UP default, pre 2.49c */
2783       make_bevel_list_3D_zup(bl);
2784       break;
2785   }
2786
2787   if (smooth_iter) {
2788     bevel_list_smooth(bl, smooth_iter);
2789   }
2790
2791   bevel_list_apply_tilt(bl);
2792 }
2793
2794 /* only for 2 points */
2795 static void make_bevel_list_segment_3D(BevList *bl)
2796 {
2797   float q[4];
2798
2799   BevPoint *bevp2 = bl->bevpoints;
2800   BevPoint *bevp1 = bevp2 + 1;
2801
2802   /* simple quat/dir */
2803   sub_v3_v3v3(bevp1->dir, bevp1->vec, bevp2->vec);
2804   normalize_v3(bevp1->dir);
2805
2806   vec_to_quat(bevp1->quat, bevp1->dir, 5, 1);
2807
2808   axis_angle_to_quat(q, bevp1->dir, bevp1->tilt);
2809   mul_qt_qtqt(bevp1->quat, q, bevp1->quat);
2810   normalize_qt(bevp1->quat);
2811   copy_v3_v3(bevp2->dir, bevp1->dir);
2812   copy_qt_qt(bevp2->quat, bevp1->quat);
2813 }
2814
2815 /* only for 2 points */
2816 static void make_bevel_list_segment_2D(BevList *bl)
2817 {
2818   BevPoint *bevp2 = bl->bevpoints;
2819   BevPoint *bevp1 = bevp2 + 1;
2820
2821   const float x1 = bevp1->vec[0] - bevp2->vec[0];
2822   const float y1 = bevp1->vec[1] - bevp2->vec[1];
2823
2824   calc_bevel_sin_cos(x1, y1, -x1, -y1, &(bevp1->sina), &(bevp1->cosa));
2825   bevp2->sina = bevp1->sina;
2826   bevp2->cosa = bevp1->cosa;
2827
2828   /* fill in dir & quat */
2829   make_bevel_list_segment_3D(bl);
2830 }
2831
2832 static void make_bevel_list_2D(BevList *bl)
2833 {
2834   /* note: bevp->dir and bevp->quat are not needed for beveling but are
2835    * used when making a path from a 2D curve, therefore they need to be set - Campbell */
2836
2837   BevPoint *bevp0, *bevp1, *bevp2;
2838   int nr;
2839
2840   if (bl->poly != -1) {
2841     bevp2 = bl->bevpoints;
2842     bevp1 = bevp2 + (bl->nr - 1);
2843     bevp0 = bevp1 - 1;
2844     nr = bl->nr;
2845   }
2846   else {
2847     bevp0 = bl->bevpoints;
2848     bevp1 = bevp0 + 1;
2849     bevp2 = bevp1 + 1;
2850
2851     nr = bl->nr - 2;
2852   }
2853
2854   while (nr--) {
2855     const float x1 = bevp1->vec[0] - bevp0->vec[0];
2856     const float x2 = bevp1->vec[0] - bevp2->vec[0];
2857     const float y1 = bevp1->vec[1] - bevp0->vec[1];
2858     const float y2 = bevp1->vec[1] - bevp2->vec[1];
2859
2860     calc_bevel_sin_cos(x1, y1, x2, y2, &(bevp1->sina), &(bevp1->cosa));
2861
2862     /* from: make_bevel_list_3D_zup, could call but avoid a second loop.
2863      * no need for tricky tilt calculation as with 3D curves */
2864     bisect_v3_v3v3v3(bevp1->dir, bevp0->vec, bevp1->vec, bevp2->vec);
2865     vec_to_quat(bevp1->quat, bevp1->dir, 5, 1);
2866     /* done with inline make_bevel_list_3D_zup */
2867
2868     bevp0 = bevp1;
2869     bevp1 = bevp2;
2870     bevp2++;
2871   }
2872
2873   /* correct non-cyclic cases */
2874   if (bl->poly == -1) {
2875     BevPoint *bevp;
2876     float angle;
2877
2878     /* first */
2879     bevp = bl->bevpoints;
2880     angle = atan2f(bevp->dir[0], bevp->dir[1]) - (float)M_PI_2;
2881     bevp->sina = sinf(angle);
2882     bevp->cosa = cosf(angle);
2883     vec_to_quat(bevp->quat, bevp->dir, 5, 1);
2884
2885     /* last */
2886     bevp = bl->bevpoints;
2887     bevp += (bl->nr - 1);
2888     angle = atan2f(bevp->dir[0], bevp->dir[1]) - (float)M_PI_2;
2889     bevp->sina = sinf(angle);
2890     bevp->cosa = cosf(angle);
2891     vec_to_quat(bevp->quat, bevp->dir, 5, 1);
2892   }
2893 }
2894
2895 static void bevlist_firstlast_direction_calc_from_bpoint(Nurb *nu, BevList *bl)
2896 {
2897   if (nu->pntsu > 1) {
2898     BPoint *first_bp = nu->bp, *last_bp = nu->bp + (nu->pntsu - 1);
2899     BevPoint *first_bevp, *last_bevp;
2900
2901     first_bevp = bl->bevpoints;
2902     last_bevp = first_bevp + (bl->nr - 1);
2903
2904     sub_v3_v3v3(first_bevp->dir, (first_bp + 1)->vec, first_bp->vec);
2905     normalize_v3(first_bevp->dir);
2906
2907     sub_v3_v3v3(last_bevp->dir, last_bp->vec, (last_bp - 1)->vec);
2908     normalize_v3(last_bevp->dir);
2909   }
2910 }
2911
2912 void BKE_curve_bevelList_free(ListBase *bev)
2913 {
2914   BevList *bl, *blnext;
2915   for (bl = bev->first; bl != NULL; bl = blnext) {
2916     blnext = bl->next;
2917     if (bl->seglen != NULL) {
2918       MEM_freeN(bl->seglen);
2919     }
2920     if (bl->segbevcount != NULL) {
2921       MEM_freeN(bl->segbevcount);
2922     }
2923     if (bl->bevpoints != NULL) {
2924       MEM_freeN(bl->bevpoints);
2925     }
2926     MEM_freeN(bl);
2927   }
2928
2929   BLI_listbase_clear(bev);
2930 }
2931
2932 void BKE_curve_bevelList_make(Object *ob, ListBase *nurbs, bool for_render)
2933 {
2934   /*
2935    * - convert all curves to polys, with indication of resol and flags for double-vertices
2936    * - possibly; do a smart vertice removal (in case Nurb)
2937    * - separate in individual blocks with BoundBox
2938    * - AutoHole detection
2939    */
2940
2941   /* this function needs an object, because of tflag and upflag */
2942   Curve *cu = ob->data;
2943   Nurb *nu;
2944   BezTriple *bezt, *prevbezt;
2945   BPoint *bp;
2946   BevList *bl, *blnew, *blnext;
2947   BevPoint *bevp2, *bevp1 = NULL, *bevp0;
2948   const float treshold = 0.00001f;
2949   float min, inp;
2950   float *seglen = NULL;
2951   struct BevelSort *sortdata, *sd, *sd1;
2952   int a, b, nr, poly, resolu = 0, len = 0, segcount;
2953   int *segbevcount;
2954   bool do_tilt, do_radius, do_weight;
2955   bool is_editmode = false;
2956   ListBase *bev;
2957
2958   /* segbevcount alsp requires seglen. */
2959   const bool need_seglen = ELEM(
2960                                cu->bevfac1_mapping, CU_BEVFAC_MAP_SEGMENT, CU_BEVFAC_MAP_SPLINE) ||
2961                            ELEM(cu->bevfac2_mapping, CU_BEVFAC_MAP_SEGMENT, CU_BEVFAC_MAP_SPLINE);
2962
2963   bev = &ob->runtime.curve_cache->bev;
2964
2965   /* do we need to calculate the radius for each point? */
2966   /* do_radius = (cu->bevobj || cu->taperobj || (cu->flag & CU_FRONT) || (cu->flag & CU_BACK)) ? 0 : 1; */
2967
2968   /* STEP 1: MAKE POLYS  */
2969
2970   BKE_curve_bevelList_free(&ob->runtime.curve_cache->bev);
2971   nu = nurbs->first;
2972   if (cu->editnurb && ob->type != OB_FONT) {
2973     is_editmode = 1;
2974   }
2975
2976   for (; nu; nu = nu->next) {
2977
2978     if (nu->hide && is_editmode) {
2979       continue;
2980     }
2981
2982     /* check if we will calculate tilt data */
2983     do_tilt = CU_DO_TILT(cu, nu);
2984     do_radius = CU_DO_RADIUS(
2985         cu, nu); /* normal display uses the radius, better just to calculate them */
2986     do_weight = true;
2987
2988     /* check we are a single point? also check we are not a surface and that the orderu is sane,
2989      * enforced in the UI but can go wrong possibly */
2990     if (!BKE_nurb_check_valid_u(nu)) {
2991       bl = MEM_callocN(sizeof(BevList), "makeBevelList1");
2992       bl->bevpoints = MEM_calloc_arrayN(1, sizeof(BevPoint), "makeBevelPoints1");
2993       BLI_addtail(bev, bl);
2994       bl->nr = 0;
2995       bl->charidx = nu->charidx;
2996     }
2997     else {
2998       BevPoint *bevp;
2999
3000       if (for_render && cu->resolu_ren != 0) {
3001         resolu = cu->resolu_ren;
3002       }
3003       else {
3004         resolu = nu->resolu;
3005       }
3006
3007       segcount = SEGMENTSU(nu);
3008
3009       if (nu->type == CU_POLY) {
3010         len = nu->pntsu;
3011         bl = MEM_callocN(sizeof(BevList), "makeBevelList2");
3012         bl->bevpoints = MEM_calloc_arrayN(len, sizeof(BevPoint), "makeBevelPoints2");
3013         if (need_seglen && (nu->flagu & CU_NURB_CYCLIC) == 0) {
3014           bl->seglen = MEM_malloc_arrayN(segcount, sizeof(float), "makeBevelList2_seglen");
3015           bl->segbevcount = MEM_malloc_arrayN(segcount, sizeof(int), "makeBevelList2_segbevcount");
3016         }
3017         BLI_addtail(bev, bl);
3018
3019         bl->poly = (nu->flagu & CU_NURB_CYCLIC) ? 0 : -1;
3020         bl->nr = len;
3021         bl->dupe_nr = 0;
3022         bl->charidx = nu->charidx;
3023         bevp = bl->bevpoints;
3024         bevp->offset = 0;
3025         bp = nu->bp;
3026         seglen = bl->seglen;
3027         segbevcount = bl->segbevcount;
3028
3029         while (len--) {
3030           copy_v3_v3(bevp->vec, bp->vec);
3031           bevp->tilt = bp->tilt;
3032           bevp->radius = bp->radius;
3033           bevp->weight = bp->weight;
3034           bevp->split_tag = true;
3035           bp++;
3036           if (seglen != NULL && len != 0) {
3037             *seglen = len_v3v3(bevp->vec, bp->vec);
3038             bevp++;
3039             bevp->offset = *seglen;
3040             if (*seglen > treshold) {
3041               *segbevcount = 1;
3042             }
3043             else {
3044               *segbevcount = 0;
3045             }
3046             seglen++;
3047             segbevcount++;
3048           }
3049           else {
3050             bevp++;
3051           }
3052         }
3053
3054         if ((nu->flagu & CU_NURB_CYCLIC) == 0) {
3055           bevlist_firstlast_direction_calc_from_bpoint(nu, bl);
3056         }
3057       }
3058       else if (nu->type == CU_BEZIER) {
3059         /* in case last point is not cyclic */
3060         len = segcount * resolu + 1;
3061
3062         bl = MEM_callocN(sizeof(BevList), "makeBevelBPoints");
3063         bl->bevpoints = MEM_calloc_arrayN(len, sizeof(BevPoint), "makeBevelBPointsPoints");
3064         if (need_seglen && (nu->flagu & CU_NURB_CYCLIC) == 0) {
3065           bl->seglen = MEM_malloc_arrayN(segcount, sizeof(float), "makeBevelBPoints_seglen");
3066           bl->segbevcount = MEM_malloc_arrayN(
3067               segcount, sizeof(int), "makeBevelBPoints_segbevcount");
3068         }
3069         BLI_addtail(bev, bl);
3070
3071         bl->poly = (nu->flagu & CU_NURB_CYCLIC) ? 0 : -1;
3072         bl->charidx = nu->charidx;
3073
3074         bevp = bl->bevpoints;
3075         seglen = bl->seglen;
3076         segbevcount = bl->segbevcount;
3077
3078         bevp->offset = 0;
3079         if (seglen != NULL) {
3080           *seglen = 0;
3081           *segbevcount = 0;
3082         }
3083
3084         a = nu->pntsu - 1;
3085         bezt = nu->bezt;
3086         if (nu->flagu & CU_NURB_CYCLIC) {
3087           a++;
3088           prevbezt = nu->bezt + (nu->pntsu - 1);
3089         }
3090         else {
3091           prevbezt = bezt;
3092           bezt++;
3093         }
3094
3095         sub_v3_v3v3(bevp->dir, prevbezt->vec[2], prevbezt->vec[1]);
3096         normalize_v3(bevp->dir);
3097
3098         BLI_assert(segcount >= a);
3099
3100         while (a--) {
3101           if (prevbezt->h2 == HD_VECT && bezt->h1 == HD_VECT) {
3102
3103             copy_v3_v3(bevp->vec, prevbezt->vec[1]);
3104             bevp->tilt = prevbezt->tilt;
3105             bevp->radius = prevbezt->radius;
3106             bevp->weight = prevbezt->weight;
3107             bevp->split_tag = true;
3108             bevp->dupe_tag = false;
3109             bevp++;
3110             bl->nr++;
3111             bl->dupe_nr = 1;
3112             if (seglen != NULL) {
3113               *seglen = len_v3v3(prevbezt->vec[1], bezt->vec[1]);
3114               bevp->offset = *seglen;
3115               seglen++;
3116               /* match segbevcount to the cleaned up bevel lists (see STEP 2) */
3117               if (bevp->offset > treshold) {
3118                 *segbevcount = 1;
3119               }
3120               segbevcount++;
3121             }
3122           }
3123           else {
3124             /* always do all three, to prevent data hanging around */
3125             int j;
3126
3127             /* BevPoint must stay aligned to 4 so sizeof(BevPoint)/sizeof(float) works */
3128             for (j = 0; j < 3; j++) {
3129               BKE_curve_forward_diff_bezier(prevbezt->vec[1][j],
3130                                             prevbezt->vec[2][j],
3131                                             bezt->vec[0][j],
3132                                             bezt->vec[1][j],
3133                                             &(bevp->vec[j]),
3134                                             resolu,
3135                                             sizeof(BevPoint));
3136             }
3137
3138             /* if both arrays are NULL do nothiong */
3139             tilt_bezpart(prevbezt,
3140                          bezt,
3141                          nu,
3142                          do_tilt ? &bevp->tilt : NULL,
3143                          do_radius ? &bevp->radius : NULL,
3144                          do_weight ? &bevp->weight : NULL,
3145                          resolu,
3146                          sizeof(BevPoint));
3147
3148             if (cu->twist_mode == CU_TWIST_TANGENT) {
3149               forward_diff_bezier_cotangent(prevbezt->vec[1],
3150                                             prevbezt->vec[2],
3151                                             bezt->vec[0],
3152                                             bezt->vec[1],
3153                                             bevp->tan,
3154                                             resolu,
3155                                             sizeof(BevPoint));
3156             }
3157
3158             /* indicate with handlecodes double points */
3159             if (prevbezt->h1 == prevbezt->h2) {
3160               if (prevbezt->h1 == 0 || prevbezt->h1 == HD_VECT) {
3161                 bevp->split_tag = true;
3162               }
3163             }
3164             else {
3165               if (prevbezt->h1 == 0 || prevbezt->h1 == HD_VECT) {
3166                 bevp->split_tag = true;
3167               }
3168               else if (prevbezt->h2 == 0 || prevbezt->h2 == HD_VECT) {
3169                 bevp->split_tag = true;
3170               }
3171             }
3172
3173             /* seglen */
3174             if (seglen != NULL) {
3175               *seglen = 0;
3176               *segbevcount = 0;
3177               for (j = 0; j < resolu; j++) {
3178                 bevp0 = bevp;
3179                 bevp++;
3180                 bevp->offset = len_v3v3(bevp0->vec, bevp->vec);
3181                 /* match seglen and segbevcount to the cleaned up bevel lists (see STEP 2) */
3182                 if (bevp->offset > treshold) {
3183                   *seglen += bevp->offset;
3184                   *segbevcount += 1;
3185                 }
3186               }
3187               seglen++;
3188               segbevcount++;
3189             }
3190             else {
3191               bevp += resolu;
3192             }
3193             bl->nr += resolu;
3194           }
3195           prevbezt = bezt;
3196           bezt++;
3197         }
3198
3199         if ((nu->flagu & CU_NURB_CYCLIC) == 0) { /* not cyclic: endpoint */
3200           copy_v3_v3(bevp->vec, prevbezt->vec[1]);
3201           bevp->tilt = prevbezt->tilt;
3202           bevp->radius = prevbezt->radius;
3203           bevp->weight = prevbezt->weight;
3204
3205           sub_v3_v3v3(bevp->dir, prevbezt->vec[1], prevbezt->vec[0]);
3206           normalize_v3(bevp->dir);
3207
3208           bl->nr++;
3209         }
3210       }
3211       else if (nu->type == CU_NURBS) {
3212         if (nu->pntsv == 1) {
3213           len = (resolu * segcount);
3214
3215           bl = MEM_callocN(sizeof(BevList), "makeBevelList3");
3216           bl->bevpoints = MEM_calloc_arrayN(len, sizeof(BevPoint), "makeBevelPoints3");
3217           if (need_seglen && (nu->flagu & CU_NURB_CYCLIC) == 0) {
3218             bl->seglen = MEM_malloc_arrayN(segcount, sizeof(float), "makeBevelList3_seglen");
3219             bl->segbevcount = MEM_malloc_arrayN(
3220                 segcount, sizeof(int), "makeBevelList3_segbevcount");
3221           }
3222           BLI_addtail(bev, bl);
3223           bl->nr = len;
3224           bl->dupe_nr = 0;
3225           bl->poly = (nu->flagu & CU_NURB_CYCLIC) ? 0 : -1;
3226           bl->charidx = nu->charidx;
3227
3228           bevp = bl->bevpoints;
3229           seglen = bl->seglen;
3230           segbevcount = bl->segbevcount;
3231
3232           BKE_nurb_makeCurve(nu,
3233                              &bevp->vec[0],
3234                              do_tilt ? &bevp->tilt : NULL,
3235                              do_radius ? &bevp->radius : NULL,
3236                              do_weight ? &bevp->weight : NULL,
3237                              resolu,
3238                              sizeof(BevPoint));
3239
3240           /* match seglen and segbevcount to the cleaned up bevel lists (see STEP 2) */
3241           if (seglen != NULL) {
3242             nr = segcount;
3243             bevp0 = bevp;
3244             bevp++;
3245             while (nr) {
3246               int j;
3247               *seglen = 0;
3248               *segbevcount = 0;
3249               /* We keep last bevel segment zero-length. */
3250               for (j = 0; j < ((nr == 1) ? (resolu - 1) : resolu); j++) {
3251                 bevp->offset = len_v3v3(bevp0->vec, bevp->vec);
3252                 if (bevp->offset > treshold) {
3253                   *seglen += bevp->offset;
3254                   *segbevcount += 1;
3255                 }
3256                 bevp0 = bevp;
3257                 bevp++;
3258               }
3259               seglen++;
3260               segbevcount++;
3261               nr--;
3262             }
3263           }
3264
3265           if ((nu->flagu & CU_NURB_CYCLIC) == 0) {
3266             bevlist_firstlast_direction_calc_from_bpoint(nu, bl);
3267           }
3268         }
3269       }
3270     }
3271   }
3272
3273   /* STEP 2: DOUBLE POINTS AND AUTOMATIC RESOLUTION, REDUCE DATABLOCKS */
3274   bl = bev->first;
3275   while (bl) {
3276     if (bl->nr) { /* null bevel items come from single points */
3277       bool is_cyclic = bl->poly != -1;
3278       nr = bl->nr;
3279       if (is_cyclic) {
3280         bevp1 = bl->bevpoints;
3281         bevp0 = bevp1 + (nr - 1);
3282       }
3283       else {
3284         bevp0 = bl->bevpoints;
3285         bevp0->offset = 0;
3286         bevp1 = bevp0 + 1;
3287       }
3288       nr--;
3289       while (nr--) {
3290         if (seglen != NULL) {
3291           if (fabsf(bevp1->offset) < treshold) {
3292             bevp0->dupe_tag = true;
3293             bl->dupe_nr++;
3294           }
3295         }
3296         else {
3297           if (fabsf(bevp0->vec[0] - bevp1->vec[0]) < 0.00001f) {
3298             if (fabsf(bevp0->vec[1] - bevp1->vec[1]) < 0.00001f) {
3299               if (fabsf(bevp0->vec[2] - bevp1->vec[2]) < 0.00001f) {
3300                 bevp0->dupe_tag = true;
3301                 bl->dupe_nr++;
3302               }
3303             }
3304           }
3305         }
3306         bevp0 = bevp1;
3307         bevp1++;
3308       }
3309     }
3310     bl = bl->next;
3311   }
3312   bl = bev->first;
3313   while (bl) {
3314     blnext = bl->next;
3315     if (bl->nr && bl->dupe_nr) {
3316       nr = bl->nr - bl->dupe_nr + 1; /* +1 because vectorbezier sets flag too */
3317       blnew = MEM_callocN(sizeof(BevList), "makeBevelList4");
3318       memcpy(blnew, bl, sizeof(BevList));
3319       blnew->bevpoints = MEM_calloc_arrayN(nr, sizeof(BevPoint), "makeBevelPoints4");
3320       if (!blnew->bevpoints) {
3321         MEM_freeN(blnew);
3322         break;
3323       }
3324       blnew->segbevcount = bl->segbevcount;
3325       blnew->seglen = bl->seglen;
3326       blnew->nr = 0;
3327       BLI_remlink(bev, bl);
3328       BLI_insertlinkbefore(bev, blnext, blnew); /* to make sure bevlijst is tuned with nurblist */
3329       bevp0 = bl->bevpoints;
3330       bevp1 = blnew->bevpoints;
3331       nr = bl->nr;
3332       while (nr--) {
3333         if (bevp0->dupe_tag == 0) {
3334           memcpy(bevp1, bevp0, sizeof(BevPoint));
3335           bevp1++;
3336           blnew->nr++;
3337         }
3338         bevp0++;
3339       }
3340       if (bl->bevpoints != NULL) {
3341         MEM_freeN(bl->bevpoints);
3342       }
3343       MEM_freeN(bl);
3344       blnew->dupe_nr = 0;
3345     }
3346     bl = blnext;
3347   }
3348
3349   /* STEP 3: POLYS COUNT AND AUTOHOLE */
3350   bl = bev->first;
3351   poly = 0;
3352   while (bl) {
3353     if (bl->nr && bl->poly >= 0) {
3354       poly++;
3355       bl->poly = poly;
3356       bl->hole = 0;
3357     }
3358     bl = bl->next;
3359   }
3360
3361   /* find extreme left points, also test (turning) direction */
3362   if (poly > 0) {
3363     sd = sortdata = MEM_malloc_arrayN(poly, sizeof(struct BevelSort), "makeBevelList5");
3364     bl = bev->first;
3365     while (bl) {
3366       if (bl->poly > 0) {
3367         BevPoint *bevp;
3368
3369         min = 300000.0;
3370         bevp = bl->bevpoints;
3371         nr = bl->nr;
3372         while (nr--) {
3373           if (min > bevp->vec[0]) {
3374             min = bevp->vec[0];
3375             bevp1 = bevp;
3376           }
3377           bevp++;
3378         }
3379         sd->bl = bl;
3380         sd->left = min;
3381
3382         bevp = bl->bevpoints;
3383         if (bevp1 == bevp) {
3384           bevp0 = bevp + (bl->nr - 1);
3385         }
3386         else {
3387           bevp0 = bevp1 - 1;
3388         }
3389         bevp = bevp + (bl->nr - 1);
3390         if (bevp1 == bevp) {
3391           bevp2 = bl->bevpoints;
3392         }
3393         else {
3394           bevp2 = bevp1 + 1;
3395         }
3396
3397         inp = ((bevp1->vec[0] - bevp0->vec[0]) * (bevp0->vec[1] - bevp2->vec[1]) +
3398                (bevp0->vec[1] - bevp1->vec[1]) * (bevp0->vec[0] - bevp2->vec[0]));
3399
3400         if (inp > 0.0f) {
3401           sd->dir = 1;
3402         }
3403         else {
3404           sd->dir = 0;
3405         }
3406
3407         sd++;
3408       }
3409
3410       bl = bl->next;
3411     }
3412     qsort(sortdata, poly, sizeof(struct BevelSort), vergxcobev);
3413
3414     sd = sortdata + 1;
3415     for (a = 1; a < poly; a++, sd++) {
3416       bl = sd->bl; /* is bl a hole? */
3417       sd1 = sortdata + (a - 1);
3418       for (b = a - 1; b >= 0; b--, sd1--) {    /* all polys to the left */
3419         if (sd1->bl->charidx == bl->charidx) { /* for text, only check matching char */
3420           if (bevelinside(sd1->bl, bl)) {
3421             bl->hole = 1 - sd1->bl->hole;
3422             break;
3423           }
3424         }
3425       }
3426     }
3427
3428     /* turning direction */
3429     if ((cu->flag & CU_3D) == 0) {
3430       sd = sortdata;
3431       for (a = 0; a < poly; a++, sd++) {
3432         if (sd->bl->hole == sd->dir) {
3433           bl = sd->bl;
3434           bevp1 = bl->bevpoints;
3435           bevp2 = bevp1 + (bl->nr - 1);
3436           nr = bl->nr / 2;
3437           while (nr--) {
3438             SWAP(BevPoint, *bevp1, *bevp2);
3439             bevp1++;
3440             bevp2--;
3441           }
3442         }
3443       }
3444     }
3445     MEM_freeN(sortdata);
3446   }
3447
3448   /* STEP 4: 2D-COSINES or 3D ORIENTATION */
3449   if ((cu->flag & CU_3D) == 0) {
3450     /* 2D Curves */
3451     for (bl = bev->first; bl; bl = bl->next) {
3452       if (bl->nr < 2) {
3453         BevPoint *bevp = bl->bevpoints;
3454         unit_qt(bevp->quat);
3455       }
3456       else if (bl->nr == 2) { /* 2 pnt, treat separate */
3457         make_bevel_list_segment_2D(bl);
3458       }
3459       else {
3460         make_bevel_list_2D(bl);
3461       }
3462     }
3463   }
3464   else {
3465     /* 3D Curves */
3466     for (bl = bev->first; bl; bl = bl->next) {
3467       if (bl->nr < 2) {
3468         BevPoint *bevp = bl->bevpoints;
3469         unit_qt(bevp->quat);
3470       }
3471       else if (bl->nr == 2) { /* 2 pnt, treat separate */
3472         make_bevel_list_segment_3D(bl);
3473       }
3474       else {
3475         make_bevel_list_3D(bl, (int)(resolu * cu->twist_smooth), cu->twist_mode);
3476       }
3477     }
3478   }
3479 }
3480
3481 /* ****************** HANDLES ************** */
3482
3483 static void calchandleNurb_intern(BezTriple *bezt,
3484                                   const BezTriple *prev,
3485                                   const BezTriple *next,
3486                                   bool is_fcurve,
3487                                   bool skip_align,
3488                                   char fcurve_smoothing)
3489 {
3490   /* defines to avoid confusion */
3491 #define p2_h1 ((p2)-3)
3492 #define p2_h2 ((p2) + 3)
3493
3494   const float *p1, *p3;
3495   float *p2;
3496   float pt[3];
3497   float dvec_a[3], dvec_b[3];
3498   float len, len_a, len_b;
3499   float len_ratio;
3500   const float eps = 1e-5;
3501
3502   /* assume normal handle until we check */
3503   bezt->f5 = HD_AUTOTYPE_NORMAL;
3504
3505   if (bezt->h1 == 0 && bezt->h2 == 0) {
3506     return;
3507   }
3508
3509   p2 = bezt->vec[1];
3510
3511   if (prev == NULL) {
3512     p3 = next->vec[1];
3513     pt[0] = 2.0f * p2[0] - p3[0];
3514     pt[1] = 2.0f * p2[1] - p3[1];
3515     pt[2] = 2.0f * p2[2] - p3[2];
3516     p1 = pt;
3517   }
3518   else {
3519     p1 = prev->vec[1];
3520   }
3521
3522   if (next == NULL) {
3523     pt[0] = 2.0f * p2[0] - p1[0];
3524     pt[1] = 2.0f * p2[1] - p1[1];
3525     pt[2] = 2.0f * p2[2] - p1[2];
3526     p3 = pt;
3527   }
3528   else {
3529     p3 = next->vec[1];
3530   }
3531
3532   sub_v3_v3v3(dvec_a, p2, p1);
3533   sub_v3_v3v3(dvec_b, p3, p2);
3534
3535   if (is_fcurve) {
3536     len_a = dvec_a[0];
3537     len_b = dvec_b[0];
3538   }
3539   else {
3540     len_a = len_v3(dvec_a);
3541     len_b = len_v3(dvec_b);
3542   }
3543
3544   if (len_a == 0.0f) {
3545     len_a = 1.0f;
3546   }
3547   if (len_b == 0.0f) {
3548     len_b = 1.0f;
3549   }
3550
3551   len_ratio = len_a / len_b;
3552
3553   if (ELEM(bezt->h1, HD_AUTO, HD_AUTO_ANIM) || ELEM(bezt->h2, HD_AUTO, HD_AUTO_ANIM)) { /* auto */
3554     float tvec[3];
3555     tvec[0] = dvec_b[0] / len_b + dvec_a[0] / len_a;
3556     tvec[1] = dvec_b[1] / len_b + dvec_a[1] / len_a;
3557     tvec[2] = dvec_b[2] / len_b + dvec_a[2] / len_a;
3558
3559     if (is_fcurve) {
3560       if (fcurve_smoothing != FCURVE_SMOOTH_NONE) {
3561         /* force the horizontal handle size to be 1/3 of the key interval so that
3562          * the X component of the parametric bezier curve is a linear spline */
3563         len = 6.0f / 2.5614f;
3564       }
3565       else {
3566         len = tvec[0];
3567       }
3568     }
3569     else {
3570       len = len_v3(tvec);
3571     }
3572     len *= 2.5614f;
3573
3574     if (len != 0.0f) {
3575       /* only for fcurves */
3576       bool leftviolate = false, rightviolate = false;
3577
3578       if (!is_fcurve || fcurve_smoothing == FCURVE_SMOOTH_NONE) {
3579         if (len_a > 5.0f * len_b) {
3580           len_a = 5.0f * len_b;
3581         }
3582         if (len_b > 5.0f * len_a) {
3583           len_b = 5.0f * len_a;
3584         }
3585       }
3586
3587       if (ELEM(bezt->h1, HD_AUTO, HD_AUTO_ANIM)) {
3588         len_a /= len;
3589         madd_v3_v3v3fl(p2_h1, p2, tvec, -len_a);
3590
3591         if ((bezt->h1 == HD_AUTO_ANIM) && next && prev) { /* keep horizontal if extrema */
3592           float ydiff1 = prev->vec[1][1] - bezt->vec[1][1];
3593           float ydiff2 = next->vec[1][1] - bezt->vec[1][1];
3594           if ((ydiff1 <= 0.0f && ydiff2 <= 0.0f) || (ydiff1 >= 0.0f && ydiff2 >= 0.0f)) {
3595             bezt->vec[0][1] = bezt->vec[1][1];
3596             bezt->f5 = HD_AUTOTYPE_SPECIAL;
3597           }
3598           else { /* handles should not be beyond y coord of two others */
3599             if (ydiff1 <= 0.0f) {
3600               if (prev->vec[1][1] > bezt->vec[0][1]) {
3601                 bezt->vec[0][1] = prev->vec[1][1];
3602                 leftviolate = 1;
3603               }
3604             }
3605             else {
3606               if (prev->vec[1][1] < bezt->vec[0][1]) {
3607                 bezt->vec[0][1] = prev->vec[1][1];
3608                 leftviolate = 1;
3609               }
3610             }
3611           }
3612         }
3613       }
3614       if (ELEM(bezt->h2, HD_AUTO, HD_AUTO_ANIM)) {
3615         len_b /= len;
3616         madd_v3_v3v3fl(p2_h2, p2, tvec, len_b);
3617
3618         if ((bezt->h2 == HD_AUTO_ANIM) && next && prev) { /* keep horizontal if extrema */
3619           float ydiff1 = prev->vec[1][1] - bezt->vec[1][1];
3620           float ydiff2 = next->vec[1][1] - bezt->vec[1][1];
3621           if ((ydiff1 <= 0.0f && ydiff2 <= 0.0f) || (ydiff1 >= 0.0f && ydiff2 >= 0.0f)) {
3622             bezt->vec[2][1] = bezt->vec[1][1];
3623             bezt->f5 = HD_AUTOTYPE_SPECIAL;
3624           }
3625           else { /* handles should not be beyond y coord of two others */
3626             if (ydiff1 <= 0.0f) {
3627               if (next->vec[1][1] < bezt->vec[2][1]) {
3628                 bezt->vec[2][1] = next->vec[1][1];
3629                 rightviolate = 1;
3630               }
3631             }
3632             else {
3633               if (next->vec[1][1] > bezt->vec[2][1]) {
3634                 bezt->vec[2][1] = next->vec[1][1];
3635                 rightviolate = 1;
3636               }
3637             }
3638           }
3639         }
3640       }
3641       if (leftviolate || rightviolate) { /* align left handle */
3642         BLI_assert(is_fcurve);
3643         /* simple 2d calculation */
3644         float h1_x = p2_h1[0] - p2[0];
3645         float h2_x = p2[0] - p2_h2[0];
3646
3647         if (leftviolate) {
3648           p2_h2[1] = p2[1] + ((p2[1] - p2_h1[1]) / h1_x) * h2_x;
3649         }
3650         else {
3651           p2_h1[1] = p2[1] + ((p2[1] - p2_h2[1]) / h2_x) * h1_x;
3652         }
3653       }
3654     }
3655   }
3656
3657   if (bezt->h1 == HD_VECT) { /* vector */
3658     madd_v3_v3v3fl(p2_h1, p2, dvec_a, -1.0f / 3.0f);
3659   }
3660   if (bezt->h2 == HD_VECT) {
3661     madd_v3_v3v3fl(p2_h2, p2, dvec_b, 1.0f / 3.0f);
3662   }
3663
3664   if (skip_align ||
3665       /* when one handle is free, alignming makes no sense, see: T35952 */
3666       (ELEM(HD_FREE, bezt->h1, bezt->h2)) ||
3667       /* also when no handles are aligned, skip this step */
3668       (!ELEM(HD_ALIGN, bezt->h1, bezt->h2) && !ELEM(HD_ALIGN_DOUBLESIDE, bezt->h1, bezt->h2))) {
3669     /* handles need to be updated during animation and applying stuff like hooks,
3670      * but in such situations it's quite difficult to distinguish in which order
3671      * align handles should be aligned so skip them for now */
3672     return;
3673   }
3674
3675   len_a = len_v3v3(p2, p2_h1);
3676   len_b = len_v3v3(p2, p2_h2);
3677
3678   if (len_a == 0.0f) {
3679     len_a = 1.0f;
3680   }
3681   if (len_b == 0.0f) {
3682     len_b = 1.0f;
3683   }
3684
3685   len_ratio = len_a / len_b;
3686
3687   if (bezt->f1 & SELECT) {                               /* order of calculation */
3688     if (ELEM(bezt->h2, HD_ALIGN, HD_ALIGN_DOUBLESIDE)) { /* aligned */
3689       if (len_a > eps) {
3690         len = 1.0f / len_ratio;
3691         p2_h2[0] = p2[0] + len * (p2[0] - p2_h1[0]);
3692         p2_h2[1] = p2[1] + len * (p2[1] - p2_h1[1]);
3693         p2_h2[2] = p2[2] + len * (p2[2] - p2_h1[2]);
3694       }
3695     }
3696     if (ELEM(bezt->h1, HD_ALIGN, HD_ALIGN_DOUBLESIDE)) {
3697       if (len_b > eps) {
3698         len = len_ratio;
3699         p2_h1[0] = p2[0] + len * (p2[0] - p2_h2[0]);
3700         p2_h1[1] = p2[1] + len * (p2[1] - p2_h2[1]);
3701         p2_h1[2] = p2[2] + len * (p2[2] - p2_h2[2]);
3702       }
3703     }
3704   }
3705   else {
3706     if (ELEM(bezt->h1, HD_ALIGN, HD_ALIGN_DOUBLESIDE)) {
3707       if (len_b > eps) {
3708         len = len_ratio;
3709         p2_h1[0] = p2[0] + len * (p2[0] - p2_h2[0]);
3710         p2_h1[1] = p2[1] + len * (p2[1] - p2_h2[1]);
3711         p2_h1[2] = p2[2] + len * (p2[2] - p2_h2[2]);
3712       }
3713     }
3714     if (ELEM(bezt->h2, HD_ALIGN, HD_ALIGN_DOUBLESIDE)) { /* aligned */
3715       if (len_a > eps) {
3716         len = 1.0f / len_ratio;
3717         p2_h2[0] = p2[0] + len * (p2[0] - p2_h1[0]);
3718         p2_h2[1] = p2[1] + len * (p2[1] - p2_h1[1]);
3719         p2_h2[2] = p2[2] + len * (p2[2] - p2_h1[2]);
3720       }
3721     }
3722   }
3723
3724 #undef p2_h1
3725 #undef p2_h2
3726 }
3727
3728 static void calchandlesNurb_intern(Nurb *nu, bool skip_align)
3729 {
3730   BezTriple *bezt, *prev, *next;
3731   int a;
3732
3733   if (n