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