Cleanup: style, use braces in source/ (include disabled blocks)
[blender.git] / source / blender / freestyle / intern / view_map / ViewMapBuilder.cpp
1 /*
2  * This program is free software; you can redistribute it and/or
3  * modify it under the terms of the GNU General Public License
4  * as published by the Free Software Foundation; either version 2
5  * of the License, or (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program; if not, write to the Free Software Foundation,
14  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
15  */
16
17 /** \file
18  * \ingroup freestyle
19  * \brief Class to build silhouette edges from a Winged-Edge structure
20  */
21
22 #include <algorithm>
23 #include <memory>
24 #include <stdexcept>
25 #include <sstream>
26
27 #include "FRS_freestyle.h"
28
29 #include "BoxGrid.h"
30 #include "CulledOccluderSource.h"
31 #include "HeuristicGridDensityProviderFactory.h"
32 #include "OccluderSource.h"
33 #include "SphericalGrid.h"
34 #include "ViewMapBuilder.h"
35
36 #include "../geometry/GridHelpers.h"
37 #include "../geometry/GeomUtils.h"
38
39 #include "../winged_edge/WFillGrid.h"
40
41 #include "BKE_global.h"
42
43 namespace Freestyle {
44
45 // XXX Grmll... G is used as template's typename parameter :/
46 static const Global &_global = G;
47
48 #define LOGGING 0
49
50 using namespace std;
51
52 template<typename G, typename I>
53 static void findOccludee(FEdge *fe,
54                          G & /*grid*/,
55                          I &occluders,
56                          real epsilon,
57                          WFace **oaWFace,
58                          Vec3r &u,
59                          Vec3r &A,
60                          Vec3r &origin,
61                          Vec3r &edgeDir,
62                          vector<WVertex *> &faceVertices)
63 {
64   WFace *face = NULL;
65   if (fe->isSmooth()) {
66     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
67     face = (WFace *)fes->face();
68   }
69   WFace *oface;
70   bool skipFace;
71
72   WVertex::incoming_edge_iterator ie;
73
74   *oaWFace = NULL;
75   if (((fe)->getNature() & Nature::SILHOUETTE) || ((fe)->getNature() & Nature::BORDER)) {
76     // we cast a ray from A in the same direction but looking behind
77     Vec3r v(-u[0], -u[1], -u[2]);
78     bool noIntersection = true;
79     real mint = FLT_MAX;
80
81     for (occluders.initAfterTarget(); occluders.validAfterTarget(); occluders.nextOccludee()) {
82 #if LOGGING
83       if (_global.debug & G_DEBUG_FREESTYLE) {
84         cout << "\t\tEvaluating intersection for occludee " << occluders.getWFace() << " and ray "
85              << A << " * " << u << endl;
86       }
87 #endif
88       oface = occluders.getWFace();
89       Polygon3r *p = occluders.getCameraSpacePolygon();
90       real d = -((p->getVertices())[0] * p->getNormal());
91       real t, t_u, t_v;
92
93       if (0 != face) {
94         skipFace = false;
95
96         if (face == oface) {
97           continue;
98         }
99
100         if (faceVertices.empty()) {
101           continue;
102         }
103
104         for (vector<WVertex *>::iterator fv = faceVertices.begin(), fvend = faceVertices.end();
105              fv != fvend;
106              ++fv) {
107           if ((*fv)->isBoundary()) {
108             continue;
109           }
110           WVertex::incoming_edge_iterator iebegin = (*fv)->incoming_edges_begin();
111           WVertex::incoming_edge_iterator ieend = (*fv)->incoming_edges_end();
112           for (ie = iebegin; ie != ieend; ++ie) {
113             if ((*ie) == 0) {
114               continue;
115             }
116
117             WFace *sface = (*ie)->GetbFace();
118             if (sface == oface) {
119               skipFace = true;
120               break;
121             }
122           }
123           if (skipFace) {
124             break;
125           }
126         }
127         if (skipFace) {
128           continue;
129         }
130       }
131       else {
132         // check whether the edge and the polygon plane are coincident:
133         //-------------------------------------------------------------
134         // first let us compute the plane equation.
135         if (GeomUtils::COINCIDENT ==
136             GeomUtils::intersectRayPlane(origin, edgeDir, p->getNormal(), d, t, epsilon)) {
137 #if LOGGING
138           if (_global.debug & G_DEBUG_FREESTYLE) {
139             cout << "\t\tRejecting occluder for target coincidence." << endl;
140           }
141 #endif
142           continue;
143         }
144       }
145
146       if (p->rayIntersect(A, v, t, t_u, t_v)) {
147 #if LOGGING
148         if (_global.debug & G_DEBUG_FREESTYLE) {
149           cout << "\t\tRay " << A << " * " << v << " intersects at time " << t << endl;
150           cout << "\t\t(v * normal) == " << (v * p->getNormal()) << " for normal "
151                << p->getNormal() << endl;
152         }
153 #endif
154         if (fabs(v * p->getNormal()) > 0.0001) {
155           if ((t > 0.0) /* && (t<1.0) */) {
156             if (t < mint) {
157               *oaWFace = oface;
158               mint = t;
159               noIntersection = false;
160               fe->setOccludeeIntersection(Vec3r(A + t * v));
161 #if LOGGING
162               if (_global.debug & G_DEBUG_FREESTYLE) {
163                 cout << "\t\tIs occludee" << endl;
164               }
165 #endif
166             }
167           }
168         }
169
170         occluders.reportDepth(A, v, t);
171       }
172     }
173
174     if (noIntersection) {
175       *oaWFace = NULL;
176     }
177   }
178 }
179
180 template<typename G, typename I>
181 static void findOccludee(FEdge *fe, G &grid, real epsilon, ViewEdge * /*ve*/, WFace **oaFace)
182 {
183   Vec3r A;
184   Vec3r edgeDir;
185   Vec3r origin;
186   A = Vec3r(((fe)->vertexA()->point3D() + (fe)->vertexB()->point3D()) / 2.0);
187   edgeDir = Vec3r((fe)->vertexB()->point3D() - (fe)->vertexA()->point3D());
188   edgeDir.normalize();
189   origin = Vec3r((fe)->vertexA()->point3D());
190   Vec3r u;
191   if (grid.orthographicProjection()) {
192     u = Vec3r(0.0, 0.0, grid.viewpoint().z() - A.z());
193   }
194   else {
195     u = Vec3r(grid.viewpoint() - A);
196   }
197   u.normalize();
198
199   vector<WVertex *> faceVertices;
200
201   WFace *face = NULL;
202   if (fe->isSmooth()) {
203     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
204     face = (WFace *)fes->face();
205   }
206
207   if (face) {
208     face->RetrieveVertexList(faceVertices);
209   }
210
211   I occluders(grid, A, epsilon);
212   findOccludee<G, I>(fe, grid, occluders, epsilon, oaFace, u, A, origin, edgeDir, faceVertices);
213 }
214
215 // computeVisibility takes a pointer to foundOccluders, instead of using a reference,
216 // so that computeVeryFastVisibility can skip the AddOccluders step with minimal overhead.
217 template<typename G, typename I>
218 static int computeVisibility(ViewMap *viewMap,
219                              FEdge *fe,
220                              G &grid,
221                              real epsilon,
222                              ViewEdge * /*ve*/,
223                              WFace **oaWFace,
224                              set<ViewShape *> *foundOccluders)
225 {
226   int qi = 0;
227
228   Vec3r center;
229   Vec3r edgeDir;
230   Vec3r origin;
231
232   center = fe->center3d();
233   edgeDir = Vec3r(fe->vertexB()->point3D() - fe->vertexA()->point3D());
234   edgeDir.normalize();
235   origin = Vec3r(fe->vertexA()->point3D());
236
237   Vec3r vp;
238   if (grid.orthographicProjection()) {
239     vp = Vec3r(center.x(), center.y(), grid.viewpoint().z());
240   }
241   else {
242     vp = Vec3r(grid.viewpoint());
243   }
244   Vec3r u(vp - center);
245   real raylength = u.norm();
246   u.normalize();
247
248   WFace *face = NULL;
249   if (fe->isSmooth()) {
250     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
251     face = (WFace *)fes->face();
252   }
253   vector<WVertex *> faceVertices;
254   WVertex::incoming_edge_iterator ie;
255
256   WFace *oface;
257   bool skipFace;
258
259   if (face) {
260     face->RetrieveVertexList(faceVertices);
261   }
262
263   I occluders(grid, center, epsilon);
264
265   for (occluders.initBeforeTarget(); occluders.validBeforeTarget(); occluders.nextOccluder()) {
266     // If we're dealing with an exact silhouette, check whether we must take care of this occluder
267     // of not. (Indeed, we don't consider the occluders that share at least one vertex with the
268     // face containing this edge).
269     //-----------
270     oface = occluders.getWFace();
271     Polygon3r *p = occluders.getCameraSpacePolygon();
272     real t, t_u, t_v;
273 #if LOGGING
274     if (_global.debug & G_DEBUG_FREESTYLE) {
275       cout << "\t\tEvaluating intersection for occluder " << (p->getVertices())[0]
276            << (p->getVertices())[1] << (p->getVertices())[2] << endl
277            << "\t\t\tand ray " << vp << " * " << u << " (center " << center << ")" << endl;
278     }
279 #endif
280
281 #if LOGGING
282     Vec3r v(vp - center);
283     real rl = v.norm();
284     v.normalize();
285     vector<Vec3r> points;
286     // Iterate over vertices, storing projections in points
287     for (vector<WOEdge *>::const_iterator woe = oface->getEdgeList().begin(),
288                                           woend = oface->getEdgeList().end();
289          woe != woend;
290          woe++) {
291       points.push_back(Vec3r((*woe)->GetaVertex()->GetVertex()));
292     }
293     Polygon3r p1(points, oface->GetNormal());
294     Vec3r v1((p1.getVertices())[0]);
295     real d = -(v1 * p->getNormal());
296     if (_global.debug & G_DEBUG_FREESTYLE) {
297       cout << "\t\tp:  " << (p->getVertices())[0] << (p->getVertices())[1] << (p->getVertices())[2]
298            << ", norm: " << p->getNormal() << endl;
299       cout << "\t\tp1: " << (p1.getVertices())[0] << (p1.getVertices())[1] << (p1.getVertices())[2]
300            << ", norm: " << p1.getNormal() << endl;
301     }
302 #else
303     real d = -((p->getVertices())[0] * p->getNormal());
304 #endif
305
306     if (face) {
307 #if LOGGING
308       if (_global.debug & G_DEBUG_FREESTYLE) {
309         cout << "\t\tDetermining face adjacency...";
310       }
311 #endif
312       skipFace = false;
313
314       if (face == oface) {
315 #if LOGGING
316         if (_global.debug & G_DEBUG_FREESTYLE) {
317           cout << "  Rejecting occluder for face concurrency." << endl;
318         }
319 #endif
320         continue;
321       }
322
323       for (vector<WVertex *>::iterator fv = faceVertices.begin(), fvend = faceVertices.end();
324            fv != fvend;
325            ++fv) {
326         if ((*fv)->isBoundary()) {
327           continue;
328         }
329
330         WVertex::incoming_edge_iterator iebegin = (*fv)->incoming_edges_begin();
331         WVertex::incoming_edge_iterator ieend = (*fv)->incoming_edges_end();
332         for (ie = iebegin; ie != ieend; ++ie) {
333           if ((*ie) == 0) {
334             continue;
335           }
336
337           WFace *sface = (*ie)->GetbFace();
338           // WFace *sfacea = (*ie)->GetaFace();
339           // if ((sface == oface) || (sfacea == oface))
340           if (sface == oface) {
341             skipFace = true;
342             break;
343           }
344         }
345         if (skipFace) {
346           break;
347         }
348       }
349       if (skipFace) {
350 #if LOGGING
351         if (_global.debug & G_DEBUG_FREESTYLE) {
352           cout << "  Rejecting occluder for face adjacency." << endl;
353         }
354 #endif
355         continue;
356       }
357     }
358     else {
359       // check whether the edge and the polygon plane are coincident:
360       //-------------------------------------------------------------
361       // first let us compute the plane equation.
362       if (GeomUtils::COINCIDENT ==
363           GeomUtils::intersectRayPlane(origin, edgeDir, p->getNormal(), d, t, epsilon)) {
364 #if LOGGING
365         if (_global.debug & G_DEBUG_FREESTYLE) {
366           cout << "\t\tRejecting occluder for target coincidence." << endl;
367         }
368 #endif
369         continue;
370       }
371     }
372
373 #if LOGGING
374     if (_global.debug & G_DEBUG_FREESTYLE) {
375       real x;
376       if (p1.rayIntersect(center, v, x, t_u, t_v)) {
377         cout << "\t\tRay should intersect at time " << (rl - x) << ". Center: " << center
378              << ", V: " << v << ", RL: " << rl << ", T:" << x << endl;
379       }
380       else {
381         cout << "\t\tRay should not intersect. Center: " << center << ", V: " << v
382              << ", RL: " << rl << endl;
383       }
384     }
385 #endif
386
387     if (p->rayIntersect(center, u, t, t_u, t_v)) {
388 #if LOGGING
389       if (_global.debug & G_DEBUG_FREESTYLE) {
390         cout << "\t\tRay " << center << " * " << u << " intersects at time " << t
391              << " (raylength is " << raylength << ")" << endl;
392         cout << "\t\t(u * normal) == " << (u * p->getNormal()) << " for normal " << p->getNormal()
393              << endl;
394       }
395 #endif
396       if (fabs(u * p->getNormal()) > 0.0001) {
397         if ((t > 0.0) && (t < raylength)) {
398 #if LOGGING
399           if (_global.debug & G_DEBUG_FREESTYLE) {
400             cout << "\t\tIs occluder" << endl;
401           }
402 #endif
403           if (foundOccluders != NULL) {
404             ViewShape *vshape = viewMap->viewShape(oface->GetVertex(0)->shape()->GetId());
405             foundOccluders->insert(vshape);
406           }
407           ++qi;
408
409           if (!grid.enableQI()) {
410             break;
411           }
412         }
413
414         occluders.reportDepth(center, u, t);
415       }
416     }
417   }
418
419   // Find occludee
420   findOccludee<G, I>(
421       fe, grid, occluders, epsilon, oaWFace, u, center, origin, edgeDir, faceVertices);
422
423   return qi;
424 }
425
426 // computeCumulativeVisibility returns the lowest x such that the majority of FEdges have QI <= x
427 //
428 // This was probably the original intention of the "normal" algorithm on which
429 // computeDetailedVisibility is based. But because the "normal" algorithm chooses the most popular
430 // QI, without considering any other values, a ViewEdge with FEdges having QIs of 0, 21, 22, 23, 24
431 // and 25 will end up having a total QI of 0, even though most of the FEdges are heavily occluded.
432 // computeCumulativeVisibility will treat this case as a QI of 22 because 3 out of 6 occluders have
433 // QI <= 22.
434
435 template<typename G, typename I>
436 static void computeCumulativeVisibility(ViewMap *ioViewMap,
437                                         G &grid,
438                                         real epsilon,
439                                         RenderMonitor *iRenderMonitor)
440 {
441   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
442
443   FEdge *fe, *festart;
444   int nSamples = 0;
445   vector<WFace *> wFaces;
446   WFace *wFace = NULL;
447   unsigned cnt = 0;
448   unsigned cntStep = (unsigned)ceil(0.01f * vedges.size());
449   unsigned tmpQI = 0;
450   unsigned qiClasses[256];
451   unsigned maxIndex, maxCard;
452   unsigned qiMajority;
453   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
454     if (iRenderMonitor) {
455       if (iRenderMonitor->testBreak()) {
456         break;
457       }
458       if (cnt % cntStep == 0) {
459         stringstream ss;
460         ss << "Freestyle: Visibility computations " << (100 * cnt / vedges.size()) << "%";
461         iRenderMonitor->setInfo(ss.str());
462         iRenderMonitor->progress((float)cnt / vedges.size());
463       }
464       cnt++;
465     }
466 #if LOGGING
467     if (_global.debug & G_DEBUG_FREESTYLE) {
468       cout << "Processing ViewEdge " << (*ve)->getId() << endl;
469     }
470 #endif
471     // Find an edge to test
472     if (!(*ve)->isInImage()) {
473       // This view edge has been proscenium culled
474       (*ve)->setQI(255);
475       (*ve)->setaShape(0);
476 #if LOGGING
477       if (_global.debug & G_DEBUG_FREESTYLE) {
478         cout << "\tCulled." << endl;
479       }
480 #endif
481       continue;
482     }
483
484     // Test edge
485     festart = (*ve)->fedgeA();
486     fe = (*ve)->fedgeA();
487     qiMajority = 0;
488     do {
489       if (fe != NULL && fe->isInImage()) {
490         qiMajority++;
491       }
492       fe = fe->nextEdge();
493     } while (fe && fe != festart);
494
495     if (qiMajority == 0) {
496       // There are no occludable FEdges on this ViewEdge
497       // This should be impossible.
498       if (_global.debug & G_DEBUG_FREESTYLE) {
499         cout << "View Edge in viewport without occludable FEdges: " << (*ve)->getId() << endl;
500       }
501       // We can recover from this error:
502       // Treat this edge as fully visible with no occludee
503       (*ve)->setQI(0);
504       (*ve)->setaShape(0);
505       continue;
506     }
507     else {
508       ++qiMajority;
509       qiMajority >>= 1;
510     }
511 #if LOGGING
512     if (_global.debug & G_DEBUG_FREESTYLE) {
513       cout << "\tqiMajority: " << qiMajority << endl;
514     }
515 #endif
516
517     tmpQI = 0;
518     maxIndex = 0;
519     maxCard = 0;
520     nSamples = 0;
521     memset(qiClasses, 0, 256 * sizeof(*qiClasses));
522     set<ViewShape *> foundOccluders;
523
524     fe = (*ve)->fedgeA();
525     do {
526       if (!fe || !fe->isInImage()) {
527         fe = fe->nextEdge();
528         continue;
529       }
530       if ((maxCard < qiMajority)) {
531         // ARB: change &wFace to wFace and use reference in called function
532         tmpQI = computeVisibility<G, I>(
533             ioViewMap, fe, grid, epsilon, *ve, &wFace, &foundOccluders);
534 #if LOGGING
535         if (_global.debug & G_DEBUG_FREESTYLE) {
536           cout << "\tFEdge: visibility " << tmpQI << endl;
537         }
538 #endif
539
540         // ARB: This is an error condition, not an alert condition.
541         // Some sort of recovery or abort is necessary.
542         if (tmpQI >= 256) {
543           cerr << "Warning: too many occluding levels" << endl;
544           // ARB: Wild guess: instead of aborting or corrupting memory, treat as tmpQI == 255
545           tmpQI = 255;
546         }
547
548         if (++qiClasses[tmpQI] > maxCard) {
549           maxCard = qiClasses[tmpQI];
550           maxIndex = tmpQI;
551         }
552       }
553       else {
554         // ARB: FindOccludee is redundant if ComputeRayCastingVisibility has been called
555         // ARB: change &wFace to wFace and use reference in called function
556         findOccludee<G, I>(fe, grid, epsilon, *ve, &wFace);
557 #if LOGGING
558         if (_global.debug & G_DEBUG_FREESTYLE) {
559           cout << "\tFEdge: occludee only (" << (wFace != NULL ? "found" : "not found") << ")"
560                << endl;
561         }
562 #endif
563       }
564
565       // Store test results
566       if (wFace) {
567         vector<Vec3r> vertices;
568         for (int i = 0, numEdges = wFace->numberOfEdges(); i < numEdges; ++i) {
569           vertices.push_back(Vec3r(wFace->GetVertex(i)->GetVertex()));
570         }
571         Polygon3r poly(vertices, wFace->GetNormal());
572         poly.userdata = (void *)wFace;
573         fe->setaFace(poly);
574         wFaces.push_back(wFace);
575         fe->setOccludeeEmpty(false);
576 #if LOGGING
577         if (_global.debug & G_DEBUG_FREESTYLE) {
578           cout << "\tFound occludee" << endl;
579         }
580 #endif
581       }
582       else {
583         fe->setOccludeeEmpty(true);
584       }
585
586       ++nSamples;
587       fe = fe->nextEdge();
588     } while ((maxCard < qiMajority) && (fe) && (fe != festart));
589
590 #if LOGGING
591     if (_global.debug & G_DEBUG_FREESTYLE) {
592       cout << "\tFinished with " << nSamples << " samples, maxCard = " << maxCard << endl;
593     }
594 #endif
595
596     // ViewEdge
597     // qi --
598     // Find the minimum value that is >= the majority of the QI
599     for (unsigned count = 0, i = 0; i < 256; ++i) {
600       count += qiClasses[i];
601       if (count >= qiMajority) {
602         (*ve)->setQI(i);
603         break;
604       }
605     }
606     // occluders --
607     // I would rather not have to go through the effort of creating this set and then copying out
608     // its contents. Is there a reason why ViewEdge::_Occluders cannot be converted to a set<>?
609     for (set<ViewShape *>::iterator o = foundOccluders.begin(), oend = foundOccluders.end();
610          o != oend;
611          ++o) {
612       (*ve)->AddOccluder((*o));
613     }
614 #if LOGGING
615     if (_global.debug & G_DEBUG_FREESTYLE) {
616       cout << "\tConclusion: QI = " << maxIndex << ", " << (*ve)->occluders_size() << " occluders."
617            << endl;
618     }
619 #else
620     (void)maxIndex;
621 #endif
622     // occludee --
623     if (!wFaces.empty()) {
624       if (wFaces.size() <= (float)nSamples / 2.0f) {
625         (*ve)->setaShape(0);
626       }
627       else {
628         ViewShape *vshape = ioViewMap->viewShape(
629             (*wFaces.begin())->GetVertex(0)->shape()->GetId());
630         (*ve)->setaShape(vshape);
631       }
632     }
633
634     wFaces.clear();
635   }
636   if (iRenderMonitor && vedges.size()) {
637     stringstream ss;
638     ss << "Freestyle: Visibility computations " << (100 * cnt / vedges.size()) << "%";
639     iRenderMonitor->setInfo(ss.str());
640     iRenderMonitor->progress((float)cnt / vedges.size());
641   }
642 }
643
644 template<typename G, typename I>
645 static void computeDetailedVisibility(ViewMap *ioViewMap,
646                                       G &grid,
647                                       real epsilon,
648                                       RenderMonitor *iRenderMonitor)
649 {
650   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
651
652   FEdge *fe, *festart;
653   int nSamples = 0;
654   vector<WFace *> wFaces;
655   WFace *wFace = NULL;
656   unsigned tmpQI = 0;
657   unsigned qiClasses[256];
658   unsigned maxIndex, maxCard;
659   unsigned qiMajority;
660   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
661     if (iRenderMonitor && iRenderMonitor->testBreak()) {
662       break;
663     }
664 #if LOGGING
665     if (_global.debug & G_DEBUG_FREESTYLE) {
666       cout << "Processing ViewEdge " << (*ve)->getId() << endl;
667     }
668 #endif
669     // Find an edge to test
670     if (!(*ve)->isInImage()) {
671       // This view edge has been proscenium culled
672       (*ve)->setQI(255);
673       (*ve)->setaShape(0);
674 #if LOGGING
675       if (_global.debug & G_DEBUG_FREESTYLE) {
676         cout << "\tCulled." << endl;
677       }
678 #endif
679       continue;
680     }
681
682     // Test edge
683     festart = (*ve)->fedgeA();
684     fe = (*ve)->fedgeA();
685     qiMajority = 0;
686     do {
687       if (fe != NULL && fe->isInImage()) {
688         qiMajority++;
689       }
690       fe = fe->nextEdge();
691     } while (fe && fe != festart);
692
693     if (qiMajority == 0) {
694       // There are no occludable FEdges on this ViewEdge
695       // This should be impossible.
696       if (_global.debug & G_DEBUG_FREESTYLE) {
697         cout << "View Edge in viewport without occludable FEdges: " << (*ve)->getId() << endl;
698       }
699       // We can recover from this error:
700       // Treat this edge as fully visible with no occludee
701       (*ve)->setQI(0);
702       (*ve)->setaShape(0);
703       continue;
704     }
705     else {
706       ++qiMajority;
707       qiMajority >>= 1;
708     }
709 #if LOGGING
710     if (_global.debug & G_DEBUG_FREESTYLE) {
711       cout << "\tqiMajority: " << qiMajority << endl;
712     }
713 #endif
714
715     tmpQI = 0;
716     maxIndex = 0;
717     maxCard = 0;
718     nSamples = 0;
719     memset(qiClasses, 0, 256 * sizeof(*qiClasses));
720     set<ViewShape *> foundOccluders;
721
722     fe = (*ve)->fedgeA();
723     do {
724       if (fe == NULL || !fe->isInImage()) {
725         fe = fe->nextEdge();
726         continue;
727       }
728       if ((maxCard < qiMajority)) {
729         // ARB: change &wFace to wFace and use reference in called function
730         tmpQI = computeVisibility<G, I>(
731             ioViewMap, fe, grid, epsilon, *ve, &wFace, &foundOccluders);
732 #if LOGGING
733         if (_global.debug & G_DEBUG_FREESTYLE) {
734           cout << "\tFEdge: visibility " << tmpQI << endl;
735         }
736 #endif
737
738         // ARB: This is an error condition, not an alert condition.
739         // Some sort of recovery or abort is necessary.
740         if (tmpQI >= 256) {
741           cerr << "Warning: too many occluding levels" << endl;
742           // ARB: Wild guess: instead of aborting or corrupting memory, treat as tmpQI == 255
743           tmpQI = 255;
744         }
745
746         if (++qiClasses[tmpQI] > maxCard) {
747           maxCard = qiClasses[tmpQI];
748           maxIndex = tmpQI;
749         }
750       }
751       else {
752         // ARB: FindOccludee is redundant if ComputeRayCastingVisibility has been called
753         // ARB: change &wFace to wFace and use reference in called function
754         findOccludee<G, I>(fe, grid, epsilon, *ve, &wFace);
755 #if LOGGING
756         if (_global.debug & G_DEBUG_FREESTYLE) {
757           cout << "\tFEdge: occludee only (" << (wFace != NULL ? "found" : "not found") << ")"
758                << endl;
759         }
760 #endif
761       }
762
763       // Store test results
764       if (wFace) {
765         vector<Vec3r> vertices;
766         for (int i = 0, numEdges = wFace->numberOfEdges(); i < numEdges; ++i) {
767           vertices.push_back(Vec3r(wFace->GetVertex(i)->GetVertex()));
768         }
769         Polygon3r poly(vertices, wFace->GetNormal());
770         poly.userdata = (void *)wFace;
771         fe->setaFace(poly);
772         wFaces.push_back(wFace);
773         fe->setOccludeeEmpty(false);
774 #if LOGGING
775         if (_global.debug & G_DEBUG_FREESTYLE) {
776           cout << "\tFound occludee" << endl;
777         }
778 #endif
779       }
780       else {
781         fe->setOccludeeEmpty(true);
782       }
783
784       ++nSamples;
785       fe = fe->nextEdge();
786     } while ((maxCard < qiMajority) && (fe) && (fe != festart));
787
788 #if LOGGING
789     if (_global.debug & G_DEBUG_FREESTYLE) {
790       cout << "\tFinished with " << nSamples << " samples, maxCard = " << maxCard << endl;
791     }
792 #endif
793
794     // ViewEdge
795     // qi --
796     (*ve)->setQI(maxIndex);
797     // occluders --
798     // I would rather not have to go through the effort of creating this this set and then copying
799     // out its contents. Is there a reason why ViewEdge::_Occluders cannot be converted to a set<>?
800     for (set<ViewShape *>::iterator o = foundOccluders.begin(), oend = foundOccluders.end();
801          o != oend;
802          ++o) {
803       (*ve)->AddOccluder((*o));
804     }
805 #if LOGGING
806     if (_global.debug & G_DEBUG_FREESTYLE) {
807       cout << "\tConclusion: QI = " << maxIndex << ", " << (*ve)->occluders_size() << " occluders."
808            << endl;
809     }
810 #endif
811     // occludee --
812     if (!wFaces.empty()) {
813       if (wFaces.size() <= (float)nSamples / 2.0f) {
814         (*ve)->setaShape(0);
815       }
816       else {
817         ViewShape *vshape = ioViewMap->viewShape(
818             (*wFaces.begin())->GetVertex(0)->shape()->GetId());
819         (*ve)->setaShape(vshape);
820       }
821     }
822
823     wFaces.clear();
824   }
825 }
826
827 template<typename G, typename I>
828 static void computeFastVisibility(ViewMap *ioViewMap, G &grid, real epsilon)
829 {
830   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
831
832   FEdge *fe, *festart;
833   unsigned nSamples = 0;
834   vector<WFace *> wFaces;
835   WFace *wFace = NULL;
836   unsigned tmpQI = 0;
837   unsigned qiClasses[256];
838   unsigned maxIndex, maxCard;
839   unsigned qiMajority;
840   bool even_test;
841   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
842     // Find an edge to test
843     if (!(*ve)->isInImage()) {
844       // This view edge has been proscenium culled
845       (*ve)->setQI(255);
846       (*ve)->setaShape(0);
847       continue;
848     }
849
850     // Test edge
851     festart = (*ve)->fedgeA();
852     fe = (*ve)->fedgeA();
853
854     even_test = true;
855     qiMajority = 0;
856     do {
857       if (even_test && fe && fe->isInImage()) {
858         qiMajority++;
859         even_test = !even_test;
860       }
861       fe = fe->nextEdge();
862     } while (fe && fe != festart);
863
864     if (qiMajority == 0) {
865       // There are no occludable FEdges on this ViewEdge
866       // This should be impossible.
867       if (_global.debug & G_DEBUG_FREESTYLE) {
868         cout << "View Edge in viewport without occludable FEdges: " << (*ve)->getId() << endl;
869       }
870       // We can recover from this error:
871       // Treat this edge as fully visible with no occludee
872       (*ve)->setQI(0);
873       (*ve)->setaShape(0);
874       continue;
875     }
876     else {
877       ++qiMajority;
878       qiMajority >>= 1;
879     }
880
881     even_test = true;
882     maxIndex = 0;
883     maxCard = 0;
884     nSamples = 0;
885     memset(qiClasses, 0, 256 * sizeof(*qiClasses));
886     set<ViewShape *> foundOccluders;
887
888     fe = (*ve)->fedgeA();
889     do {
890       if (!fe || !fe->isInImage()) {
891         fe = fe->nextEdge();
892         continue;
893       }
894       if (even_test) {
895         if ((maxCard < qiMajority)) {
896           // ARB: change &wFace to wFace and use reference in called function
897           tmpQI = computeVisibility<G, I>(
898               ioViewMap, fe, grid, epsilon, *ve, &wFace, &foundOccluders);
899
900           // ARB: This is an error condition, not an alert condition.
901           // Some sort of recovery or abort is necessary.
902           if (tmpQI >= 256) {
903             cerr << "Warning: too many occluding levels" << endl;
904             // ARB: Wild guess: instead of aborting or corrupting memory, treat as tmpQI == 255
905             tmpQI = 255;
906           }
907
908           if (++qiClasses[tmpQI] > maxCard) {
909             maxCard = qiClasses[tmpQI];
910             maxIndex = tmpQI;
911           }
912         }
913         else {
914           // ARB: FindOccludee is redundant if ComputeRayCastingVisibility has been called
915           // ARB: change &wFace to wFace and use reference in called function
916           findOccludee<G, I>(fe, grid, epsilon, *ve, &wFace);
917         }
918
919         if (wFace) {
920           vector<Vec3r> vertices;
921           for (int i = 0, numEdges = wFace->numberOfEdges(); i < numEdges; ++i) {
922             vertices.push_back(Vec3r(wFace->GetVertex(i)->GetVertex()));
923           }
924           Polygon3r poly(vertices, wFace->GetNormal());
925           poly.userdata = (void *)wFace;
926           fe->setaFace(poly);
927           wFaces.push_back(wFace);
928         }
929         ++nSamples;
930       }
931
932       even_test = !even_test;
933       fe = fe->nextEdge();
934     } while ((maxCard < qiMajority) && (fe) && (fe != festart));
935
936     // qi --
937     (*ve)->setQI(maxIndex);
938
939     // occluders --
940     for (set<ViewShape *>::iterator o = foundOccluders.begin(), oend = foundOccluders.end();
941          o != oend;
942          ++o) {
943       (*ve)->AddOccluder((*o));
944     }
945
946     // occludee --
947     if (!wFaces.empty()) {
948       if (wFaces.size() < nSamples / 2) {
949         (*ve)->setaShape(0);
950       }
951       else {
952         ViewShape *vshape = ioViewMap->viewShape(
953             (*wFaces.begin())->GetVertex(0)->shape()->GetId());
954         (*ve)->setaShape(vshape);
955       }
956     }
957
958     wFaces.clear();
959   }
960 }
961
962 template<typename G, typename I>
963 static void computeVeryFastVisibility(ViewMap *ioViewMap, G &grid, real epsilon)
964 {
965   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
966
967   FEdge *fe;
968   unsigned qi = 0;
969   WFace *wFace = 0;
970
971   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
972     // Find an edge to test
973     if (!(*ve)->isInImage()) {
974       // This view edge has been proscenium culled
975       (*ve)->setQI(255);
976       (*ve)->setaShape(0);
977       continue;
978     }
979     fe = (*ve)->fedgeA();
980     // Find a FEdge inside the occluder proscenium to test for visibility
981     FEdge *festart = fe;
982     while (fe && !fe->isInImage() && fe != festart) {
983       fe = fe->nextEdge();
984     }
985
986     // Test edge
987     if (!fe || !fe->isInImage()) {
988       // There are no occludable FEdges on this ViewEdge
989       // This should be impossible.
990       if (_global.debug & G_DEBUG_FREESTYLE) {
991         cout << "View Edge in viewport without occludable FEdges: " << (*ve)->getId() << endl;
992       }
993       // We can recover from this error:
994       // Treat this edge as fully visible with no occludee
995       qi = 0;
996       wFace = NULL;
997     }
998     else {
999       qi = computeVisibility<G, I>(ioViewMap, fe, grid, epsilon, *ve, &wFace, NULL);
1000     }
1001
1002     // Store test results
1003     if (wFace) {
1004       vector<Vec3r> vertices;
1005       for (int i = 0, numEdges = wFace->numberOfEdges(); i < numEdges; ++i) {
1006         vertices.push_back(Vec3r(wFace->GetVertex(i)->GetVertex()));
1007       }
1008       Polygon3r poly(vertices, wFace->GetNormal());
1009       poly.userdata = (void *)wFace;
1010       fe->setaFace(poly);  // This works because setaFace *copies* the polygon
1011       ViewShape *vshape = ioViewMap->viewShape(wFace->GetVertex(0)->shape()->GetId());
1012       (*ve)->setaShape(vshape);
1013     }
1014     else {
1015       (*ve)->setaShape(0);
1016     }
1017     (*ve)->setQI(qi);
1018   }
1019 }
1020
1021 void ViewMapBuilder::BuildGrid(WingedEdge &we, const BBox<Vec3r> &bbox, unsigned int sceneNumFaces)
1022 {
1023   _Grid->clear();
1024   Vec3r size;
1025   for (unsigned int i = 0; i < 3; i++) {
1026     size[i] = fabs(bbox.getMax()[i] - bbox.getMin()[i]);
1027     // let make the grid 1/10 bigger to avoid numerical errors while computing triangles/cells
1028     // intersections.
1029     size[i] += size[i] / 10.0;
1030     if (size[i] == 0) {
1031       if (_global.debug & G_DEBUG_FREESTYLE) {
1032         cout << "Warning: the bbox size is 0 in dimension " << i << endl;
1033       }
1034     }
1035   }
1036   _Grid->configure(Vec3r(bbox.getMin() - size / 20.0), size, sceneNumFaces);
1037
1038   // Fill in the grid:
1039   WFillGrid fillGridRenderer(_Grid, &we);
1040   fillGridRenderer.fillGrid();
1041
1042   // DEBUG
1043   _Grid->displayDebug();
1044 }
1045
1046 ViewMap *ViewMapBuilder::BuildViewMap(WingedEdge &we,
1047                                       visibility_algo iAlgo,
1048                                       real epsilon,
1049                                       const BBox<Vec3r> &bbox,
1050                                       unsigned int sceneNumFaces)
1051 {
1052   _ViewMap = new ViewMap;
1053   _currentId = 1;
1054   _currentFId = 0;
1055   _currentSVertexId = 0;
1056
1057   // Builds initial view edges
1058   computeInitialViewEdges(we);
1059
1060   // Detects cusps
1061   computeCusps(_ViewMap);
1062
1063   // Compute intersections
1064   ComputeIntersections(_ViewMap, sweep_line, epsilon);
1065
1066   // Compute visibility
1067   ComputeEdgesVisibility(_ViewMap, we, bbox, sceneNumFaces, iAlgo, epsilon);
1068
1069   return _ViewMap;
1070 }
1071
1072 static inline real distance2D(const Vec3r &point, const real origin[2])
1073 {
1074   return ::hypot((point[0] - origin[0]), (point[1] - origin[1]));
1075 }
1076
1077 static inline bool crossesProscenium(real proscenium[4], FEdge *fe)
1078 {
1079   Vec2r min(proscenium[0], proscenium[2]);
1080   Vec2r max(proscenium[1], proscenium[3]);
1081   Vec2r A(fe->vertexA()->getProjectedX(), fe->vertexA()->getProjectedY());
1082   Vec2r B(fe->vertexB()->getProjectedX(), fe->vertexB()->getProjectedY());
1083
1084   return GeomUtils::intersect2dSeg2dArea(min, max, A, B);
1085 }
1086
1087 static inline bool insideProscenium(real proscenium[4], const Vec3r &point)
1088 {
1089   return !(point[0] < proscenium[0] || point[0] > proscenium[1] || point[1] < proscenium[2] ||
1090            point[1] > proscenium[3]);
1091 }
1092
1093 void ViewMapBuilder::CullViewEdges(ViewMap *ioViewMap,
1094                                    real viewProscenium[4],
1095                                    real occluderProscenium[4],
1096                                    bool extensiveFEdgeSearch)
1097 {
1098   // Cull view edges by marking them as non-displayable.
1099   // This avoids the complications of trying to delete edges from the ViewMap.
1100
1101   // Non-displayable view edges will be skipped over during visibility calculation.
1102
1103   // View edges will be culled according to their position w.r.t. the viewport proscenium (viewport
1104   // + 5% border, or some such).
1105
1106   // Get proscenium boundary for culling
1107   GridHelpers::getDefaultViewProscenium(viewProscenium);
1108   real prosceniumOrigin[2];
1109   prosceniumOrigin[0] = (viewProscenium[1] - viewProscenium[0]) / 2.0;
1110   prosceniumOrigin[1] = (viewProscenium[3] - viewProscenium[2]) / 2.0;
1111   if (_global.debug & G_DEBUG_FREESTYLE) {
1112     cout << "Proscenium culling:" << endl;
1113     cout << "Proscenium: [" << viewProscenium[0] << ", " << viewProscenium[1] << ", "
1114          << viewProscenium[2] << ", " << viewProscenium[3] << "]" << endl;
1115     cout << "Origin: [" << prosceniumOrigin[0] << ", " << prosceniumOrigin[1] << "]" << endl;
1116   }
1117
1118   // A separate occluder proscenium will also be maintained, starting out the same as the viewport
1119   // proscenium, and expanding as necessary so that it encompasses the center point of at least one
1120   // feature edge in each retained view edge. The occluder proscenium will be used later to cull
1121   // occluding triangles before they are inserted into the Grid. The occluder proscenium starts out
1122   // the same size as the view proscenium
1123   GridHelpers::getDefaultViewProscenium(occluderProscenium);
1124
1125   // N.B. Freestyle is inconsistent in its use of ViewMap::viewedges_container and
1126   // vector<ViewEdge*>::iterator. Probably all occurences of vector<ViewEdge*>::iterator should be
1127   // replaced ViewMap::viewedges_container throughout the code. For each view edge
1128   ViewMap::viewedges_container::iterator ve, veend;
1129
1130   for (ve = ioViewMap->ViewEdges().begin(), veend = ioViewMap->ViewEdges().end(); ve != veend;
1131        ve++) {
1132     // Overview:
1133     //    Search for a visible feature edge
1134     //    If none: mark view edge as non-displayable
1135     //    Otherwise:
1136     //        Find a feature edge with center point inside occluder proscenium.
1137     //        If none exists, find the feature edge with center point closest to viewport origin.
1138     //            Expand occluder proscenium to enclose center point.
1139
1140     // For each feature edge, while bestOccluderTarget not found and view edge not visibile
1141     bool bestOccluderTargetFound = false;
1142     FEdge *bestOccluderTarget = NULL;
1143     real bestOccluderDistance = 0.0;
1144     FEdge *festart = (*ve)->fedgeA();
1145     FEdge *fe = festart;
1146     // All ViewEdges start culled
1147     (*ve)->setIsInImage(false);
1148
1149     // For simple visibility calculation: mark a feature edge that is known to have a center point
1150     // inside the occluder proscenium. Cull all other feature edges.
1151     do {
1152       // All FEdges start culled
1153       fe->setIsInImage(false);
1154
1155       // Look for the visible edge that can most easily be included in the occluder proscenium.
1156       if (!bestOccluderTargetFound) {
1157         // If center point is inside occluder proscenium,
1158         if (insideProscenium(occluderProscenium, fe->center2d())) {
1159           // Use this feature edge for visibility deterimination
1160           fe->setIsInImage(true);
1161           // Mark bestOccluderTarget as found
1162           bestOccluderTargetFound = true;
1163           bestOccluderTarget = fe;
1164         }
1165         else {
1166           real d = distance2D(fe->center2d(), prosceniumOrigin);
1167           // If center point is closer to viewport origin than current target
1168           if (bestOccluderTarget == NULL || d < bestOccluderDistance) {
1169             // Then store as bestOccluderTarget
1170             bestOccluderDistance = d;
1171             bestOccluderTarget = fe;
1172           }
1173         }
1174       }
1175
1176       // If feature edge crosses the view proscenium
1177       if (!(*ve)->isInImage() && crossesProscenium(viewProscenium, fe)) {
1178         // Then the view edge will be included in the image
1179         (*ve)->setIsInImage(true);
1180       }
1181       fe = fe->nextEdge();
1182     } while (fe && fe != festart && !(bestOccluderTargetFound && (*ve)->isInImage()));
1183
1184     // Either we have run out of FEdges, or we already have the one edge we need to determine
1185     // visibility Cull all remaining edges.
1186     while (fe && fe != festart) {
1187       fe->setIsInImage(false);
1188       fe = fe->nextEdge();
1189     }
1190
1191     // If bestOccluderTarget was not found inside the occluder proscenium, we need to expand the
1192     // occluder proscenium to include it.
1193     if ((*ve)->isInImage() && bestOccluderTarget != NULL && !bestOccluderTargetFound) {
1194       // Expand occluder proscenium to enclose bestOccluderTarget
1195       Vec3r point = bestOccluderTarget->center2d();
1196       if (point[0] < occluderProscenium[0]) {
1197         occluderProscenium[0] = point[0];
1198       }
1199       else if (point[0] > occluderProscenium[1]) {
1200         occluderProscenium[1] = point[0];
1201       }
1202       if (point[1] < occluderProscenium[2]) {
1203         occluderProscenium[2] = point[1];
1204       }
1205       else if (point[1] > occluderProscenium[3]) {
1206         occluderProscenium[3] = point[1];
1207       }
1208       // Use bestOccluderTarget for visibility determination
1209       bestOccluderTarget->setIsInImage(true);
1210     }
1211   }
1212
1213   // We are done calculating the occluder proscenium.
1214   // Expand the occluder proscenium by an epsilon to avoid rounding errors.
1215   const real epsilon = 1.0e-6;
1216   occluderProscenium[0] -= epsilon;
1217   occluderProscenium[1] += epsilon;
1218   occluderProscenium[2] -= epsilon;
1219   occluderProscenium[3] += epsilon;
1220
1221   // For "Normal" or "Fast" style visibility computation only:
1222
1223   // For more detailed visibility calculation, make a second pass through the view map, marking all
1224   // feature edges with center points inside the final occluder proscenium. All of these feature
1225   // edges can be considered during visibility calculation.
1226
1227   // So far we have only found one FEdge per ViewEdge.  The "Normal" and "Fast" styles of
1228   // visibility computation want to consider many FEdges for each ViewEdge. Here we re-scan the
1229   // view map to find any usable FEdges that we skipped on the first pass, or that have become
1230   // usable because the occluder proscenium has been expanded since the edge was visited on the
1231   // first pass.
1232   if (extensiveFEdgeSearch) {
1233     // For each view edge,
1234     for (ve = ioViewMap->ViewEdges().begin(), veend = ioViewMap->ViewEdges().end(); ve != veend;
1235          ve++) {
1236       if (!(*ve)->isInImage()) {
1237         continue;
1238       }
1239       // For each feature edge,
1240       FEdge *festart = (*ve)->fedgeA();
1241       FEdge *fe = festart;
1242       do {
1243         // If not (already) visible and center point inside occluder proscenium,
1244         if (!fe->isInImage() && insideProscenium(occluderProscenium, fe->center2d())) {
1245           // Use the feature edge for visibility determination
1246           fe->setIsInImage(true);
1247         }
1248         fe = fe->nextEdge();
1249       } while (fe && fe != festart);
1250     }
1251   }
1252 }
1253
1254 void ViewMapBuilder::computeInitialViewEdges(WingedEdge &we)
1255 {
1256   vector<WShape *> wshapes = we.getWShapes();
1257   SShape *psShape;
1258
1259   for (vector<WShape *>::const_iterator it = wshapes.begin(); it != wshapes.end(); it++) {
1260     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
1261       break;
1262     }
1263
1264     // create the embedding
1265     psShape = new SShape;
1266     psShape->setId((*it)->GetId());
1267     psShape->setName((*it)->getName());
1268     psShape->setLibraryPath((*it)->getLibraryPath());
1269     psShape->setFrsMaterials((*it)->frs_materials());  // FIXME
1270
1271     // create the view shape
1272     ViewShape *vshape = new ViewShape(psShape);
1273     // add this view shape to the view map:
1274     _ViewMap->AddViewShape(vshape);
1275
1276     // we want to number the view edges in a unique way for the while scene.
1277     _pViewEdgeBuilder->setCurrentViewId(_currentId);
1278     // we want to number the feature edges in a unique way for the while scene.
1279     _pViewEdgeBuilder->setCurrentFId(_currentFId);
1280     // we want to number the SVertex in a unique way for the while scene.
1281     _pViewEdgeBuilder->setCurrentSVertexId(_currentFId);
1282     _pViewEdgeBuilder->BuildViewEdges(dynamic_cast<WXShape *>(*it),
1283                                       vshape,
1284                                       _ViewMap->ViewEdges(),
1285                                       _ViewMap->ViewVertices(),
1286                                       _ViewMap->FEdges(),
1287                                       _ViewMap->SVertices());
1288
1289     _currentId = _pViewEdgeBuilder->currentViewId() + 1;
1290     _currentFId = _pViewEdgeBuilder->currentFId() + 1;
1291     _currentSVertexId = _pViewEdgeBuilder->currentSVertexId() + 1;
1292
1293     psShape->ComputeBBox();
1294   }
1295 }
1296
1297 void ViewMapBuilder::computeCusps(ViewMap *ioViewMap)
1298 {
1299   vector<ViewVertex *> newVVertices;
1300   vector<ViewEdge *> newVEdges;
1301   ViewMap::viewedges_container &vedges = ioViewMap->ViewEdges();
1302   ViewMap::viewedges_container::iterator ve = vedges.begin(), veend = vedges.end();
1303   for (; ve != veend; ++ve) {
1304     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
1305       break;
1306     }
1307     if ((!((*ve)->getNature() & Nature::SILHOUETTE)) || (!((*ve)->fedgeA()->isSmooth()))) {
1308       continue;
1309     }
1310     FEdge *fe = (*ve)->fedgeA();
1311     FEdge *fefirst = fe;
1312     bool first = true;
1313     bool positive = true;
1314     do {
1315       FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
1316       Vec3r A((fes)->vertexA()->point3d());
1317       Vec3r B((fes)->vertexB()->point3d());
1318       Vec3r AB(B - A);
1319       AB.normalize();
1320       Vec3r m((A + B) / 2.0);
1321       Vec3r crossP(AB ^ (fes)->normal());
1322       crossP.normalize();
1323       Vec3r viewvector;
1324       if (_orthographicProjection) {
1325         viewvector = Vec3r(0.0, 0.0, m.z() - _viewpoint.z());
1326       }
1327       else {
1328         viewvector = Vec3r(m - _viewpoint);
1329       }
1330       viewvector.normalize();
1331       if (first) {
1332         if (((crossP) * (viewvector)) > 0) {
1333           positive = true;
1334         }
1335         else {
1336           positive = false;
1337         }
1338         first = false;
1339       }
1340       // If we're in a positive part, we need a stronger negative value to change
1341       NonTVertex *cusp = NULL;
1342       if (positive) {
1343         if (((crossP) * (viewvector)) < -0.1) {
1344           // state changes
1345           positive = false;
1346           // creates and insert cusp
1347           cusp = dynamic_cast<NonTVertex *>(
1348               ioViewMap->InsertViewVertex(fes->vertexA(), newVEdges));
1349           if (cusp) {
1350             cusp->setNature(cusp->getNature() | Nature::CUSP);
1351           }
1352         }
1353       }
1354       else {
1355         // If we're in a negative part, we need a stronger negative value to change
1356         if (((crossP) * (viewvector)) > 0.1) {
1357           positive = true;
1358           cusp = dynamic_cast<NonTVertex *>(
1359               ioViewMap->InsertViewVertex(fes->vertexA(), newVEdges));
1360           if (cusp) {
1361             cusp->setNature(cusp->getNature() | Nature::CUSP);
1362           }
1363         }
1364       }
1365       fe = fe->nextEdge();
1366     } while (fe && fe != fefirst);
1367   }
1368   for (ve = newVEdges.begin(), veend = newVEdges.end(); ve != veend; ++ve) {
1369     (*ve)->viewShape()->AddEdge(*ve);
1370     vedges.push_back(*ve);
1371   }
1372 }
1373
1374 void ViewMapBuilder::ComputeCumulativeVisibility(ViewMap *ioViewMap,
1375                                                  WingedEdge &we,
1376                                                  const BBox<Vec3r> &bbox,
1377                                                  real epsilon,
1378                                                  bool cull,
1379                                                  GridDensityProviderFactory &factory)
1380 {
1381   AutoPtr<GridHelpers::Transform> transform;
1382   AutoPtr<OccluderSource> source;
1383
1384   if (_orthographicProjection) {
1385     transform.reset(new BoxGrid::Transform);
1386   }
1387   else {
1388     transform.reset(new SphericalGrid::Transform);
1389   }
1390
1391   if (cull) {
1392     source.reset(new CulledOccluderSource(*transform, we, *ioViewMap, true));
1393   }
1394   else {
1395     source.reset(new OccluderSource(*transform, we));
1396   }
1397
1398   AutoPtr<GridDensityProvider> density(factory.newGridDensityProvider(*source, bbox, *transform));
1399
1400   if (_orthographicProjection) {
1401     BoxGrid grid(*source, *density, ioViewMap, _viewpoint, _EnableQI);
1402     computeCumulativeVisibility<BoxGrid, BoxGrid::Iterator>(
1403         ioViewMap, grid, epsilon, _pRenderMonitor);
1404   }
1405   else {
1406     SphericalGrid grid(*source, *density, ioViewMap, _viewpoint, _EnableQI);
1407     computeCumulativeVisibility<SphericalGrid, SphericalGrid::Iterator>(
1408         ioViewMap, grid, epsilon, _pRenderMonitor);
1409   }
1410 }
1411
1412 void ViewMapBuilder::ComputeDetailedVisibility(ViewMap *ioViewMap,
1413                                                WingedEdge &we,
1414                                                const BBox<Vec3r> &bbox,
1415                                                real epsilon,
1416                                                bool cull,
1417                                                GridDensityProviderFactory &factory)
1418 {
1419   AutoPtr<GridHelpers::Transform> transform;
1420   AutoPtr<OccluderSource> source;
1421
1422   if (_orthographicProjection) {
1423     transform.reset(new BoxGrid::Transform);
1424   }
1425   else {
1426     transform.reset(new SphericalGrid::Transform);
1427   }
1428
1429   if (cull) {
1430     source.reset(new CulledOccluderSource(*transform, we, *ioViewMap, true));
1431   }
1432   else {
1433     source.reset(new OccluderSource(*transform, we));
1434   }
1435
1436   AutoPtr<GridDensityProvider> density(factory.newGridDensityProvider(*source, bbox, *transform));
1437
1438   if (_orthographicProjection) {
1439     BoxGrid grid(*source, *density, ioViewMap, _viewpoint, _EnableQI);
1440     computeDetailedVisibility<BoxGrid, BoxGrid::Iterator>(
1441         ioViewMap, grid, epsilon, _pRenderMonitor);
1442   }
1443   else {
1444     SphericalGrid grid(*source, *density, ioViewMap, _viewpoint, _EnableQI);
1445     computeDetailedVisibility<SphericalGrid, SphericalGrid::Iterator>(
1446         ioViewMap, grid, epsilon, _pRenderMonitor);
1447   }
1448 }
1449
1450 void ViewMapBuilder::ComputeEdgesVisibility(ViewMap *ioViewMap,
1451                                             WingedEdge &we,
1452                                             const BBox<Vec3r> &bbox,
1453                                             unsigned int sceneNumFaces,
1454                                             visibility_algo iAlgo,
1455                                             real epsilon)
1456 {
1457 #if 0
1458   iAlgo = ray_casting;  // for testing algorithms equivalence
1459 #endif
1460   switch (iAlgo) {
1461     case ray_casting:
1462       if (_global.debug & G_DEBUG_FREESTYLE) {
1463         cout << "Using ordinary ray casting" << endl;
1464       }
1465       BuildGrid(we, bbox, sceneNumFaces);
1466       ComputeRayCastingVisibility(ioViewMap, epsilon);
1467       break;
1468     case ray_casting_fast:
1469       if (_global.debug & G_DEBUG_FREESTYLE) {
1470         cout << "Using fast ray casting" << endl;
1471       }
1472       BuildGrid(we, bbox, sceneNumFaces);
1473       ComputeFastRayCastingVisibility(ioViewMap, epsilon);
1474       break;
1475     case ray_casting_very_fast:
1476       if (_global.debug & G_DEBUG_FREESTYLE) {
1477         cout << "Using very fast ray casting" << endl;
1478       }
1479       BuildGrid(we, bbox, sceneNumFaces);
1480       ComputeVeryFastRayCastingVisibility(ioViewMap, epsilon);
1481       break;
1482     case ray_casting_culled_adaptive_traditional:
1483       if (_global.debug & G_DEBUG_FREESTYLE) {
1484         cout << "Using culled adaptive grid with heuristic density and traditional QI calculation"
1485              << endl;
1486       }
1487       try {
1488         HeuristicGridDensityProviderFactory factory(0.5f, sceneNumFaces);
1489         ComputeDetailedVisibility(ioViewMap, we, bbox, epsilon, true, factory);
1490       }
1491       catch (...) {
1492         // Last resort catch to make sure RAII semantics hold for OptimizedGrid. Can be replaced
1493         // with try...catch block around main() if the program as a whole is converted to RAII
1494
1495         // This is the little-mentioned caveat of RAII: RAII does not work unless destructors are
1496         // always called, but destructors are only called if all exceptions are caught (or
1497         // std::terminate() is replaced).
1498
1499         // We don't actually handle the exception here, so re-throw it now that our destructors
1500         // have had a chance to run.
1501         throw;
1502       }
1503       break;
1504     case ray_casting_adaptive_traditional:
1505       if (_global.debug & G_DEBUG_FREESTYLE) {
1506         cout
1507             << "Using unculled adaptive grid with heuristic density and traditional QI calculation"
1508             << endl;
1509       }
1510       try {
1511         HeuristicGridDensityProviderFactory factory(0.5f, sceneNumFaces);
1512         ComputeDetailedVisibility(ioViewMap, we, bbox, epsilon, false, factory);
1513       }
1514       catch (...) {
1515         throw;
1516       }
1517       break;
1518     case ray_casting_culled_adaptive_cumulative:
1519       if (_global.debug & G_DEBUG_FREESTYLE) {
1520         cout << "Using culled adaptive grid with heuristic density and cumulative QI calculation"
1521              << endl;
1522       }
1523       try {
1524         HeuristicGridDensityProviderFactory factory(0.5f, sceneNumFaces);
1525         ComputeCumulativeVisibility(ioViewMap, we, bbox, epsilon, true, factory);
1526       }
1527       catch (...) {
1528         throw;
1529       }
1530       break;
1531     case ray_casting_adaptive_cumulative:
1532       if (_global.debug & G_DEBUG_FREESTYLE) {
1533         cout << "Using unculled adaptive grid with heuristic density and cumulative QI calculation"
1534              << endl;
1535       }
1536       try {
1537         HeuristicGridDensityProviderFactory factory(0.5f, sceneNumFaces);
1538         ComputeCumulativeVisibility(ioViewMap, we, bbox, epsilon, false, factory);
1539       }
1540       catch (...) {
1541         throw;
1542       }
1543       break;
1544     default:
1545       break;
1546   }
1547 }
1548
1549 static const unsigned gProgressBarMaxSteps = 10;
1550 static const unsigned gProgressBarMinSize = 2000;
1551
1552 void ViewMapBuilder::ComputeRayCastingVisibility(ViewMap *ioViewMap, real epsilon)
1553 {
1554   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
1555   bool progressBarDisplay = false;
1556   unsigned progressBarStep = 0;
1557   unsigned vEdgesSize = vedges.size();
1558   unsigned fEdgesSize = ioViewMap->FEdges().size();
1559
1560   if (_pProgressBar != NULL && fEdgesSize > gProgressBarMinSize) {
1561     unsigned progressBarSteps = min(gProgressBarMaxSteps, vEdgesSize);
1562     progressBarStep = vEdgesSize / progressBarSteps;
1563     _pProgressBar->reset();
1564     _pProgressBar->setLabelText("Computing Ray casting Visibility");
1565     _pProgressBar->setTotalSteps(progressBarSteps);
1566     _pProgressBar->setProgress(0);
1567     progressBarDisplay = true;
1568   }
1569
1570   unsigned counter = progressBarStep;
1571   FEdge *fe, *festart;
1572   int nSamples = 0;
1573   vector<Polygon3r *> aFaces;
1574   Polygon3r *aFace = NULL;
1575   unsigned tmpQI = 0;
1576   unsigned qiClasses[256];
1577   unsigned maxIndex, maxCard;
1578   unsigned qiMajority;
1579   static unsigned timestamp = 1;
1580   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
1581     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
1582       break;
1583     }
1584 #if LOGGING
1585     if (_global.debug & G_DEBUG_FREESTYLE) {
1586       cout << "Processing ViewEdge " << (*ve)->getId() << endl;
1587     }
1588 #endif
1589     festart = (*ve)->fedgeA();
1590     fe = (*ve)->fedgeA();
1591     qiMajority = 1;
1592     do {
1593       qiMajority++;
1594       fe = fe->nextEdge();
1595     } while (fe && fe != festart);
1596     qiMajority >>= 1;
1597 #if LOGGING
1598     if (_global.debug & G_DEBUG_FREESTYLE) {
1599       cout << "\tqiMajority: " << qiMajority << endl;
1600     }
1601 #endif
1602
1603     tmpQI = 0;
1604     maxIndex = 0;
1605     maxCard = 0;
1606     nSamples = 0;
1607     fe = (*ve)->fedgeA();
1608     memset(qiClasses, 0, 256 * sizeof(*qiClasses));
1609     set<ViewShape *> occluders;
1610     do {
1611       if ((maxCard < qiMajority)) {
1612         tmpQI = ComputeRayCastingVisibility(fe, _Grid, epsilon, occluders, &aFace, timestamp++);
1613
1614 #if LOGGING
1615         if (_global.debug & G_DEBUG_FREESTYLE) {
1616           cout << "\tFEdge: visibility " << tmpQI << endl;
1617         }
1618 #endif
1619         // ARB: This is an error condition, not an alert condition.
1620         // Some sort of recovery or abort is necessary.
1621         if (tmpQI >= 256) {
1622           cerr << "Warning: too many occluding levels" << endl;
1623           // ARB: Wild guess: instead of aborting or corrupting memory, treat as tmpQI == 255
1624           tmpQI = 255;
1625         }
1626
1627         if (++qiClasses[tmpQI] > maxCard) {
1628           maxCard = qiClasses[tmpQI];
1629           maxIndex = tmpQI;
1630         }
1631       }
1632       else {
1633         // ARB: FindOccludee is redundant if ComputeRayCastingVisibility has been called
1634         FindOccludee(fe, _Grid, epsilon, &aFace, timestamp++);
1635 #if LOGGING
1636         if (_global.debug & G_DEBUG_FREESTYLE) {
1637           cout << "\tFEdge: occludee only (" << (aFace != NULL ? "found" : "not found") << ")"
1638                << endl;
1639         }
1640 #endif
1641       }
1642
1643       if (aFace) {
1644         fe->setaFace(*aFace);
1645         aFaces.push_back(aFace);
1646         fe->setOccludeeEmpty(false);
1647 #if LOGGING
1648         if (_global.debug & G_DEBUG_FREESTYLE) {
1649           cout << "\tFound occludee" << endl;
1650         }
1651 #endif
1652       }
1653       else {
1654         // ARB: We are arbitrarily using the last observed value for occludee (almost always the
1655         // value observed
1656         //     for the edge before festart). Is that meaningful?
1657         // ...in fact, _occludeeEmpty seems to be unused.
1658         fe->setOccludeeEmpty(true);
1659       }
1660
1661       ++nSamples;
1662       fe = fe->nextEdge();
1663     } while ((maxCard < qiMajority) && (fe) && (fe != festart));
1664 #if LOGGING
1665     if (_global.debug & G_DEBUG_FREESTYLE) {
1666       cout << "\tFinished with " << nSamples << " samples, maxCard = " << maxCard << endl;
1667     }
1668 #endif
1669
1670     // ViewEdge
1671     // qi --
1672     (*ve)->setQI(maxIndex);
1673     // occluders --
1674     for (set<ViewShape *>::iterator o = occluders.begin(), oend = occluders.end(); o != oend;
1675          ++o) {
1676       (*ve)->AddOccluder((*o));
1677     }
1678 #if LOGGING
1679     if (_global.debug & G_DEBUG_FREESTYLE) {
1680       cout << "\tConclusion: QI = " << maxIndex << ", " << (*ve)->occluders_size() << " occluders."
1681            << endl;
1682     }
1683 #endif
1684     // occludee --
1685     if (!aFaces.empty()) {
1686       if (aFaces.size() <= (float)nSamples / 2.0f) {
1687         (*ve)->setaShape(0);
1688       }
1689       else {
1690         vector<Polygon3r *>::iterator p = aFaces.begin();
1691         WFace *wface = (WFace *)((*p)->userdata);
1692         ViewShape *vshape = ioViewMap->viewShape(wface->GetVertex(0)->shape()->GetId());
1693         ++p;
1694         (*ve)->setaShape(vshape);
1695       }
1696     }
1697
1698     if (progressBarDisplay) {
1699       counter--;
1700       if (counter <= 0) {
1701         counter = progressBarStep;
1702         _pProgressBar->setProgress(_pProgressBar->getProgress() + 1);
1703       }
1704     }
1705     aFaces.clear();
1706   }
1707 }
1708
1709 void ViewMapBuilder::ComputeFastRayCastingVisibility(ViewMap *ioViewMap, real epsilon)
1710 {
1711   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
1712   bool progressBarDisplay = false;
1713   unsigned progressBarStep = 0;
1714   unsigned vEdgesSize = vedges.size();
1715   unsigned fEdgesSize = ioViewMap->FEdges().size();
1716
1717   if (_pProgressBar != NULL && fEdgesSize > gProgressBarMinSize) {
1718     unsigned progressBarSteps = min(gProgressBarMaxSteps, vEdgesSize);
1719     progressBarStep = vEdgesSize / progressBarSteps;
1720     _pProgressBar->reset();
1721     _pProgressBar->setLabelText("Computing Ray casting Visibility");
1722     _pProgressBar->setTotalSteps(progressBarSteps);
1723     _pProgressBar->setProgress(0);
1724     progressBarDisplay = true;
1725   }
1726
1727   unsigned counter = progressBarStep;
1728   FEdge *fe, *festart;
1729   unsigned nSamples = 0;
1730   vector<Polygon3r *> aFaces;
1731   Polygon3r *aFace = NULL;
1732   unsigned tmpQI = 0;
1733   unsigned qiClasses[256];
1734   unsigned maxIndex, maxCard;
1735   unsigned qiMajority;
1736   static unsigned timestamp = 1;
1737   bool even_test;
1738   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
1739     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
1740       break;
1741     }
1742
1743     festart = (*ve)->fedgeA();
1744     fe = (*ve)->fedgeA();
1745     qiMajority = 1;
1746     do {
1747       qiMajority++;
1748       fe = fe->nextEdge();
1749     } while (fe && fe != festart);
1750     if (qiMajority >= 4) {
1751       qiMajority >>= 2;
1752     }
1753     else {
1754       qiMajority = 1;
1755     }
1756
1757     set<ViewShape *> occluders;
1758
1759     even_test = true;
1760     maxIndex = 0;
1761     maxCard = 0;
1762     nSamples = 0;
1763     memset(qiClasses, 0, 256 * sizeof(*qiClasses));
1764     fe = (*ve)->fedgeA();
1765     do {
1766       if (even_test) {
1767         if ((maxCard < qiMajority)) {
1768           tmpQI = ComputeRayCastingVisibility(fe, _Grid, epsilon, occluders, &aFace, timestamp++);
1769
1770           // ARB: This is an error condition, not an alert condition.
1771           // Some sort of recovery or abort is necessary.
1772           if (tmpQI >= 256) {
1773             cerr << "Warning: too many occluding levels" << endl;
1774             // ARB: Wild guess: instead of aborting or corrupting memory, treat as tmpQI == 255
1775             tmpQI = 255;
1776           }
1777
1778           if (++qiClasses[tmpQI] > maxCard) {
1779             maxCard = qiClasses[tmpQI];
1780             maxIndex = tmpQI;
1781           }
1782         }
1783         else {
1784           // ARB: FindOccludee is redundant if ComputeRayCastingVisibility has been called
1785           FindOccludee(fe, _Grid, epsilon, &aFace, timestamp++);
1786         }
1787
1788         if (aFace) {
1789           fe->setaFace(*aFace);
1790           aFaces.push_back(aFace);
1791         }
1792         ++nSamples;
1793         even_test = false;
1794       }
1795       else {
1796         even_test = true;
1797       }
1798       fe = fe->nextEdge();
1799     } while ((maxCard < qiMajority) && (fe) && (fe != festart));
1800
1801     (*ve)->setQI(maxIndex);
1802
1803     if (!aFaces.empty()) {
1804       if (aFaces.size() < nSamples / 2) {
1805         (*ve)->setaShape(0);
1806       }
1807       else {
1808         vector<Polygon3r *>::iterator p = aFaces.begin();
1809         WFace *wface = (WFace *)((*p)->userdata);
1810         ViewShape *vshape = ioViewMap->viewShape(wface->GetVertex(0)->shape()->GetId());
1811         ++p;
1812 #if 0
1813         for (; p != pend; ++p) {
1814           WFace *f = (WFace *)((*p)->userdata);
1815           ViewShape *vs = ioViewMap->viewShape(f->GetVertex(0)->shape()->GetId());
1816           if (vs != vshape) {
1817             sameShape = false;
1818             break;
1819           }
1820         }
1821         if (sameShape)
1822 #endif
1823         {
1824           (*ve)->setaShape(vshape);
1825         }
1826       }
1827     }
1828
1829     //(*ve)->setaFace(aFace);
1830
1831     if (progressBarDisplay) {
1832       counter--;
1833       if (counter <= 0) {
1834         counter = progressBarStep;
1835         _pProgressBar->setProgress(_pProgressBar->getProgress() + 1);
1836       }
1837     }
1838     aFaces.clear();
1839   }
1840 }
1841
1842 void ViewMapBuilder::ComputeVeryFastRayCastingVisibility(ViewMap *ioViewMap, real epsilon)
1843 {
1844   vector<ViewEdge *> &vedges = ioViewMap->ViewEdges();
1845   bool progressBarDisplay = false;
1846   unsigned progressBarStep = 0;
1847   unsigned vEdgesSize = vedges.size();
1848   unsigned fEdgesSize = ioViewMap->FEdges().size();
1849
1850   if (_pProgressBar != NULL && fEdgesSize > gProgressBarMinSize) {
1851     unsigned progressBarSteps = min(gProgressBarMaxSteps, vEdgesSize);
1852     progressBarStep = vEdgesSize / progressBarSteps;
1853     _pProgressBar->reset();
1854     _pProgressBar->setLabelText("Computing Ray casting Visibility");
1855     _pProgressBar->setTotalSteps(progressBarSteps);
1856     _pProgressBar->setProgress(0);
1857     progressBarDisplay = true;
1858   }
1859
1860   unsigned counter = progressBarStep;
1861   FEdge *fe;
1862   unsigned qi = 0;
1863   Polygon3r *aFace = NULL;
1864   static unsigned timestamp = 1;
1865   for (vector<ViewEdge *>::iterator ve = vedges.begin(), veend = vedges.end(); ve != veend; ve++) {
1866     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
1867       break;
1868     }
1869
1870     set<ViewShape *> occluders;
1871
1872     fe = (*ve)->fedgeA();
1873     qi = ComputeRayCastingVisibility(fe, _Grid, epsilon, occluders, &aFace, timestamp++);
1874     if (aFace) {
1875       fe->setaFace(*aFace);
1876       WFace *wface = (WFace *)(aFace->userdata);
1877       ViewShape *vshape = ioViewMap->viewShape(wface->GetVertex(0)->shape()->GetId());
1878       (*ve)->setaShape(vshape);
1879     }
1880     else {
1881       (*ve)->setaShape(0);
1882     }
1883
1884     (*ve)->setQI(qi);
1885
1886     if (progressBarDisplay) {
1887       counter--;
1888       if (counter <= 0) {
1889         counter = progressBarStep;
1890         _pProgressBar->setProgress(_pProgressBar->getProgress() + 1);
1891       }
1892     }
1893   }
1894 }
1895
1896 void ViewMapBuilder::FindOccludee(FEdge *fe,
1897                                   Grid *iGrid,
1898                                   real epsilon,
1899                                   Polygon3r **oaPolygon,
1900                                   unsigned timestamp,
1901                                   Vec3r &u,
1902                                   Vec3r &A,
1903                                   Vec3r &origin,
1904                                   Vec3r &edgeDir,
1905                                   vector<WVertex *> &faceVertices)
1906 {
1907   WFace *face = NULL;
1908   if (fe->isSmooth()) {
1909     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
1910     face = (WFace *)fes->face();
1911   }
1912   OccludersSet occluders;
1913   WFace *oface;
1914   bool skipFace;
1915
1916   WVertex::incoming_edge_iterator ie;
1917   OccludersSet::iterator p, pend;
1918
1919   *oaPolygon = NULL;
1920   if (((fe)->getNature() & Nature::SILHOUETTE) || ((fe)->getNature() & Nature::BORDER)) {
1921     occluders.clear();
1922     // we cast a ray from A in the same direction but looking behind
1923     Vec3r v(-u[0], -u[1], -u[2]);
1924     iGrid->castInfiniteRay(A, v, occluders, timestamp);
1925
1926     bool noIntersection = true;
1927     real mint = FLT_MAX;
1928     // we met some occluders, let us fill the aShape field with the first intersected occluder
1929     for (p = occluders.begin(), pend = occluders.end(); p != pend; p++) {
1930       // check whether the edge and the polygon plane are coincident:
1931       //-------------------------------------------------------------
1932       // first let us compute the plane equation.
1933       oface = (WFace *)(*p)->userdata;
1934       Vec3r v1(((*p)->getVertices())[0]);
1935       Vec3r normal((*p)->getNormal());
1936       real d = -(v1 * normal);
1937       real t, t_u, t_v;
1938
1939       if (face) {
1940         skipFace = false;
1941
1942         if (face == oface) {
1943           continue;
1944         }
1945
1946         if (faceVertices.empty()) {
1947           continue;
1948         }
1949
1950         for (vector<WVertex *>::iterator fv = faceVertices.begin(), fvend = faceVertices.end();
1951              fv != fvend;
1952              ++fv) {
1953           if ((*fv)->isBoundary()) {
1954             continue;
1955           }
1956           WVertex::incoming_edge_iterator iebegin = (*fv)->incoming_edges_begin();
1957           WVertex::incoming_edge_iterator ieend = (*fv)->incoming_edges_end();
1958           for (ie = iebegin; ie != ieend; ++ie) {
1959             if ((*ie) == 0) {
1960               continue;
1961             }
1962
1963             WFace *sface = (*ie)->GetbFace();
1964             if (sface == oface) {
1965               skipFace = true;
1966               break;
1967             }
1968           }
1969           if (skipFace) {
1970             break;
1971           }
1972         }
1973         if (skipFace) {
1974           continue;
1975         }
1976       }
1977       else {
1978         if (GeomUtils::COINCIDENT ==
1979             GeomUtils::intersectRayPlane(origin, edgeDir, normal, d, t, epsilon)) {
1980           continue;
1981         }
1982       }
1983       if ((*p)->rayIntersect(A, v, t, t_u, t_v)) {
1984         if (fabs(v * normal) > 0.0001) {
1985           if (t > 0.0) {  // && t < 1.0) {
1986             if (t < mint) {
1987               *oaPolygon = (*p);
1988               mint = t;
1989               noIntersection = false;
1990               fe->setOccludeeIntersection(Vec3r(A + t * v));
1991             }
1992           }
1993         }
1994       }
1995     }
1996
1997     if (noIntersection) {
1998       *oaPolygon = NULL;
1999     }
2000   }
2001 }
2002
2003 void ViewMapBuilder::FindOccludee(
2004     FEdge *fe, Grid *iGrid, real epsilon, Polygon3r **oaPolygon, unsigned timestamp)
2005 {
2006   OccludersSet occluders;
2007
2008   Vec3r A;
2009   Vec3r edgeDir;
2010   Vec3r origin;
2011   A = Vec3r(((fe)->vertexA()->point3D() + (fe)->vertexB()->point3D()) / 2.0);
2012   edgeDir = Vec3r((fe)->vertexB()->point3D() - (fe)->vertexA()->point3D());
2013   edgeDir.normalize();
2014   origin = Vec3r((fe)->vertexA()->point3D());
2015   Vec3r u;
2016   if (_orthographicProjection) {
2017     u = Vec3r(0.0, 0.0, _viewpoint.z() - A.z());
2018   }
2019   else {
2020     u = Vec3r(_viewpoint - A);
2021   }
2022   u.normalize();
2023   if (A < iGrid->getOrigin()) {
2024     cerr << "Warning: point is out of the grid for fedge " << fe->getId().getFirst() << "-"
2025          << fe->getId().getSecond() << endl;
2026   }
2027
2028   vector<WVertex *> faceVertices;
2029
2030   WFace *face = NULL;
2031   if (fe->isSmooth()) {
2032     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
2033     face = (WFace *)fes->face();
2034   }
2035   if (face) {
2036     face->RetrieveVertexList(faceVertices);
2037   }
2038
2039   return FindOccludee(
2040       fe, iGrid, epsilon, oaPolygon, timestamp, u, A, origin, edgeDir, faceVertices);
2041 }
2042
2043 int ViewMapBuilder::ComputeRayCastingVisibility(FEdge *fe,
2044                                                 Grid *iGrid,
2045                                                 real epsilon,
2046                                                 set<ViewShape *> &oOccluders,
2047                                                 Polygon3r **oaPolygon,
2048                                                 unsigned timestamp)
2049 {
2050   OccludersSet occluders;
2051   int qi = 0;
2052
2053   Vec3r center;
2054   Vec3r edgeDir;
2055   Vec3r origin;
2056
2057   center = fe->center3d();
2058   edgeDir = Vec3r(fe->vertexB()->point3D() - fe->vertexA()->point3D());
2059   edgeDir.normalize();
2060   origin = Vec3r(fe->vertexA()->point3D());
2061   // Is the edge outside the view frustum ?
2062   Vec3r gridOrigin(iGrid->getOrigin());
2063   Vec3r gridExtremity(iGrid->getOrigin() + iGrid->gridSize());
2064
2065   if ((center.x() < gridOrigin.x()) || (center.y() < gridOrigin.y()) ||
2066       (center.z() < gridOrigin.z()) || (center.x() > gridExtremity.x()) ||
2067       (center.y() > gridExtremity.y()) || (center.z() > gridExtremity.z())) {
2068     cerr << "Warning: point is out of the grid for fedge " << fe->getId() << endl;
2069     // return 0;
2070   }
2071
2072 #if 0
2073   Vec3r A(fe->vertexA()->point2d());
2074   Vec3r B(fe->vertexB()->point2d());
2075   int viewport[4];
2076   SilhouetteGeomEngine::retrieveViewport(viewport);
2077   if ((A.x() < viewport[0]) || (A.x() > viewport[2]) || (A.y() < viewport[1]) ||
2078       (A.y() > viewport[3]) || (B.x() < viewport[0]) || (B.x() > viewport[2]) ||
2079       (B.y() < viewport[1]) || (B.y() > viewport[3])) {
2080     cerr << "Warning: point is out of the grid for fedge " << fe->getId() << endl;
2081     //return 0;
2082   }
2083 #endif
2084
2085   Vec3r vp;
2086   if (_orthographicProjection) {
2087     vp = Vec3r(center.x(), center.y(), _viewpoint.z());
2088   }
2089   else {
2090     vp = Vec3r(_viewpoint);
2091   }
2092   Vec3r u(vp - center);
2093   real raylength = u.norm();
2094   u.normalize();
2095 #if 0
2096   if (_global.debug & G_DEBUG_FREESTYLE) {
2097     cout << "grid origin " << iGrid->getOrigin().x() << "," << iGrid->getOrigin().y() << ","
2098          << iGrid->getOrigin().z() << endl;
2099     cout << "center " << center.x() << "," << center.y() << "," << center.z() << endl;
2100   }
2101 #endif
2102
2103   iGrid->castRay(center, vp, occluders, timestamp);
2104
2105   WFace *face = NULL;
2106   if (fe->isSmooth()) {
2107     FEdgeSmooth *fes = dynamic_cast<FEdgeSmooth *>(fe);
2108     face = (WFace *)fes->face();
2109   }
2110   vector<WVertex *> faceVertices;
2111   WVertex::incoming_edge_iterator ie;
2112
2113   WFace *oface;
2114   bool skipFace;
2115   OccludersSet::iterator p, pend;
2116   if (face) {
2117     face->RetrieveVertexList(faceVertices);
2118   }
2119
2120   for (p = occluders.begin(), pend = occluders.end(); p != pend; p++) {
2121     // If we're dealing with an exact silhouette, check whether we must take care of this occluder
2122     // of not. (Indeed, we don't consider the occluders that share at least one vertex with the
2123     // face containing this edge).
2124     //-----------
2125     oface = (WFace *)(*p)->userdata;
2126 #if LOGGING
2127     if (_global.debug & G_DEBUG_FREESTYLE) {
2128       cout << "\t\tEvaluating intersection for occluder " << ((*p)->getVertices())[0]
2129            << ((*p)->getVertices())[1] << ((*p)->getVertices())[2] << endl
2130            << "\t\t\tand ray " << vp << " * " << u << " (center " << center << ")" << endl;
2131     }
2132 #endif
2133     Vec3r v1(((*p)->getVertices())[0]);
2134     Vec3r normal((*p)->getNormal());
2135     real d = -(v1 * normal);
2136     real t, t_u, t_v;
2137
2138 #if LOGGING
2139     if (_global.debug & G_DEBUG_FREESTYLE) {
2140       cout << "\t\tp:  " << ((*p)->getVertices())[0] << ((*p)->getVertices())[1]
2141            << ((*p)->getVertices())[2] << ", norm: " << (*p)->getNormal() << endl;
2142     }
2143 #endif
2144
2145     if (face) {
2146 #if LOGGING
2147       if (_global.debug & G_DEBUG_FREESTYLE) {
2148         cout << "\t\tDetermining face adjacency...";
2149       }
2150 #endif
2151       skipFace = false;
2152
2153       if (face == oface) {
2154 #if LOGGING
2155         if (_global.debug & G_DEBUG_FREESTYLE) {
2156           cout << "  Rejecting occluder for face concurrency." << endl;
2157         }
2158 #endif
2159         continue;
2160       }
2161
2162       for (vector<WVertex *>::iterator fv = faceVertices.begin(), fvend = faceVertices.end();
2163            fv != fvend;
2164            ++fv) {
2165         if ((*fv)->isBoundary()) {
2166           continue;
2167         }
2168
2169         WVertex::incoming_edge_iterator iebegin = (*fv)->incoming_edges_begin();
2170         WVertex::incoming_edge_iterator ieend = (*fv)->incoming_edges_end();
2171         for (ie = iebegin; ie != ieend; ++ie) {
2172           if ((*ie) == 0) {
2173             continue;
2174           }
2175
2176           WFace *sface = (*ie)->GetbFace();
2177           // WFace *sfacea = (*ie)->GetaFace();
2178           // if ((sface == oface) || (sfacea == oface)) {
2179           if (sface == oface) {
2180             skipFace = true;
2181             break;
2182           }
2183         }
2184         if (skipFace) {
2185           break;
2186         }
2187       }
2188       if (skipFace) {
2189 #if LOGGING
2190         if (_global.debug & G_DEBUG_FREESTYLE) {
2191           cout << "  Rejecting occluder for face adjacency." << endl;
2192         }
2193 #endif
2194         continue;
2195       }
2196     }
2197     else {
2198       // check whether the edge and the polygon plane are coincident:
2199       //-------------------------------------------------------------
2200       // first let us compute the plane equation.
2201
2202       if (GeomUtils::COINCIDENT ==
2203           GeomUtils::intersectRayPlane(origin, edgeDir, normal, d, t, epsilon)) {
2204 #if LOGGING
2205         if (_global.debug & G_DEBUG_FREESTYLE) {
2206           cout << "\t\tRejecting occluder for target coincidence." << endl;
2207         }
2208 #endif
2209         continue;
2210       }
2211     }
2212
2213     if ((*p)->rayIntersect(center, u, t, t_u, t_v)) {
2214 #if LOGGING
2215       if (_global.debug & G_DEBUG_FREESTYLE) {
2216         cout << "\t\tRay " << vp << " * " << u << " intersects at time " << t << " (raylength is "
2217              << raylength << ")" << endl;
2218         cout << "\t\t(u * normal) == " << (u * normal) << " for normal " << normal << endl;
2219       }
2220 #endif
2221       if (fabs(u * normal) > 0.0001) {
2222         if ((t > 0.0) && (t < raylength)) {
2223 #if LOGGING
2224           if (_global.debug & G_DEBUG_FREESTYLE) {
2225             cout << "\t\tIs occluder" << endl;
2226           }
2227 #endif
2228           WFace *f = (WFace *)((*p)->userdata);
2229           ViewShape *vshape = _ViewMap->viewShape(f->GetVertex(0)->shape()->GetId());
2230           oOccluders.insert(vshape);
2231           ++qi;
2232           if (!_EnableQI) {
2233             break;
2234           }
2235         }
2236       }
2237     }
2238   }
2239
2240   // Find occludee
2241   FindOccludee(fe, iGrid, epsilon, oaPolygon, timestamp, u, center, origin, edgeDir, faceVertices);
2242
2243   return qi;
2244 }
2245
2246 void ViewMapBuilder::ComputeIntersections(ViewMap *ioViewMap,
2247                                           intersection_algo iAlgo,
2248                                           real epsilon)
2249 {
2250   switch (iAlgo) {
2251     case sweep_line:
2252       ComputeSweepLineIntersections(ioViewMap, epsilon);
2253       break;
2254     default:
2255       break;
2256   }
2257 #if 0
2258   if (_global.debug & G_DEBUG_FREESTYLE) {
2259     ViewMap::viewvertices_container &vvertices = ioViewMap->ViewVertices();
2260     for (ViewMap::viewvertices_container::iterator vv = vvertices.begin(), vvend = vvertices.end();
2261          vv != vvend;
2262          ++vv) {
2263       if ((*vv)->getNature() == Nature::T_VERTEX) {
2264         TVertex *tvertex = (TVertex *)(*vv);
2265         cout << "TVertex " << tvertex->getId() << " has :" << endl;
2266         cout << "FrontEdgeA: " << tvertex->frontEdgeA().first << endl;
2267         cout << "FrontEdgeB: " << tvertex->frontEdgeB().first << endl;
2268         cout << "BackEdgeA: " << tvertex->backEdgeA().first << endl;
2269         cout << "BackEdgeB: " << tvertex->backEdgeB().first << endl << endl;
2270       }
2271     }
2272   }
2273 #endif
2274 }
2275
2276 struct less_SVertex2D : public binary_function<SVertex *, SVertex *, bool> {
2277   real epsilon;
2278
2279   less_SVertex2D(real eps) : binary_function<SVertex *, SVertex *, bool>()
2280   {
2281     epsilon = eps;
2282   }
2283
2284   bool operator()(SVertex *x, SVertex *y)
2285   {
2286     Vec3r A = x->point2D();
2287     Vec3r B = y->point2D();
2288     for (unsigned int i = 0; i < 3; i++) {
2289       if ((fabs(A[i] - B[i])) < epsilon) {
2290         continue;
2291       }
2292       if (A[i] < B[i]) {
2293         return true;
2294       }
2295       if (A[i] > B[i]) {
2296         return false;
2297       }
2298     }
2299     return false;
2300   }
2301 };
2302
2303 typedef Segment<FEdge *, Vec3r> segment;
2304 typedef Intersection<segment> intersection;
2305
2306 struct less_Intersection : public binary_function<intersection *, intersection *, bool> {
2307   segment *edge;
2308
2309   less_Intersection(segment *iEdge) : binary_function<intersection *, intersection *, bool>()
2310   {
2311     edge = iEdge;
2312   }
2313
2314   bool operator()(intersection *x, intersection *y)
2315   {
2316     real tx = x->getParameter(edge);
2317     real ty = y->getParameter(edge);
2318     if (tx > ty) {
2319       return true;
2320     }
2321     return false;
2322   }
2323 };
2324
2325 struct silhouette_binary_rule : public binary_rule<segment, segment> {
2326   silhouette_binary_rule() : binary_rule<segment, segment>()
2327   {
2328   }
2329
2330   virtual bool operator()(segment &s1, segment &s2)
2331   {
2332     FEdge *f1 = s1.edge();
2333     FEdge *f2 = s2.edge();
2334
2335     if ((!(((f1)->getNature() & Nature::SILHOUETTE) || ((f1)->getNature() & Nature::BORDER))) &&
2336         (!(((f2)->getNature() & Nature::SILHOUETTE) || ((f2)->getNature() & Nature::BORDER)))) {
2337       return false;
2338     }
2339
2340     return true;
2341   }
2342 };
2343
2344 void ViewMapBuilder::ComputeSweepLineIntersections(ViewMap *ioViewMap, real epsilon)
2345 {
2346   vector<SVertex *> &svertices = ioViewMap->SVertices();
2347   bool progressBarDisplay = false;
2348   unsigned sVerticesSize = svertices.size();
2349   unsigned fEdgesSize = ioViewMap->FEdges().size();
2350 #if 0
2351   if (_global.debug & G_DEBUG_FREESTYLE) {
2352     ViewMap::fedges_container &fedges = ioViewMap->FEdges();
2353     for (ViewMap::fedges_container::const_iterator f = fedges.begin(), end = fedges.end();
2354          f != end;
2355          ++f) {
2356       cout << (*f)->aMaterialIndex() << "-" << (*f)->bMaterialIndex() << endl;
2357     }
2358   }
2359 #endif
2360   unsigned progressBarStep = 0;
2361
2362   if (_pProgressBar != NULL && fEdgesSize > gProgressBarMinSize) {
2363     unsigned progressBarSteps = min(gProgressBarMaxSteps, sVerticesSize);
2364     progressBarStep = sVerticesSize / progressBarSteps;
2365     _pProgressBar->reset();
2366     _pProgressBar->setLabelText("Computing Sweep Line Intersections");
2367     _pProgressBar->setTotalSteps(progressBarSteps);
2368     _pProgressBar->setProgress(0);
2369     progressBarDisplay = true;
2370   }
2371
2372   unsigned counter = progressBarStep;
2373
2374   sort(svertices.begin(), svertices.end(), less_SVertex2D(epsilon));
2375
2376   SweepLine<FEdge *, Vec3r> SL;
2377
2378   vector<FEdge *> &ioEdges = ioViewMap->FEdges();
2379
2380   vector<segment *> segments;
2381
2382   vector<FEdge *>::iterator fe, fend;
2383
2384   for (fe = ioEdges.begin(), fend = ioEdges.end(); fe != fend; fe++) {
2385     segment *s = new segment((*fe), (*fe)->vertexA()->point2D(), (*fe)->vertexB()->point2D());
2386     (*fe)->userdata = s;
2387     segments.push_back(s);
2388   }
2389
2390   vector<segment *> vsegments;
2391   for (vector<SVertex *>::iterator sv = svertices.begin(), svend = svertices.end(); sv != svend;
2392        sv++) {
2393     if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
2394       break;
2395     }
2396
2397     const vector<FEdge *> &vedges = (*sv)->fedges();
2398
2399     for (vector<FEdge *>::const_iterator sve = vedges.begin(), sveend = vedges.end();
2400          sve != sveend;
2401          sve++) {
2402       vsegments.push_back((segment *)((*sve)->userdata));
2403     }
2404
2405     Vec3r evt((*sv)->point2D());
2406     silhouette_binary_rule sbr;
2407     SL.process(evt, vsegments, sbr, epsilon);
2408
2409     if (progressBarDisplay) {
2410       counter--;
2411       if (counter <= 0) {
2412         counter = progressBarStep;
2413         _pProgressBar->setProgress(_pProgressBar->getProgress() + 1);
2414       }
2415     }
2416     vsegments.clear();
2417   }
2418
2419   if (_pRenderMonitor && _pRenderMonitor->testBreak()) {
2420     // delete segments
2421     if (!segments.empty()) {
2422       vector<segment *>::iterator s, send;
2423       for (s = segments.begin(), send = segments.end(); s != send; s++) {
2424         delete *s;
2425       }
2426     }
2427     return;
2428   }
2429
2430   // reset userdata:
2431   for (fe = ioEdges.begin(), fend = ioEdges.end(); fe != fend; fe++) {
2432     (*fe)->userdata = NULL;
2433   }
2434
2435   // list containing the new edges resulting from splitting operations.
2436   vector<FEdge *> newEdges;
2437
2438   // retrieve the intersected edges:
2439   vector<segment *> &iedges = SL.intersectedEdges();
2440   // retrieve the intersections:
2441   vector<intersection *> &intersections = SL.intersections();
2442
2443   int id = 0;
2444   // create a view vertex for each intersection and linked this one with the intersection object
2445   vector<intersection *>::iterator i, iend;
2446   for (i = intersections.begin(), iend = intersections.end(); i != iend; i++) {
2447     FEdge *fA = (*i)->EdgeA->edge();
2448     FEdge *fB = (*i)->EdgeB->edge();
2449
2450     Vec3r A1 = fA->vertexA()->point3D();
2451     Vec3r A2 = fA->vertexB()->point3D();
2452     Vec3r B1 = fB->vertexA()->point3D();
2453     Vec3r B2 = fB->vertexB()->point3D();
2454
2455     Vec3r a1 = fA->vertexA()->point2D();
2456     Vec3r a2 = fA->vertexB()->point2D();
2457     Vec3r b1 = fB->vertexA()->point2D();
2458     Vec3r b2 = fB->vertexB()->point2D();
2459
2460     real ta = (*i)->tA;
2461     real tb = (*i)->tB;
2462
2463     if ((ta < -epsilon) || (ta > 1 + epsilon)) {
2464       cerr << "Warning: 2D intersection out of range for edge " << fA->vertexA()->getId() << " - "
2465            << fA->vertexB()->getId() << endl;
2466     }
2467
2468     if ((tb < -epsilon) || (tb > 1 + epsilon)) {
2469       cerr << "Warning: 2D intersection out of range for edge " << fB->vertexA()->getId() << " - "
2470            << fB->vertexB()->getId() << endl;
2471     }
2472
2473     real Ta = SilhouetteGeomEngine::ImageToWorldParameter(fA, ta);
2474     real Tb = SilhouetteGeomEngine::ImageToWorldParameter(fB, tb);
2475
2476     if ((Ta < -epsilon) || (Ta > 1 + epsilon)) {
2477       cerr << "Warning: 3D intersection out of range for edge " << fA->vertexA()->getId() << " - "
2478            << fA->vertexB()->getId() << endl;
2479     }
2480
2481     if ((Tb < -epsilon) || (Tb > 1 + epsilon)) {
2482       cerr << "Warning: 3D intersection out of range for edge " << fB->vertexA()->getId() << " - "
2483            << fB->vertexB()->getId() << endl;
2484     }
2485
2486 #if 0
2487     if (_global.debug & G_DEBUG_FREESTYLE) {
2488       if ((Ta < -epsilon) || (Ta > 1 + epsilon) || (Tb < -epsilon) || (Tb > 1 + epsilon)) {
2489         printf("ta %.12e\n", ta);
2490         printf("tb %.12e\n", tb);
2491         printf("a1 %e, %e -- a2 %e, %e\n", a1[0], a1[1], a2[0], a2[1]);
2492         printf("b1 %e, %e -- b2 %e, %e\n", b1[0], b1[1], b2[0], b2[1]);
2493         //printf("line([%e, %e], [%e, %e]);\n", a1[0], a2[0], a1[1], a2[1]);
2494         //printf("line([%e, %e], [%e, %e]);\n", b1[0], b2[0], b1[1], b2[1]);
2495         if ((Ta < -epsilon) || (Ta > 1 + epsilon)) {
2496           printf("Ta %.12e\n", Ta);
2497         }
2498         if ((Tb < -epsilon) || (Tb > 1 + epsilon)) {
2499           printf("Tb %.12e\n", Tb);
2500         }
2501         printf("A1 %e, %e, %e -- A2 %e, %e, %e\n", A1[0], A1[1], A1[2], A2[0], A2[1], A2[2]);
2502         printf("B1 %e, %e, %e -- B2 %e, %e, %e\n", B1[0], B1[1], B1[2], B2[0], B2[1], B2[2]);
2503       }
2504     }
2505 #endif
2506
2507     TVertex *tvertex = ioViewMap->CreateTVertex(Vec3r(A1 + Ta * (A2 - A1)),
2508                                                 Vec3r(a1 + ta * (a2 - a1)),
2509                                                 fA,
2510                                                 Vec3r(B1 + Tb * (B2 - B1)),
2511                                                 Vec3r(b1 + tb * (b2 - b1)),
2512                                                 fB,
2513                                                 id);
2514
2515     (*i)->userdata = tvertex;
2516     ++id;
2517   }
2518
2519   progressBarStep = 0;
2520
2521   if (progressBarDisplay) {
2522     unsigned iEdgesSize = iedges.size();
2523     unsigned progressBarSteps = min(gProgressBarMaxSteps, iEdgesSize);
2524     progressBarStep = iEdgesSize / progressBarSteps;
2525     _pProgressBar->reset();
2526     _pProgressBar->setLabelText("Splitting intersected edges");
2527     _pProgressBar->setTotalSteps(progressBarSteps);
2528     _pProgressBar->setProgress(0);
2529   }
2530
2531   counter = progressBarStep;
2532
2533   vector<TVertex *> edgeVVertices;
2534   vector<ViewEdge *> newVEdges;
2535   vector<segment *>::iterator s, send;
2536   for (s = iedges.begin(), send = iedges.end(); s != send; s++) {
2537     edgeVVertices.clear();
2538     newEdges.clear();
2539     newVEdges.clear();
2540
2541     FEdge *fedge = (*s)->edge();
2542     ViewEdge *vEdge = fedge->viewedge();
2543     ViewShape *shape = vEdge->viewShape();
2544
2545     vector<intersection *> &eIntersections = (*s)->intersections();
2546     // we first need to sort these intersections from farther to closer to A
2547     sort(eIntersections.begin(), eIntersections.end(), less_Intersection(*s));
2548     for (i = eIntersections.begin(), iend = eIntersections.end(); i != iend; i++) {
2549       edgeVVertices.push_back((TVertex *)(*i)->userdata);
2550     }
2551
2552     shape->SplitEdge(fedge, edgeVVertices, ioViewMap->FEdges(), ioViewMap->ViewEdges());
2553
2554     if (progressBarDisplay) {
2555       counter--;
2556       if (counter <= 0) {
2557         counter = progressBarStep;
2558         _pProgressBar->setProgress(_pProgressBar->getProgress() + 1);
2559       }
2560     }
2561   }
2562
2563   // reset userdata:
2564   for (fe = ioEdges.begin(), fend = ioEdges.end(); fe != fend; fe++) {
2565     (*fe)->userdata = NULL;
2566   }
2567
2568   // delete segments
2569   if (!segments.empty()) {
2570     for (s = segments.begin(), send = segments.end(); s != send; s++) {
2571       delete *s;
2572     }
2573   }
2574 }
2575
2576 } /* namespace Freestyle */