Cleanup: comment blocks
[blender.git] / source / blender / bmesh / operators / bmo_connect_pair.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  * Contributor(s): Campbell Barton.
19  *
20  * ***** END GPL LICENSE BLOCK *****
21  */
22
23 /** \file blender/bmesh/operators/bmo_connect_pair.c
24  *  \ingroup bmesh
25  *
26  * Connect vertex pair across multiple faces (splits faces).
27  */
28
29 #include "MEM_guardedalloc.h"
30
31 #include "BLI_math.h"
32 #include "BLI_utildefines.h"
33 #include "BLI_heap.h"
34
35 #include "bmesh.h"
36
37 #include "intern/bmesh_operators_private.h" /* own include */
38
39 #include "BLI_mempool.h"
40 #include "BLI_listbase.h"
41
42 /**
43  * Method for connecting across many faces.
44  *
45  * - use the line between both verts and their normal average to construct a matrix.
46  * - using the matrix, we can find all intersecting verts/edges.
47  * - walk the connected data and find the shortest path.
48  *   - store a heap of paths which are being scanned (#PathContext.states).
49  *   - continuously search the shortest path in the heap.
50  *   - never step over the same element twice (tag elements as #ELE_TOUCHED).
51  *     this avoids going into an eternal loop of there are many possible branches (see T45582).
52  *   - when running into a branch, create a new #PathLinkState state and add to the heap.
53  *   - when the target is reached, finish - since none of the other paths can be shorter then the one just found.
54  * - if the connection can't be found - fail.
55  * - with the connection found, split all edges tagging verts (or tag verts that sit on the intersection).
56  * - run the standard connect operator.
57  */
58
59 #define CONNECT_EPS 0.0001f
60 #define VERT_OUT 1
61 #define VERT_EXCLUDE 2
62
63 /* typically hidden faces */
64 #define FACE_EXCLUDE 2
65
66 /* any element we've walked over (only do it once!) */
67 #define ELE_TOUCHED 4
68
69 #define FACE_WALK_TEST(f)  (CHECK_TYPE_INLINE(f, BMFace *), \
70         BMO_face_flag_test(pc->bm_bmoflag, f, FACE_EXCLUDE) == 0)
71 #define VERT_WALK_TEST(v)  (CHECK_TYPE_INLINE(v, BMVert *), \
72         BMO_vert_flag_test(pc->bm_bmoflag, v, VERT_EXCLUDE) == 0)
73
74 #if 0
75 #define ELE_TOUCH_TEST(e) \
76         (CHECK_TYPE_ANY(e, BMVert *, BMEdge *, BMElem *, BMElemF *), \
77          BMO_elem_flag_test(pc->bm_bmoflag, (BMElemF *)e, ELE_TOUCHED))
78 #endif
79 #define ELE_TOUCH_MARK(e) \
80         { CHECK_TYPE_ANY(e, BMVert *, BMEdge *, BMElem *, BMElemF *); \
81           BMO_elem_flag_enable(pc->bm_bmoflag, (BMElemF *)e, ELE_TOUCHED); } ((void)0)
82
83
84 #define ELE_TOUCH_TEST_VERT(v) BMO_vert_flag_test(pc->bm_bmoflag, v, ELE_TOUCHED)
85 // #define ELE_TOUCH_MARK_VERT(v) BMO_vert_flag_enable(pc->bm_bmoflag, (BMElemF *)e, ELE_TOUCHED)
86
87 // #define ELE_TOUCH_TEST_EDGE(v) BMO_edge_flag_test(pc->bm_bmoflag, v, ELE_TOUCHED)
88 // #define ELE_TOUCH_MARK_EDGE(v) BMO_edge_flag_enable(pc->bm_bmoflag, (BMElemF *)e, ELE_TOUCHED)
89
90 // #define ELE_TOUCH_TEST_FACE(v) BMO_face_flag_test(pc->bm_bmoflag, v, ELE_TOUCHED)
91 // #define ELE_TOUCH_MARK_FACE(v) BMO_face_flag_enable(pc->bm_bmoflag, (BMElemF *)e, ELE_TOUCHED)
92
93 // #define DEBUG_PRINT
94
95 typedef struct PathContext {
96         Heap *states;
97         float matrix[3][3];
98         float axis_sep;
99
100         /* only to access BMO flags */
101         BMesh *bm_bmoflag;
102
103         BMVert *v_a, *v_b;
104
105         BLI_mempool *link_pool;
106 } PathContext;
107
108 /**
109  * Single linked list where each item contains state and points to previous path item.
110  */
111 typedef struct PathLink {
112         struct PathLink *next;
113         BMElem *ele;       /* edge or vert */
114         BMElem *ele_from;  /* edge or face we came from (not 'next->ele') */
115 } PathLink;
116
117 typedef struct PathLinkState {
118         /* chain of links */
119         struct PathLink *link_last;
120
121         /* length along links */
122         float dist;
123         float co_prev[3];
124 } PathLinkState;
125
126 /**
127   \name Min Dist Dir Util
128  *
129  * Simply getting the closest intersecting vert/edge is _not_ good enough. see T43792
130  * we need to get the closest in both directions since the absolute closest may be a dead-end.
131  *
132  * Logic is simple:
133  *
134  * - first intersection, store the direction.
135  * - successive intersections will update the first distance if its aligned with the first hit.
136  *   otherwise update the opposite distance.
137  * - caller stores best outcome in both directions.
138  *
139  * \{ */
140
141 typedef struct MinDistDir {
142         /* distance in both directions (FLT_MAX == uninitialized) */
143         float dist_min[2];
144         /* direction of the first intersection found */
145         float dir[3];
146 } MinDistDir;
147
148 #define MIN_DIST_DIR_INIT {{FLT_MAX, FLT_MAX}}
149
150 static int min_dist_dir_test(MinDistDir *mddir, const float dist_dir[3], const float dist_sq)
151 {
152
153         if (mddir->dist_min[0] == FLT_MAX) {
154                 return 0;
155         }
156         else {
157                 if (dot_v3v3(dist_dir, mddir->dir) > 0.0f) {
158                         if (dist_sq < mddir->dist_min[0]) {
159                                 return 0;
160                         }
161                 }
162                 else {
163                         if (dist_sq < mddir->dist_min[1]) {
164                                 return 1;
165                         }
166                 }
167         }
168
169         return -1;
170 }
171
172 static void min_dist_dir_update(MinDistDir *dist, const float dist_dir[3])
173 {
174         if (dist->dist_min[0] == FLT_MAX) {
175                 copy_v3_v3(dist->dir, dist_dir);
176         }
177 }
178
179 /** \} */
180
181
182 static int state_isect_co_pair(
183         const PathContext *pc,
184         const float co_a[3], const float co_b[3])
185 {
186         const float diff_a = dot_m3_v3_row_x((float (*)[3])pc->matrix, co_a) - pc->axis_sep;
187         const float diff_b = dot_m3_v3_row_x((float (*)[3])pc->matrix, co_b) - pc->axis_sep;
188
189         const int test_a = (fabsf(diff_a) < CONNECT_EPS) ? 0 : (diff_a < 0.0f) ? -1 : 1;
190         const int test_b = (fabsf(diff_b) < CONNECT_EPS) ? 0 : (diff_b < 0.0f) ? -1 : 1;
191
192         if ((test_a && test_b) && (test_a != test_b)) {
193                 return 1;  /* on either side */
194         }
195         else {
196                 return 0;
197         }
198 }
199
200 static int state_isect_co_exact(
201         const PathContext *pc,
202         const float co[3])
203 {
204         const float diff = dot_m3_v3_row_x((float (*)[3])pc->matrix, co) - pc->axis_sep;
205         return (fabsf(diff) <= CONNECT_EPS);
206 }
207
208 static float state_calc_co_pair_fac(
209         const PathContext *pc,
210         const float co_a[3], const float co_b[3])
211 {
212         float diff_a, diff_b, diff_tot;
213
214         diff_a = fabsf(dot_m3_v3_row_x((float (*)[3])pc->matrix, co_a) - pc->axis_sep);
215         diff_b = fabsf(dot_m3_v3_row_x((float (*)[3])pc->matrix, co_b) - pc->axis_sep);
216         diff_tot = (diff_a + diff_b);
217         return (diff_tot > FLT_EPSILON) ? (diff_a / diff_tot) : 0.5f;
218 }
219
220 static void state_calc_co_pair(
221         const PathContext *pc,
222         const float co_a[3], const float co_b[3],
223         float r_co[3])
224 {
225         const float fac = state_calc_co_pair_fac(pc, co_a, co_b);
226         interp_v3_v3v3(r_co, co_a, co_b, fac);
227 }
228
229 #ifndef NDEBUG
230 /**
231  * Ideally we wouldn't need this and for most cases we don't.
232  * But when a face has vertices that are on the boundary more than once this becomes tricky.
233  */
234 static bool state_link_find(const PathLinkState *state, BMElem *ele)
235 {
236         PathLink *link = state->link_last;
237         BLI_assert(ELEM(ele->head.htype, BM_VERT, BM_EDGE, BM_FACE));
238         if (link) {
239                 do {
240                         if (link->ele == ele) {
241                                 return true;
242                         }
243                 } while ((link = link->next));
244         }
245         return false;
246 }
247 #endif
248
249 static void state_link_add(
250         PathContext *pc, PathLinkState *state,
251         BMElem *ele, BMElem *ele_from)
252 {
253         PathLink *step_new = BLI_mempool_alloc(pc->link_pool);
254         BLI_assert(ele != ele_from);
255         BLI_assert(state_link_find(state, ele) == false);
256
257         /* never walk onto this again */
258         ELE_TOUCH_MARK(ele);
259
260 #ifdef DEBUG_PRINT
261         printf("%s: adding to state %p, %.4f - ", __func__, state, state->dist);
262         if (ele->head.htype == BM_VERT) {
263                 printf("vert %d, ", BM_elem_index_get(ele));
264         }
265         else if (ele->head.htype == BM_EDGE) {
266                 printf("edge %d, ", BM_elem_index_get(ele));
267         }
268         else {
269                 BLI_assert(0);
270         }
271
272         if (ele_from == NULL) {
273                 printf("from NULL\n");
274         }
275         else if (ele_from->head.htype == BM_EDGE) {
276                 printf("from edge %d\n", BM_elem_index_get(ele_from));
277         }
278         else if (ele_from->head.htype == BM_FACE) {
279                 printf("from face %d\n", BM_elem_index_get(ele_from));
280         }
281         else {
282                 BLI_assert(0);
283         }
284 #endif
285
286         /* track distance */
287         {
288                 float co[3];
289                 if (ele->head.htype == BM_VERT) {
290                         copy_v3_v3(co, ((BMVert *)ele)->co);
291                 }
292                 else if (ele->head.htype == BM_EDGE) {
293                         state_calc_co_pair(pc, ((BMEdge *)ele)->v1->co, ((BMEdge *)ele)->v2->co, co);
294                 }
295                 else {
296                         BLI_assert(0);
297                 }
298
299                 /* tally distance */
300                 if (ele_from) {
301                         state->dist += len_v3v3(state->co_prev, co);
302                 }
303                 copy_v3_v3(state->co_prev, co);
304         }
305
306         step_new->ele = ele;
307         step_new->ele_from = ele_from;
308         step_new->next = state->link_last;
309         state->link_last = step_new;
310 }
311
312 static PathLinkState *state_dupe_add(
313         PathLinkState *state, const PathLinkState *state_orig)
314 {
315         state = MEM_mallocN(sizeof(*state), __func__);
316         *state = *state_orig;
317         return state;
318 }
319
320 static PathLinkState *state_link_add_test(
321         PathContext *pc, PathLinkState *state, const PathLinkState *state_orig,
322         BMElem *ele, BMElem *ele_from)
323 {
324         const bool is_new = (state_orig->link_last != state->link_last);
325         if (is_new) {
326                 state = state_dupe_add(state, state_orig);
327         }
328
329         state_link_add(pc, state, ele, ele_from);
330
331         /* after adding a link so we use the updated 'state->dist' */
332         if (is_new) {
333                 BLI_heap_insert(pc->states, state->dist, state);
334         }
335
336         return state;
337 }
338
339 /* walk around the face edges */
340 static PathLinkState *state_step__face_edges(
341         PathContext *pc,
342         PathLinkState *state, const PathLinkState *state_orig,
343         BMLoop *l_iter, BMLoop *l_last,
344         MinDistDir *mddir)
345 {
346
347         BMLoop *l_iter_best[2] = {NULL, NULL};
348         int i;
349
350         do {
351                 if (state_isect_co_pair(pc, l_iter->v->co, l_iter->next->v->co)) {
352                         float dist_test;
353                         float co_isect[3];
354                         float dist_dir[3];
355                         int index;
356
357                         state_calc_co_pair(pc, l_iter->v->co, l_iter->next->v->co, co_isect);
358
359                         sub_v3_v3v3(dist_dir, co_isect, state_orig->co_prev);
360                         dist_test = len_squared_v3(dist_dir);
361                         if ((index = min_dist_dir_test(mddir, dist_dir, dist_test)) != -1) {
362                                 BMElem *ele_next      = (BMElem *)l_iter->e;
363                                 BMElem *ele_next_from = (BMElem *)l_iter->f;
364
365                                 if (FACE_WALK_TEST((BMFace *)ele_next_from) &&
366                                     (ELE_TOUCH_TEST_VERT((BMVert *)ele_next) == false))
367                                 {
368                                         min_dist_dir_update(mddir, dist_dir);
369                                         mddir->dist_min[index] = dist_test;
370                                         l_iter_best[index] = l_iter;
371                                 }
372                         }
373                 }
374         } while ((l_iter = l_iter->next) != l_last);
375
376         for (i = 0; i < 2; i++) {
377                 if ((l_iter = l_iter_best[i])) {
378                         BMElem *ele_next      = (BMElem *)l_iter->e;
379                         BMElem *ele_next_from = (BMElem *)l_iter->f;
380                         state = state_link_add_test(pc, state, state_orig, ele_next, ele_next_from);
381                 }
382         }
383
384         return state;
385 }
386
387 /* walk around the face verts */
388 static PathLinkState *state_step__face_verts(
389         PathContext *pc,
390         PathLinkState *state, const PathLinkState *state_orig,
391         BMLoop *l_iter, BMLoop *l_last,
392         MinDistDir *mddir)
393 {
394         BMLoop *l_iter_best[2] = {NULL, NULL};
395         int i;
396
397         do {
398                 if (state_isect_co_exact(pc, l_iter->v->co)) {
399                         float dist_test;
400                         const float *co_isect = l_iter->v->co;
401                         float dist_dir[3];
402                         int index;
403
404                         sub_v3_v3v3(dist_dir, co_isect, state_orig->co_prev);
405                         dist_test = len_squared_v3(dist_dir);
406                         if ((index = min_dist_dir_test(mddir, dist_dir, dist_test)) != -1) {
407                                 BMElem *ele_next      = (BMElem *)l_iter->v;
408                                 BMElem *ele_next_from = (BMElem *)l_iter->f;
409
410                                 if (FACE_WALK_TEST((BMFace *)ele_next_from) &&
411                                     (ELE_TOUCH_TEST_VERT((BMVert *)ele_next) == false))
412                                 {
413                                         min_dist_dir_update(mddir, dist_dir);
414                                         mddir->dist_min[index] = dist_test;
415                                         l_iter_best[index] = l_iter;
416                                 }
417                         }
418                 }
419         } while ((l_iter = l_iter->next) != l_last);
420
421         for (i = 0; i < 2; i++) {
422                 if ((l_iter = l_iter_best[i])) {
423                         BMElem *ele_next      = (BMElem *)l_iter->v;
424                         BMElem *ele_next_from = (BMElem *)l_iter->f;
425                         state = state_link_add_test(pc, state, state_orig, ele_next, ele_next_from);
426                 }
427         }
428
429         return state;
430 }
431
432 static bool state_step(PathContext *pc, PathLinkState *state)
433 {
434         PathLinkState state_orig = *state;
435         BMElem *ele = state->link_last->ele;
436         const void *ele_from = state->link_last->ele_from;
437
438         if (ele->head.htype == BM_EDGE) {
439                 BMEdge *e = (BMEdge *)ele;
440
441                 BMIter liter;
442                 BMLoop *l_start;
443
444                 BM_ITER_ELEM (l_start, &liter, e, BM_LOOPS_OF_EDGE) {
445                         if ((l_start->f != ele_from) &&
446                             FACE_WALK_TEST(l_start->f))
447                         {
448                                 MinDistDir mddir = MIN_DIST_DIR_INIT;
449                                 /* very similar to block below */
450                                 state = state_step__face_edges(pc, state, &state_orig,
451                                                                l_start->next, l_start, &mddir);
452                                 state = state_step__face_verts(pc, state, &state_orig,
453                                                                l_start->next->next, l_start, &mddir);
454                         }
455                 }
456         }
457         else if (ele->head.htype == BM_VERT) {
458                 BMVert *v = (BMVert *)ele;
459
460                 /* vert loops */
461                 {
462                         BMIter liter;
463                         BMLoop *l_start;
464
465                         BM_ITER_ELEM (l_start, &liter, v, BM_LOOPS_OF_VERT) {
466                                 if ((l_start->f != ele_from) &&
467                                     FACE_WALK_TEST(l_start->f))
468                                 {
469                                         MinDistDir mddir = MIN_DIST_DIR_INIT;
470                                         /* very similar to block above */
471                                         state = state_step__face_edges(pc, state, &state_orig,
472                                                                        l_start->next, l_start->prev, &mddir);
473                                         if (l_start->f->len > 3) {
474                                                 /* adjacent verts are handled in state_step__vert_edges */
475                                                 state = state_step__face_verts(pc, state, &state_orig,
476                                                                                l_start->next->next, l_start->prev, &mddir);
477                                         }
478                                 }
479                         }
480                 }
481
482                 /* vert edges  */
483                 {
484                         BMIter eiter;
485                         BMEdge *e;
486                         BM_ITER_ELEM (e, &eiter, v, BM_EDGES_OF_VERT) {
487                                 BMVert *v_other = BM_edge_other_vert(e, v);
488                                 if (((BMElem *)e != ele_from) &&
489                                     VERT_WALK_TEST(v_other))
490                                 {
491                                         if (state_isect_co_exact(pc, v_other->co)) {
492                                                 BMElem *ele_next      = (BMElem *)v_other;
493                                                 BMElem *ele_next_from = (BMElem *)e;
494                                                 if (ELE_TOUCH_TEST_VERT((BMVert *)ele_next) == false) {
495                                                         state = state_link_add_test(pc, state, &state_orig, ele_next, ele_next_from);
496                                                 }
497                                         }
498                                 }
499                         }
500                 }
501         }
502         else {
503                 BLI_assert(0);
504         }
505         return (state_orig.link_last != state->link_last);
506 }
507
508 /**
509  * Get a orientation matrix from 2 vertices.
510  */
511 static void bm_vert_pair_to_matrix(BMVert *v_pair[2], float r_unit_mat[3][3])
512 {
513         const float eps = 1e-8f;
514
515         float basis_dir[3];
516         float basis_tmp[3];
517         float basis_nor[3];
518
519         sub_v3_v3v3(basis_dir, v_pair[0]->co, v_pair[1]->co);
520         normalize_v3(basis_dir);
521
522 #if 0
523         add_v3_v3v3(basis_nor, v_pair[0]->no, v_pair[1]->no);
524         cross_v3_v3v3(basis_tmp, basis_nor, basis_dir);
525         cross_v3_v3v3(basis_nor, basis_tmp, basis_dir);
526 #else
527         /* align both normals to the directions before combining */
528         {
529                 float basis_nor_a[3];
530                 float basis_nor_b[3];
531
532                 /* align normal to direction */
533                 project_plane_v3_v3v3(basis_nor_a, v_pair[0]->no, basis_dir);
534                 project_plane_v3_v3v3(basis_nor_b, v_pair[1]->no, basis_dir);
535
536                 /* don't normalize before combining so as normals approach the direction, they have less effect (T46784). */
537
538                 /* combine the normals */
539                 /* for flipped faces */
540                 if (dot_v3v3(basis_nor_a, basis_nor_b) < 0.0f) {
541                         negate_v3(basis_nor_b);
542                 }
543                 add_v3_v3v3(basis_nor, basis_nor_a, basis_nor_b);
544         }
545 #endif
546
547         /* get third axis */
548         normalize_v3(basis_nor);
549         cross_v3_v3v3(basis_tmp, basis_dir, basis_nor);
550
551
552         /* Try get the axis from surrounding faces, fallback to 'ortho_v3_v3' */
553         if (UNLIKELY(normalize_v3(basis_tmp) < eps)) {
554                 /* vertex normals are directly opposite */
555
556                 /* find the loop with the lowest angle */
557                 struct { float nor[3]; float angle_cos; } axis_pair[2];
558                 int i;
559
560                 for (i = 0; i < 2; i++) {
561                         BMIter liter;
562                         BMLoop *l;
563
564                         zero_v2(axis_pair[i].nor);
565                         axis_pair[i].angle_cos = -FLT_MAX;
566
567                         BM_ITER_ELEM (l, &liter, v_pair[i], BM_LOOPS_OF_VERT) {
568                                 float basis_dir_proj[3];
569                                 float angle_cos_test;
570
571                                 /* project basis dir onto the normal to find its closest angle */
572                                 project_plane_v3_v3v3(basis_dir_proj, basis_dir, l->f->no);
573
574                                 if (normalize_v3(basis_dir_proj) > eps) {
575                                         angle_cos_test = dot_v3v3(basis_dir_proj, basis_dir);
576
577                                         if (angle_cos_test > axis_pair[i].angle_cos) {
578                                                 axis_pair[i].angle_cos = angle_cos_test;
579                                                 copy_v3_v3(axis_pair[i].nor, basis_dir_proj);
580                                         }
581                                 }
582                         }
583                 }
584
585                 /* create a new 'basis_nor' from the best direction.
586                  * note: we could add the directions,
587                  * but this more often gives 45d rotated matrix, so just use the best one. */
588                 copy_v3_v3(basis_nor, axis_pair[axis_pair[0].angle_cos < axis_pair[1].angle_cos].nor);
589                 project_plane_v3_v3v3(basis_nor, basis_nor, basis_dir);
590
591                 cross_v3_v3v3(basis_tmp, basis_dir, basis_nor);
592
593                 /* last resort, pick _any_ ortho axis */
594                 if (UNLIKELY(normalize_v3(basis_tmp) < eps)) {
595                         ortho_v3_v3(basis_nor, basis_dir);
596                         normalize_v3(basis_nor);
597                         cross_v3_v3v3(basis_tmp, basis_dir, basis_nor);
598                         normalize_v3(basis_tmp);
599                 }
600         }
601
602         copy_v3_v3(r_unit_mat[0], basis_tmp);
603         copy_v3_v3(r_unit_mat[1], basis_dir);
604         copy_v3_v3(r_unit_mat[2], basis_nor);
605         if (invert_m3(r_unit_mat) == false) {
606                 unit_m3(r_unit_mat);
607         }
608 }
609
610 void bmo_connect_vert_pair_exec(BMesh *bm, BMOperator *op)
611 {
612         BMOpSlot *op_verts_slot = BMO_slot_get(op->slots_in, "verts");
613
614         PathContext pc;
615         PathLinkState state_best = {NULL};
616
617         if (op_verts_slot->len != 2) {
618                 /* fail! */
619                 return;
620         }
621
622         pc.bm_bmoflag = bm;
623         pc.v_a = ((BMVert **)op_verts_slot->data.p)[0];
624         pc.v_b = ((BMVert **)op_verts_slot->data.p)[1];
625
626         /* fail! */
627         if (!(pc.v_a && pc.v_b)) {
628                 return;
629         }
630
631 #ifdef DEBUG_PRINT
632         printf("%s: v_a: %d\n", __func__, BM_elem_index_get(pc.v_a));
633         printf("%s: v_b: %d\n", __func__, BM_elem_index_get(pc.v_b));
634 #endif
635
636         /* tag so we won't touch ever (typically hidden faces) */
637         BMO_slot_buffer_flag_enable(bm, op->slots_in, "faces_exclude", BM_FACE, FACE_EXCLUDE);
638         BMO_slot_buffer_flag_enable(bm, op->slots_in, "verts_exclude", BM_VERT, VERT_EXCLUDE);
639
640         /* setup context */
641         {
642                 pc.states = BLI_heap_new();
643                 pc.link_pool = BLI_mempool_create(sizeof(PathLink), 0, 512, BLI_MEMPOOL_NOP);
644         }
645
646         /* calculate matrix */
647         {
648                 bm_vert_pair_to_matrix(&pc.v_a, pc.matrix);
649                 pc.axis_sep = dot_m3_v3_row_x(pc.matrix, pc.v_a->co);
650         }
651
652         /* add first vertex */
653         {
654                 PathLinkState *state;
655                 state = MEM_callocN(sizeof(*state), __func__);
656                 state_link_add(&pc, state, (BMElem *)pc.v_a, NULL);
657                 BLI_heap_insert(pc.states, state->dist, state);
658         }
659
660
661         while (!BLI_heap_is_empty(pc.states)) {
662
663 #ifdef DEBUG_PRINT
664                 printf("\n%s: stepping %d\n", __func__, BLI_heap_size(pc.states));
665 #endif
666
667                 while (!BLI_heap_is_empty(pc.states)) {
668                         PathLinkState *state = BLI_heap_popmin(pc.states);
669
670                         /* either we insert this into 'pc.states' or its freed */
671                         bool continue_search;
672
673                         if (state->link_last->ele == (BMElem *)pc.v_b) {
674                                 /* pass, wait until all are found */
675 #ifdef DEBUG_PRINT
676                                 printf("%s: state %p loop found %.4f\n", __func__, state, state->dist);
677 #endif
678                                 state_best = *state;
679
680                                 /* we're done, exit all loops */
681                                 BLI_heap_clear(pc.states, MEM_freeN);
682                                 continue_search = false;
683                         }
684                         else if (state_step(&pc, state)) {
685                                 continue_search = true;
686                         }
687                         else {
688                                 /* didn't reach the end, remove it,
689                                  * links are shared between states so just free the link_pool at the end */
690
691 #ifdef DEBUG_PRINT
692                                 printf("%s: state %p removed\n", __func__, state);
693 #endif
694                                 continue_search = false;
695                         }
696
697                         if (continue_search) {
698                                 BLI_heap_insert(pc.states, state->dist, state);
699                         }
700                         else {
701                                 MEM_freeN(state);
702                         }
703                 }
704         }
705
706         if (state_best.link_last) {
707                 PathLink *link;
708
709                 /* find the best state */
710                 link = state_best.link_last;
711                 do {
712                         if (link->ele->head.htype == BM_EDGE) {
713                                 BMEdge *e = (BMEdge *)link->ele;
714                                 BMVert *v_new;
715                                 float e_fac = state_calc_co_pair_fac(&pc, e->v1->co, e->v2->co);
716                                 v_new = BM_edge_split(bm, e, e->v1, NULL, e_fac);
717                                 BMO_vert_flag_enable(bm, v_new, VERT_OUT);
718                         }
719                         else if (link->ele->head.htype == BM_VERT) {
720                                 BMVert *v = (BMVert *)link->ele;
721                                 BMO_vert_flag_enable(bm, v, VERT_OUT);
722                         }
723                         else {
724                                 BLI_assert(0);
725                         }
726                 } while ((link = link->next));
727         }
728
729         BMO_vert_flag_enable(bm, pc.v_a, VERT_OUT);
730         BMO_vert_flag_enable(bm, pc.v_b, VERT_OUT);
731
732         BLI_mempool_destroy(pc.link_pool);
733
734         BLI_heap_free(pc.states, MEM_freeN);
735
736 #if 1
737         if (state_best.link_last) {
738                 BMOperator op_sub;
739                 BMO_op_initf(bm, &op_sub, 0,
740                              "connect_verts verts=%fv faces_exclude=%s check_degenerate=%b",
741                              VERT_OUT, op, "faces_exclude", true);
742                 BMO_op_exec(bm, &op_sub);
743                 BMO_slot_copy(&op_sub, slots_out, "edges.out",
744                               op,      slots_out, "edges.out");
745                 BMO_op_finish(bm, &op_sub);
746         }
747 #endif
748 }