diff --git a/Code/FBSD/BlurredSegment/blurredsegment.cpp b/Code/FBSD/BlurredSegment/blurredsegment.cpp
index 205f4cf19178cedcfad70b56fa1e7aa9b4a9e0a7..29c869d476c370a6e8bddd4415d40357f01fb1f3 100755
--- a/Code/FBSD/BlurredSegment/blurredsegment.cpp
+++ b/Code/FBSD/BlurredSegment/blurredsegment.cpp
@@ -51,6 +51,14 @@ vector<Pt2i> *BlurredSegment::getAllLeft () const
 }
 
 
+DigitalStraightSegment *BlurredSegment::holdSegment ()
+{
+  DigitalStraightSegment *tmp = dss;
+  dss = NULL;
+  return tmp;
+}
+
+
 int BlurredSegment::size () const
 {
   return (plist->size ());
@@ -117,7 +125,7 @@ vector <vector <Pt2i> > BlurredSegment::connectedComponents () const
     vector <Pt2i> cc;
     bool started = false;
     vector <Pt2i>::const_iterator it = pts.begin ();
-    Pt2i pix = *it++;
+    Pt2i pix (*it++);
     while (it != pts.end ())
     {
       if (it->isConnectedTo (pix))
@@ -138,7 +146,7 @@ vector <vector <Pt2i> > BlurredSegment::connectedComponents () const
           started = false;;
         }
       }
-      pix = *it++;
+      pix.set (*it++);
     }
   }
   return ccs;
@@ -153,7 +161,7 @@ int BlurredSegment::countOfConnectedPoints () const
   {
     bool started = false;
     vector <Pt2i>::const_iterator it = pts.begin ();
-    Pt2i pix = *it++;
+    Pt2i pix (*it++);
     while (it != pts.end ())
     {
       if (it->isConnectedTo (pix))
@@ -169,7 +177,7 @@ int BlurredSegment::countOfConnectedPoints () const
       {
         if (started) started = false;;
       }
-      pix = *it++;
+      pix.set (*it++);
     }
   }
   return count;
@@ -184,7 +192,7 @@ int BlurredSegment::countOfConnectedComponents () const
   {
     bool started = false;
     vector <Pt2i>::const_iterator it = pts.begin ();
-    Pt2i pix = *it++;
+    Pt2i pix (*it++);
     while (it != pts.end ())
     {
       if (it->isConnectedTo (pix))
@@ -199,7 +207,7 @@ int BlurredSegment::countOfConnectedComponents () const
       {
         if (started) started = false;;
       }
-      pix = *it++;
+      pix.set (*it++);
     }
   }
   return count;
@@ -214,7 +222,7 @@ int BlurredSegment::countOfConnectedPoints (int min) const
   {
     int cpt = 1;
     vector <Pt2i>::const_iterator it = pts.begin ();
-    Pt2i pix = *it++;
+    Pt2i pix (*it++);
     while (it != pts.end ())
     {
       if (it->isConnectedTo (pix))
@@ -223,7 +231,7 @@ int BlurredSegment::countOfConnectedPoints (int min) const
         else if (cpt > min) count ++;
       }
       else cpt = 1;
-      pix = *it++;
+      pix.set (*it++);
     }
   }
   return count;
@@ -238,7 +246,7 @@ int BlurredSegment::countOfConnectedComponents (int min) const
   {
     int cpt = 1;
     vector <Pt2i>::const_iterator it = pts.begin ();
-    Pt2i pix = *it++;
+    Pt2i pix (*it++);
     while (it != pts.end ())
     {
       if (it->isConnectedTo (pix))
@@ -246,7 +254,7 @@ int BlurredSegment::countOfConnectedComponents (int min) const
         if (++cpt == min) count ++;
       }
       else cpt = 1;;
-      pix = *it++;
+      pix.set (*it++);
     }
   }
   return count;
@@ -271,7 +279,7 @@ BlurredSegment::getConnectedComponents () const
       {
         lres.push_back (pix);
         compose = bit->isConnectedTo (pix);
-        if (compose) pix = *bit++;
+        if (compose) pix.set (*bit++);
       }
       while (compose && bit != eit);
       if (compose) lres.push_back (pix);
diff --git a/Code/FBSD/BlurredSegment/blurredsegment.h b/Code/FBSD/BlurredSegment/blurredsegment.h
index cf65f490a3c5ba36eba0680c4ea66b88fd2ee3a3..ffb5c5c364ff2826f33231dca072f2a21a15aad6 100755
--- a/Code/FBSD/BlurredSegment/blurredsegment.h
+++ b/Code/FBSD/BlurredSegment/blurredsegment.h
@@ -55,6 +55,11 @@ public:
    */
   inline DigitalStraightSegment *getSegment () { return dss; }
 
+  /**
+   * \brief Returns a pointer to a copy of the enclosing DSS.
+   */
+  DigitalStraightSegment *holdSegment ();
+
   /**
    * \brief Returns the count of points of the blurred segment.
    */
@@ -189,6 +194,13 @@ public:
    */
   vector <vector <Pt2i> > getConnectedComponents () const;
 
+  /**
+   * \brief Checks if given point is one of the three antipodal points.
+   * @param pt Given point.
+   */
+  inline bool isAntipodal (Pt2i pt) const {
+    return (pt.equals (laps) || pt.equals (lape) || pt.equals (lapv)); }
+
 
 protected:
 
diff --git a/Code/FBSD/BlurredSegment/bsdetector.cpp b/Code/FBSD/BlurredSegment/bsdetector.cpp
index 60f6d5687906d974dff10997d7c5f45f2da45a30..f1cfa36a2fbdeb7e0a5507aeb34b9e623eb050e0 100755
--- a/Code/FBSD/BlurredSegment/bsdetector.cpp
+++ b/Code/FBSD/BlurredSegment/bsdetector.cpp
@@ -93,6 +93,7 @@ BSDetector::~BSDetector ()
   if (staticDetOn) delete bstStatic;
   delete bst1;
   delete bst2;
+  
   if (lsf1 != NULL) delete lsf1;
   if (lsf2 != NULL) delete lsf2;
   if (bsini != NULL) delete bsini;
@@ -102,6 +103,34 @@ BSDetector::~BSDetector ()
 }
 
 
+void BSDetector::clearAll ()
+{
+  if (lsf1 != NULL)
+  {
+    delete lsf1;
+    lsf1 = NULL;
+  }
+  if (lsf2 != NULL)
+  {
+    delete lsf2;
+    lsf2 = NULL;
+  }
+  if (bsini != NULL)
+  {
+    delete bsini;
+    bsini = NULL;
+  }
+  if (bsf != NULL)
+  {
+    delete bsf;
+    bsf = NULL;
+  }
+  vector <BlurredSegment *>::iterator it = mbsf.begin ();
+  while (it != mbsf.end ()) delete (*it++);
+  mbsf.clear ();
+}
+
+
 void BSDetector::setGradientMap (VMap *data)
 {
   gMap = data;
@@ -230,7 +259,7 @@ bool BSDetector::runMultiDetection (const Pt2i &p1, const Pt2i &p2)
 {
   vector<Pt2i> pts;
   p1.draw (pts, p2);
-  int locmax[pts.size ()];
+  int *locmax = new int[pts.size ()];
   int nlm = gMap->localMax (locmax, pts);
   bool isnext = true;
   for (int i = 0; isnext && i < nlm; i++)
@@ -313,8 +342,8 @@ int BSDetector::detect (const Pt2i &p1, const Pt2i &p2,
       if (detw < PRELIM_MIN_HALF_WIDTH) detw = PRELIM_MIN_HALF_WIDTH;
       int dx = (int) ((v0.y () * detw) / l);
       int dy = (int) (- (v0.x () * detw) / l);
-      inip1 = Pt2i (pcentral.x () + dx, pcentral.y () + dy);
-      inip2 = Pt2i (pcentral.x () - dx, pcentral.y () - dy);
+      inip1.set (pcentral.x () + dx, pcentral.y () + dy);
+      inip2.set (pcentral.x () - dx, pcentral.y () - dy);
       iniwidth = 0;
     }
   }
@@ -377,7 +406,7 @@ int BSDetector::detect (const Pt2i &p1, const Pt2i &p2,
   // Scan recentering and fitting
   //-----------------------------
   if (recenteringOn)
-    pCenter = bsini->getSegment()->centerOfIntersection (inip1, inip2);
+    pCenter.set (bsini->getSegment()->centerOfIntersection (inip1, inip2));
 
   // Finer detection based on gradient maxima with orientation constraint
   //---------------------------------------------------------------------
@@ -500,7 +529,7 @@ int BSDetector::staticDetect (const Pt2i &p1, const Pt2i &p2,
   // Scan recentering and fitting
   //-----------------------------
   if (recenteringOn)
-    pCenter = bsini->getSegment()->centerOfIntersection (inip1, inip2);
+    pCenter.set (bsini->getSegment()->centerOfIntersection (inip1, inip2));
 
   // Finer detection based on gradient maxima with orientation constraint
   //---------------------------------------------------------------------
@@ -511,9 +540,9 @@ int BSDetector::staticDetect (const Pt2i &p1, const Pt2i &p2,
 
   // Scan recentering and fitting
   //-----------------------------
-  pCenter = bsini->getCenter ();
+  pCenter.set (bsini->getCenter ());
   if (recenteringOn)
-    pCenter = bsf->getSegment()->centerOfIntersection (inip1, inip2);
+    pCenter.set (bsf->getSegment()->centerOfIntersection (inip1, inip2));
 
   // Third detection based on gradient maxima with orientation constraint
   //---------------------------------------------------------------------
@@ -584,6 +613,66 @@ BlurredSegment *BSDetector::getBlurredSegment (int step) const
 }
 
 
+int BSDetector::getDigitalStraightSegments (
+                        vector<DigitalStraightSegment *> &dss) const
+{
+  if (mbsf.empty ())
+  {
+    if (bsf == NULL) return 0;
+    DigitalStraightSegment *ds = bsf->getSegment ();
+    if (ds == NULL) return (0);
+    dss.push_back (ds);
+    return (1);
+  }
+  int nb = 0;
+  vector<BlurredSegment *>::const_iterator it = mbsf.begin ();
+  while (it != mbsf.end ())
+  {
+    BlurredSegment *bs = (*it++);
+    if (bs != NULL)
+    {
+      DigitalStraightSegment *ds = bs->getSegment ();
+      if (ds != NULL)
+      {
+        dss.push_back (ds);
+        nb ++;
+      }
+    }
+  }
+  return nb;
+}
+
+
+int BSDetector::copyDigitalStraightSegments (
+                             vector<DigitalStraightSegment> &dss) const
+{
+  if (mbsf.empty ())
+  {
+    if (bsf == NULL) return 0;
+    DigitalStraightSegment *ds = bsf->getSegment ();
+    if (ds == NULL) return 0;
+    dss.push_back (DigitalStraightSegment (ds));
+    return 1;
+  }
+  int nb = 0;
+  vector<BlurredSegment *>::const_iterator it = mbsf.begin ();
+  while (it != mbsf.end ())
+  {
+    BlurredSegment *bs = (*it++);
+    if (bs != NULL)
+    {
+      DigitalStraightSegment *ds = bs->getSegment ();
+      if (ds != NULL)
+      {
+        dss.push_back (DigitalStraightSegment (ds));
+        nb ++;
+      }
+    }
+  }
+  return nb;
+}
+
+
 void BSDetector::incMaxDetections (bool dir)
 {
   maxtrials = maxtrials + (dir ? -1 : 1);
diff --git a/Code/FBSD/BlurredSegment/bsdetector.h b/Code/FBSD/BlurredSegment/bsdetector.h
index 4340fb1991eb2a6fff03112b6357e6793480aa9a..0a53f3e8a060d5426fec4f63c9f9af3490dfb24a 100755
--- a/Code/FBSD/BlurredSegment/bsdetector.h
+++ b/Code/FBSD/BlurredSegment/bsdetector.h
@@ -71,6 +71,11 @@ public:
    */
   ~BSDetector ();
 
+  /**
+   * \brief Delete all detected structures.
+   */
+  void clearAll ();
+
   /**
    * \brief Returns the minimal vertical or horizontal width.
    */
@@ -167,6 +172,18 @@ public:
    */
   inline void preserveFormerBlurredSegments () { mbsf.clear (); }
 
+  /**
+   * \brief Returns the digital straight segments in provided list.
+   * @param dss List to fill in with digital straight segments.
+   */
+  int getDigitalStraightSegments (vector<DigitalStraightSegment *> &dss) const;
+
+  /**
+   * \brief Returns a copy of the digital straight segments in provided list.
+   * @param dss List to fill in with the copy of digital straight segments.
+   */
+  int copyDigitalStraightSegments (vector<DigitalStraightSegment> &dss) const;
+
   /**
    * \brief Returns the assigned maximal thickness to detector.
    */
diff --git a/Code/FBSD/BlurredSegment/bsfilter.cpp b/Code/FBSD/BlurredSegment/bsfilter.cpp
index ba081d3f9e06e2ff8ec6dd5d5d937e8120ca7b75..94695c86e6dbd659018ebaa579e6f7a3a738e850 100755
--- a/Code/FBSD/BlurredSegment/bsfilter.cpp
+++ b/Code/FBSD/BlurredSegment/bsfilter.cpp
@@ -1,5 +1,5 @@
 #include "bsfilter.h"
-#include "blurredsegmentproto.h"
+#include "bsproto.h"
 
 
 
@@ -37,7 +37,7 @@ BlurredSegment *BSFilter::filter (BlurredSegment *bs)
   int mw = 1;
   AbsRat sw = bs->minimalWidth ();
   if (sw.denominator () != 0) mw += sw.numerator () / sw.denominator ();
-  BlurredSegmentProto resbs (mw, bs->getCenter (), leftIn, rightIn);
+  BSProto resbs (mw, bs->getCenter (), leftIn, rightIn);
   bsFinalSize = resbs.size ();
   return (resbs.endOfBirth ());
 }
diff --git a/Code/FBSD/BlurredSegment/bsproto.cpp b/Code/FBSD/BlurredSegment/bsproto.cpp
new file mode 100755
index 0000000000000000000000000000000000000000..6581e80d33822d76324fc529b2fc55ac66dddf6d
--- /dev/null
+++ b/Code/FBSD/BlurredSegment/bsproto.cpp
@@ -0,0 +1,318 @@
+#include "bsproto.h"
+
+
+BSProto::BSProto (int maxWidth, Pt2i pix)
+{
+  this->maxWidth.set (maxWidth);
+  plist = new BiPtList (pix);
+  leftOK = false;
+  rightOK = false;
+  bsFlat = false;
+  bsOK = false;
+  convexhull = NULL;
+  chChanged = false;
+  dss = NULL;
+}
+
+
+BSProto::BSProto (int maxWidth, Pt2i center,
+                  const vector<Pt2i> &leftPts, const vector<Pt2i> &rightPts)
+{
+  this->maxWidth.set (maxWidth);
+  plist = new BiPtList (center);
+  leftOK = false;
+  rightOK = false;
+  bsFlat = false;
+  bsOK = false;
+  convexhull = NULL;
+  chChanged = false;
+  dss = NULL;
+
+  vector<Pt2i>::const_iterator itr = rightPts.begin ();
+  vector<Pt2i>::const_iterator itl = leftPts.begin ();
+  bool scanningRight = true;
+  bool scanningLeft = true;
+  while (scanningRight || scanningLeft)
+  {
+    if (scanningRight)
+    {
+      if (itr == rightPts.end ()) scanningRight = false;
+      else
+        if (! addRight (*itr++)) scanningRight = false;
+    }
+    if (scanningLeft)
+    {
+      if (itl == leftPts.end ()) scanningLeft = false;
+      else
+        if (! addLeft (*itl++)) scanningLeft = false;
+    }
+  }
+}
+
+
+BSProto::~BSProto ()
+{
+  if (convexhull != NULL) delete convexhull;
+}
+
+
+AbsRat BSProto::strictThickness () const
+{
+  return (convexhull != NULL ? convexhull->thickness () : AbsRat (0, 1));
+}
+
+
+AbsRat BSProto::digitalThickness () const
+{
+  if (bsOK)
+  {
+    Pt2i s, e, v;
+    convexhull->antipodalEdgeAndVertex (s, e, v);
+    DigitalStraightLine l (s, e, v);
+    return (AbsRat (l.width (), l.period ()));
+  }
+  return (AbsRat (1, 1));
+}
+
+
+DigitalStraightLine *BSProto::getLine () const
+{
+  if (bsOK)
+  {
+    Pt2i s, e, v;
+    convexhull->antipodalEdgeAndVertex (s, e, v);
+    return (new DigitalStraightLine (s, e, v));
+  }
+  if (bsFlat || leftOK || rightOK)
+    return (new DigitalStraightLine (getLastLeft (), getLastRight (),
+                                     DigitalStraightLine::DSL_THIN));
+  return (NULL); // No line if only one point
+}
+
+
+bool BSProto::addLeftSorted (Pt2i pix)
+{
+  if (pix.equals (plist->frontPoint ()))
+  { 
+    plist->addFront (pix);
+    chChanged = false;
+    return true;
+  }
+  return addLeft (pix);
+}
+
+
+bool BSProto::addRightSorted (Pt2i pix)
+{
+  if (pix.equals (plist->backPoint ()))
+  {
+    plist->addBack (pix);
+    chChanged = false;
+    return true;
+  }
+  return addRight (pix);
+}
+
+
+bool BSProto::addLeft (Pt2i pix)
+{
+  if (bsOK)     // convexhull defined
+  {
+    bool res = addPoint (pix, true);
+    return (res);
+  }
+
+  else if (bsFlat)    // leftOK && rightOK
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () != 0)
+    {
+      convexhull = new ConvexHull (pix, plist->frontPoint (),
+                                        plist->backPoint ());
+      bsOK = true;
+    }
+    plist->addFront (pix);
+  }
+
+  else if (leftOK)    // thus ! rightOK
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () == 0) bsFlat = true;
+    else
+    {
+      convexhull = new ConvexHull (pix, plist->frontPoint (),
+                                        plist->backPoint ());
+      bsOK = true;
+    }
+    plist->addFront (pix);
+  }
+
+  else if (rightOK)
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () == 0) bsFlat = true;
+    else
+    {
+      convexhull = new ConvexHull (pix, plist->frontPoint (),
+                                        plist->backPoint ());
+      bsOK = true;
+    }
+    plist->addFront (pix);
+  }
+
+  else   // only the central point is ok
+  {
+    plist->addFront (pix);
+    leftOK = true;
+  }
+
+  chChanged = true;
+  return true;
+}
+
+
+bool BSProto::addRight (Pt2i pix)
+{
+  if (bsOK)     // bs != NULL
+  {
+    bool res = addPoint (pix, false);
+    return (res);
+  }
+
+  else if (bsFlat)    // leftOK && rightOK
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () != 0)
+    {
+      convexhull = new ConvexHull (plist->frontPoint (),
+                                   plist->backPoint (), pix);
+      bsOK = true;
+    }
+    plist->addBack (pix);
+  }
+
+  else if (rightOK)    // thus ! leftOK
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () == 0) bsFlat = true;
+    else
+    {
+      convexhull = new ConvexHull (plist->frontPoint (),
+                                   plist->backPoint (), pix);
+      bsOK = true;
+    }
+    plist->addBack (pix);
+  }
+
+  else if (leftOK)
+  {
+    AbsRat height = plist->heightToEnds (pix);
+    if (height.greaterThan (maxWidth)) return false;
+    if (height.numerator () == 0) bsFlat = true;
+    else
+    {
+      convexhull = new ConvexHull (plist->frontPoint (),
+                                   plist->backPoint (), pix);
+      bsOK = true;
+    }
+    plist->addBack (pix);
+  }
+
+  else   // only the central point is ok
+  {
+    plist->addBack (pix);
+    rightOK = true;
+  }
+
+  chChanged = true;
+  return true;
+}
+
+
+bool BSProto::addPoint (Pt2i p, bool onleft)
+{
+  bool inserted = convexhull->addPointDS (p, onleft);
+  if ((strictThickness ()).greaterThan (maxWidth))
+  {
+    if (inserted) convexhull->restore ();
+    return false;
+  }
+ 
+  if (onleft) plist->addFront (p);
+  else plist->addBack (p);
+  chChanged = true;
+  return true;
+}
+
+
+void BSProto::removeLeft (int n)
+{
+  if (bsOK) plist->removeFront (n);
+}
+
+
+void BSProto::removeRight (int n)
+{
+  if (bsOK) plist->removeBack (n);
+}
+ 
+
+Vr2i BSProto::getSupportVector ()
+{
+  if (bsOK)
+  {
+    Pt2i s, e, v;
+    convexhull->antipodalEdgeAndVertex (s, e, v);
+    return (s.vectorTo (e));
+  }
+  if (bsFlat || leftOK || rightOK)
+    return (getLastLeft().vectorTo (getLastRight ()));
+  return (Vr2i (1, 0));    // hardly better with only one point !
+}
+
+
+BlurredSegment *BSProto::endOfBirth ()
+{
+  DigitalStraightSegment *seg = NULL;
+  if (bsOK)
+  {
+    int xmin, ymin, xmax, ymax;
+    plist->findExtrema (xmin, ymin, xmax, ymax);
+    Pt2i s, e, v;
+    convexhull->antipodalEdgeAndVertex (s, e, v);
+    seg = new DigitalStraightSegment (s, e, v, xmin, ymin, xmax, ymax);
+  }
+  else if (bsFlat || rightOK || leftOK)
+  {
+    Pt2i llast = plist->frontPoint ();
+    Pt2i rlast = plist->backPoint ();
+    if (llast.equals (rlast)) // Strange, should not be flat, rightok or leftok
+    {
+      plist = NULL;
+      return (NULL);
+    }
+    int xmin = llast.x ();
+    if (rlast.x () < llast.x ()) xmin = rlast.x ();
+    int ymin = llast.y ();
+    if (rlast.y () < llast.y ()) ymin = rlast.y ();
+    int xmax = llast.x ();
+    if (rlast.x () > llast.x ()) xmax = rlast.x ();
+    int ymax = llast.y ();
+    if (rlast.y () > llast.y ()) ymax = rlast.y ();
+    seg = new DigitalStraightSegment (llast, rlast,
+                                      DigitalStraightLine::DSL_THIN,
+                                      xmin, ymin, xmax, ymax);
+  }
+  else return (NULL);
+  Pt2i aps (-1, -1), ape (-1, -1), apv (-1, -1);
+  if (convexhull != NULL)
+    convexhull->antipodalEdgeAndVertex (aps, ape, apv);
+  BlurredSegment *bbs = new BlurredSegment (plist, seg, aps, ape, apv);
+  plist = NULL;  // NECESSARY TO AVOID CONTENTS CLEARANCE !!!
+  return (bbs);
+}
diff --git a/Code/FBSD/BlurredSegment/bsproto.h b/Code/FBSD/BlurredSegment/bsproto.h
new file mode 100755
index 0000000000000000000000000000000000000000..15bb7c9cec92be04de0df863b67baffdffb8917c
--- /dev/null
+++ b/Code/FBSD/BlurredSegment/bsproto.h
@@ -0,0 +1,161 @@
+#ifndef BLURRED_SEGMENT_PROTO_H
+#define BLURRED_SEGMENT_PROTO_H
+
+
+#include "blurredsegment.h"
+
+using namespace std;
+
+
+/** 
+ * @class BSProto bsproto.h
+ * \brief A prototype of blurred segment, untill complete specification.
+ * It is mostly based on a evolving list of points, its convex hull and
+ * the successive states of the blurred segment construction.
+ * \author {P. Even}
+ */
+class BSProto : public BlurredSegment
+{
+public:
+
+  /**
+   * Creates a blurred segment prototype.
+   * @param maxWidth Maximal width of the blurred segment to build
+   * @param pix Central point of the blurred segment to build
+   */
+  BSProto (int maxWidth, Pt2i pix);
+
+  /**
+   * Creates a blurred segment prototype with lists of points.
+   * @param maxWidth Maximal width of the blurred segment to build.
+   * @param center Central point of the blurred segment to build.
+   * @param leftPts Points to add on left side
+   * @param rightPts Points to add on right side.
+   */
+  BSProto (int maxWidth, Pt2i center,
+                  const vector<Pt2i> &leftPts, const vector<Pt2i> &rightPts);
+
+  /**
+   * \brief Deletes the blurred segment prototype.
+   */
+  ~BSProto ();
+
+  /**
+   * \brief Checks if the blurred segment has at least two points.
+   */
+  inline bool isLineable () const {
+    return (bsOK || bsFlat || leftOK || rightOK); }
+
+  /**
+   * \brief Checks if the blurred segment is not flat (true BS).
+   */
+  inline bool isNotFlat () const { return (bsOK); }
+
+  /**
+   * \brief Returns the built-in blurred segment strict thickness.
+   * The strict thickness is the distance between bounding lines, ie (nu-1)/p.
+   */
+  AbsRat strictThickness () const;
+
+  /**
+   * \brief Returns the built-in blurred segment digital thickness.
+   * The digital thickness is the width of the digital straight line, ie nu/p.
+   */
+  AbsRat digitalThickness () const;
+
+  /**
+   * \brief Returns the requested max width of the segment.
+   */
+  inline const AbsRat getMaxWidth () const { return maxWidth; }
+
+  /**
+   * \brief Sets the requested max width of the segment.
+   * @param maxwidth The new vaue for the requested max width.
+   */
+  inline void setMaxWidth (const AbsRat &val) { maxWidth.set (val); }
+
+  /**
+   * \brief Returns the underlying digital straight line.
+   */
+  DigitalStraightLine *getLine () const;
+
+  /**
+   * Adds a new sorted point to the left (possibly identical to previous one).
+   * @param pix point to be added.
+   * Returns true if the point is accepted.
+   */
+  bool addLeftSorted (Pt2i pix);
+
+  /**
+   * Adds a new sorted point to the right (possibly identical to previous one).
+   * @param pix point to be added.
+   * Returns true if the point is accepted.
+   */
+  bool addRightSorted (Pt2i pix);
+
+  /**
+   * Adds a new point on the left.
+   * @param pix point to be added.
+   * Returns true if the point is accepted.
+   */
+  bool addLeft (Pt2i pix);
+
+  /**
+   * Adds a new point on the right.
+   * @param pix point to be added.
+   * Returns true if the point is accepted.
+   */
+  bool addRight (Pt2i pix);
+
+  /**
+   * \brief Remove last points on the left side.
+   * @param n Number of points to remove.
+   */
+  void removeLeft (int n);
+
+  /**
+   * \brief Remove last points on the right side.
+   * @param n Number of points to remove.
+   */
+  void removeRight (int n);
+
+  /**
+   * \brief Returns the support vector of the blurred segment.
+   */
+  Vr2i getSupportVector ();
+
+  /**
+   * \brief Returns a stable blurred segment from available information.
+   */
+  BlurredSegment *endOfBirth ();
+
+
+protected:
+
+  /** Maximal width of the blurred segment. */
+  AbsRat maxWidth;
+
+  /** Maintained convex hull of the blurred segment. */
+  ConvexHull *convexhull;
+
+  /** Indicates if the blurred segment is constructed. */
+  bool bsOK;
+  /** Indicates if the points are aligned. */
+  bool bsFlat;
+  /** Indicates if the left point is defined. */
+  bool leftOK;
+  /** Indicates if the right point is defined. */
+  bool rightOK;
+
+  /** Flag indicating if the convex hull changed since last DSS extraction. */
+  bool chChanged;
+
+
+  /**
+   * \brief Submits a new point to extend the blurred segment.
+   * @param p Submitted point.
+   * @param onleft Adding direction (true for LEFT, false for RIGHT).
+   */
+  bool addPoint (Pt2i p, bool onleft);
+};
+#endif
diff --git a/Code/FBSD/BlurredSegment/bstracker.cpp b/Code/FBSD/BlurredSegment/bstracker.cpp
index 63abef529758c97e2385e53528613be9eddf4c63..f0e10da8d4d885ec6f1bd7ec9e7e58f5778d877e 100755
--- a/Code/FBSD/BlurredSegment/bstracker.cpp
+++ b/Code/FBSD/BlurredSegment/bstracker.cpp
@@ -1,5 +1,5 @@
 #include "bstracker.h"
-#include "blurredsegmentproto.h"
+#include "bsproto.h"
 
 
 
@@ -111,7 +111,7 @@ BlurredSegment *BSTracker::fastTrack (int bsMaxWidth,
     pfirst.set (pix.at (candide));
   }
 
-  BlurredSegmentProto bs (bsMaxWidth, pfirst);
+  BSProto bs (bsMaxWidth, pfirst);
   Pt2i lastLeft (pfirst);
   Pt2i lastRight (pfirst);
   
@@ -147,7 +147,7 @@ BlurredSegment *BSTracker::fastTrack (int bsMaxWidth,
         }
         if (added)
         {
-          lastRight = pix.at (candide);
+          lastRight.set (pix.at (candide));
           if (rstop != 0) rstop = 0;
         }
         else if (++rstop > acceptedLacks)
@@ -176,7 +176,7 @@ BlurredSegment *BSTracker::fastTrack (int bsMaxWidth,
         }
         if (added)
         {
-          lastLeft = pix.at (candide);
+          lastLeft.set (pix.at (candide));
           if (lstop != 0) lstop = 0;
         }
         else if (++lstop > acceptedLacks)
@@ -234,7 +234,7 @@ BlurredSegment *BSTracker::fineTrack (int bsMaxWidth,
     return NULL;
   }
 
-  BlurredSegmentProto bs (bsMaxWidth, pix[cand[0]]);
+  BSProto bs (bsMaxWidth, pix[cand[0]]);
 
   // Handles assigned thickness control
   bool atcOn = assignedThicknessControlOn;
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero1.cpp b/Code/FBSD/DirectionalScanner/adaptivescannero1.cpp
index db8f5d34b8622613228f8cb6e2319a756a3ea57c..26dc811d073589f8b57c1f666286971995defc3f 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero1.cpp
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero1.cpp
@@ -113,7 +113,23 @@ AdaptiveScannerO1::AdaptiveScannerO1 (
 }
 
 
-int AdaptiveScannerO1::first (vector<Pt2i> &scan)
+AdaptiveScannerO1::AdaptiveScannerO1 (AdaptiveScannerO1 *ds)
+                 : DirectionalScanner (ds)
+{
+  templ_a = ds->templ_a;
+  templ_b = ds->templ_b;
+  templ_nu = ds->templ_nu;
+  dlc1 = ds->dlc1;
+}
+
+
+DirectionalScanner *AdaptiveScannerO1::getCopy ()
+{
+  return (new AdaptiveScannerO1 (this));
+}
+
+
+int AdaptiveScannerO1::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -138,7 +154,7 @@ int AdaptiveScannerO1::first (vector<Pt2i> &scan)
 int AdaptiveScannerO1::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcx --;
   // Whenever the control line changed
   while (lcy < ymax - 1 && lcx >= xmin && dla * lcx + dlb * lcy > dlc1)
@@ -178,7 +194,7 @@ int AdaptiveScannerO1::nextOnLeft (vector<Pt2i> &scan)
 int AdaptiveScannerO1::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcx ++;
   while (rcy < ymax - 1 && rcx >= xmin && dla * rcx + dlb * rcy > dlc1)
   {
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero1.h b/Code/FBSD/DirectionalScanner/adaptivescannero1.h
index c90c1ebb71d5e69c676019fccc9b8d789de0f316..da2e38bbd8acd5204e9dea053b1217db4c782037 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero1.h
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero1.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class Adaptive dynamicalscannero1.h
  * \brief Adaptive directional scanner for the 1st octant.
- * \author {P. Even}
  */
 class AdaptiveScannerO1 : public DirectionalScanner
 {
@@ -88,13 +87,19 @@ public:
                      int nbs, bool *steps,
                      int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -132,6 +137,15 @@ protected :
 
   /** Shift coefficient of the discrete lower support line */
   int dlc1;
+
+
+  /**
+   * @fn AdaptiveScannerO1(AdaptiveScannerO1 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  AdaptiveScannerO1 (AdaptiveScannerO1 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero2.cpp b/Code/FBSD/DirectionalScanner/adaptivescannero2.cpp
index 839875e0434a487e4942cf40cec7957d44373b9f..559c5d4bbab253eb86acbea4395ca2810f227f15 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero2.cpp
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero2.cpp
@@ -113,7 +113,23 @@ AdaptiveScannerO2::AdaptiveScannerO2 (
 }
 
 
-int AdaptiveScannerO2::first (vector<Pt2i> &scan)
+AdaptiveScannerO2::AdaptiveScannerO2 (AdaptiveScannerO2 *ds)
+                 : DirectionalScanner (ds)
+{
+  templ_a = ds->templ_a;
+  templ_b = ds->templ_b;
+  templ_nu = ds->templ_nu;
+  dlc1 = ds->dlc1;
+}
+
+
+DirectionalScanner *AdaptiveScannerO2::getCopy ()
+{
+  return (new AdaptiveScannerO2 (this));
+}
+
+
+int AdaptiveScannerO2::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -138,7 +154,7 @@ int AdaptiveScannerO2::first (vector<Pt2i> &scan)
 int AdaptiveScannerO2::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcy --;
   // Whenever the control line changed
   while (lcx > xmin && lcy < ymax && dla * lcx + dlb * lcy > dlc1)
@@ -178,7 +194,7 @@ int AdaptiveScannerO2::nextOnLeft (vector<Pt2i> &scan)
 int AdaptiveScannerO2::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcy ++;
   while (rcx > xmin && rcy < ymax && dla * rcx + dlb * rcy > dlc1)
   {
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero2.h b/Code/FBSD/DirectionalScanner/adaptivescannero2.h
index d7588bbea9dd1e0e3ea4f0198efe3b146718852a..11ec2d2fbbc18c4c467cd115142f3acc25281e73 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero2.h
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero2.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class AdaptiveScannerO2 adaptivescannero2.h
  * \brief Adaptive directional scanner for the 2nd octant.
- * \author {P. Even}
  */
 class AdaptiveScannerO2 : public DirectionalScanner
 {
@@ -88,13 +87,19 @@ public:
                      int nbs, bool *steps,
                      int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -132,6 +137,15 @@ protected :
 
   /** Shift coefficient of the discrete lower support line */
   int dlc1;
+
+
+  /**
+   * @fn AdaptiveScannerO2(AdaptiveScannerO2 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  AdaptiveScannerO2 (AdaptiveScannerO2 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero7.cpp b/Code/FBSD/DirectionalScanner/adaptivescannero7.cpp
index 557a1e47add48d7c4620a5705739fe14f6f063e3..0192faa0617d5f0924396d354ec2e74739c7d7f3 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero7.cpp
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero7.cpp
@@ -113,8 +113,23 @@ AdaptiveScannerO7::AdaptiveScannerO7 (
 }
 
 
+AdaptiveScannerO7::AdaptiveScannerO7 (AdaptiveScannerO7 *ds)
+                 : DirectionalScanner (ds)
+{
+  templ_a = ds->templ_a;
+  templ_b = ds->templ_b;
+  templ_nu = ds->templ_nu;
+  dlc1 = ds->dlc1;
+}
+
+
+DirectionalScanner *AdaptiveScannerO7::getCopy ()
+{
+  return (new AdaptiveScannerO7 (this));
+}
+
 
-int AdaptiveScannerO7::first (vector<Pt2i> &scan)
+int AdaptiveScannerO7::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -139,7 +154,7 @@ int AdaptiveScannerO7::first (vector<Pt2i> &scan)
 int AdaptiveScannerO7::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcy ++;
   while (lcx < xmax - 1 && lcy < ymax && dla * lcx + dlb * lcy < dlc1)
   {
@@ -178,7 +193,7 @@ int AdaptiveScannerO7::nextOnLeft (vector<Pt2i> &scan)
 int AdaptiveScannerO7::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcy --;
   // Whenever the control corridor changed
   while (rcx < xmax - 1 && rcy < ymax && dla * rcx + dlb * rcy < dlc1)
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero7.h b/Code/FBSD/DirectionalScanner/adaptivescannero7.h
index ffb2bfa93b094ca65ad1bdabc7a3687eecb7f71d..6d0a5a591e45d2295a5dea3df56ac76639a782a1 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero7.h
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero7.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class AdaptiveScannerO7 adaptivescannero7.h
  * \brief Adaptive directional scanner for the 7th octant.
- * \author {P. Even}
  */
 class AdaptiveScannerO7 : public DirectionalScanner
 {
@@ -88,13 +87,19 @@ public:
                      int nbs, bool *steps,
                      int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -132,6 +137,15 @@ protected :
 
   /** Shift coefficient of the discrete lower support line */
   int dlc1;
+
+
+  /**
+   * @fn AdaptiveScannerO7(AdaptiveScannerO7 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  AdaptiveScannerO7 (AdaptiveScannerO7 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero8.cpp b/Code/FBSD/DirectionalScanner/adaptivescannero8.cpp
index 83ec6e6bde93bef4c92a3765ee532734d18b898f..7e0ccc2b05bb3110a4a28f5dfdaabc3bdbf7768d 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero8.cpp
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero8.cpp
@@ -113,8 +113,23 @@ AdaptiveScannerO8::AdaptiveScannerO8 (
 }
 
 
+AdaptiveScannerO8::AdaptiveScannerO8 (AdaptiveScannerO8 *ds)
+                 : DirectionalScanner (ds)
+{
+  templ_a = ds->templ_a;
+  templ_b = ds->templ_b;
+  templ_nu = ds->templ_nu;
+  dlc1 = ds->dlc1;
+}
+
+
+DirectionalScanner *AdaptiveScannerO8::getCopy ()
+{
+  return (new AdaptiveScannerO8 (this));
+}
+
 
-int AdaptiveScannerO8::first (vector<Pt2i> &scan)
+int AdaptiveScannerO8::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -139,7 +154,7 @@ int AdaptiveScannerO8::first (vector<Pt2i> &scan)
 int AdaptiveScannerO8::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcx --;
   while (lcy < ymax - 1 && lcx < xmax && dla * lcx + dlb * lcy < dlc1)
   {
@@ -178,7 +193,7 @@ int AdaptiveScannerO8::nextOnLeft (vector<Pt2i> &scan)
 int AdaptiveScannerO8::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcx ++;
   // Whenever the control corridor changed
   while (rcy < ymax - 1 && rcx < xmax && dla * rcx + dlb * rcy < dlc1)
diff --git a/Code/FBSD/DirectionalScanner/adaptivescannero8.h b/Code/FBSD/DirectionalScanner/adaptivescannero8.h
index fe057cc82a548474c6e8a4857b83f9b2b6dfbacb..d66d7bfeb15e6ebaaa85b374cc1437cb6181f3ee 100755
--- a/Code/FBSD/DirectionalScanner/adaptivescannero8.h
+++ b/Code/FBSD/DirectionalScanner/adaptivescannero8.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class AdaptiveScannerO8 adaptivescannero8.h
  * \brief Adaptive directional scanner for the 8th octant.
- * \author {P. Even}
  */
 class AdaptiveScannerO8 : public DirectionalScanner
 {
@@ -87,13 +86,19 @@ public:
                      int nbs, bool *steps,
                      int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -131,6 +136,15 @@ protected :
 
   /** Shift coefficient of the discrete lower support line */
   int dlc1;
+
+
+  /**
+   * @fn AdaptiveScannerO8(AdaptiveScannerO8 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  AdaptiveScannerO8 (AdaptiveScannerO8 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/directionalscanner.h b/Code/FBSD/DirectionalScanner/directionalscanner.h
index 38595af49bfd9b440e4ec9c09240eb4032b74c39..f212c0e4499ebdaea2ad7e812ea6da2336ccb794 100755
--- a/Code/FBSD/DirectionalScanner/directionalscanner.h
+++ b/Code/FBSD/DirectionalScanner/directionalscanner.h
@@ -13,7 +13,6 @@ using namespace std;
  * @class DirectionalScanner directionalscanner.h
  * \brief Incremental directional scanner.
  * This scanner iterately provides parallel scan lines.
- * \author {P. Even and B. Kerautret}
  */
 class DirectionalScanner
 {
@@ -26,12 +25,18 @@ public:
   virtual ~DirectionalScanner ();
 
   /**
-   * @fn first(vector<Pt2i> &scan)
+   * @fn virtual DirectionalScanner *getCopy () = 0;
+   * \brief Returns a copy of the directional scanner.
+   */
+  virtual DirectionalScanner *getCopy () = 0;
+
+  /**
+   * @fn first(vector<Pt2i> &scan) const
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  virtual int first (vector<Pt2i> &scan) = 0;
+  virtual int first (vector<Pt2i> &scan) const = 0;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -61,12 +66,18 @@ public:
 
   /**
    * @fn Pt2i locate (const Pt2i &pt) const
-   * \brief Returns the scanner coordinates of the givent point.
+   * \brief Returns the scanner coordinates of given point.
    * Scanner coordinates are the stripe index and the position in the stripe.
    * @param pt Image coordinates of the point.
    */
   virtual Pt2i locate (const Pt2i &pt) const;
 
+  /**
+   * @fn void releaseClearance ()
+   * \brief Releases output vector clearance before filling.
+   */
+  inline void releaseClearance () { clearance = false; }
+
 
 protected:
 
@@ -88,6 +99,11 @@ protected:
   /** Current step in scan direction. */
   bool *lst2, *rst2;
 
+  /** Flag indicating if the output vector should be cleared before filling.
+      Set to true by default. */
+  bool clearance;
+
+
   DirectionalScanner () { }
 
   /**
@@ -110,7 +126,21 @@ protected:
                       int nb, bool* st, int sx, int sy)
              : xmin (xmini), ymin (ymini), xmax (xmaxi), ymax (ymaxi),
                nbs (nb), steps (st),
-               ccx (sx), ccy (sy), lcx (sx), lcy (sy), rcx (sx), rcy (sy) { }
+               ccx (sx), ccy (sy), lcx (sx), lcy (sy), rcx (sx), rcy (sy),
+               clearance (true) { }
+
+  /**
+   * @fn DirectionalScanner(DirectionalScanner *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  DirectionalScanner (DirectionalScanner *ds)
+         : xmin (ds->xmin), ymin (ds->ymin), xmax (ds->xmax), ymax (ds->ymax),
+           dla (ds->dla), dlb (ds->dlb), dlc2 (ds->dlc2),
+           nbs (ds->nbs), steps (ds->steps), fs (ds->fs),
+           ccx (ds->ccx), ccy (ds->ccy),
+           lcx (ds->lcx), lcy (ds->lcy), rcx (ds->rcx), rcy (ds->rcy),
+           lst2 (ds->lst2), rst2 (ds->rst2), clearance (ds->clearance) { }
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero1.cpp b/Code/FBSD/DirectionalScanner/directionalscannero1.cpp
index 56cea0c0b5c97fd987e2a19d7f228840605811c0..1d3d3c32fe8966f9b5ae490dd35d437a844e9091 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero1.cpp
+++ b/Code/FBSD/DirectionalScanner/directionalscannero1.cpp
@@ -114,7 +114,23 @@ ccy = lcy;
 }
 
 
-int DirectionalScannerO1::first (vector<Pt2i> &scan)
+DirectionalScannerO1::DirectionalScannerO1 (DirectionalScannerO1 *ds)
+                    : DirectionalScanner (ds)
+{
+  lst1 = ds->lst1;
+  rst1 = ds->rst1;
+  ltransition = ds->ltransition;
+  rtransition = ds->rtransition;
+}
+
+
+DirectionalScanner *DirectionalScannerO1::getCopy ()
+{
+  return (new DirectionalScannerO1 (this));
+}
+
+
+int DirectionalScannerO1::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -139,7 +155,7 @@ int DirectionalScannerO1::first (vector<Pt2i> &scan)
 int DirectionalScannerO1::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (ltransition)
   {
     lcy --;
@@ -187,7 +203,7 @@ int DirectionalScannerO1::nextOnLeft (vector<Pt2i> &scan)
 int DirectionalScannerO1::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (rtransition)
   {
     rcx ++;
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero1.h b/Code/FBSD/DirectionalScanner/directionalscannero1.h
index e30df677167ed6d46760a1b9d853bd572fccbd6d..73129956c34fd34c1f7f3f0ff906069b8b47af5b 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero1.h
+++ b/Code/FBSD/DirectionalScanner/directionalscannero1.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class DirectionalScannerO1 directionalscannero1.h
  * \brief Incremental directional scanner for the 1st octant.
- * \author {P. Even and B. Kerautret}
  */
 class DirectionalScannerO1 : public DirectionalScanner
 {
@@ -87,13 +86,19 @@ public:
                         int nbs, bool *steps,
                         int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -127,6 +132,15 @@ private:
 
   /** State of the scan */
   bool ltransition, rtransition;
+
+
+  /**
+   * @fn DirectionalScannerO1(DirectionalScannerO1 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  DirectionalScannerO1 (DirectionalScannerO1 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero2.cpp b/Code/FBSD/DirectionalScanner/directionalscannero2.cpp
index 7fe62fcdd843a5a14dafd1fb5f175446d35b4679..d5277a755f77ef10ffd438eb7f8346d5dd9f73ed 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero2.cpp
+++ b/Code/FBSD/DirectionalScanner/directionalscannero2.cpp
@@ -114,7 +114,23 @@ ccy = lcy;
 }
 
 
-int DirectionalScannerO2::first (vector<Pt2i> &scan)
+DirectionalScannerO2::DirectionalScannerO2 (DirectionalScannerO2 *ds)
+                    : DirectionalScanner (ds)
+{
+  lst1 = ds->lst1;
+  rst1 = ds->rst1;
+  ltransition = ds->ltransition;
+  rtransition = ds->rtransition;
+}
+
+
+DirectionalScanner *DirectionalScannerO2::getCopy ()
+{
+  return (new DirectionalScannerO2 (this));
+}
+
+
+int DirectionalScannerO2::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -139,7 +155,7 @@ int DirectionalScannerO2::first (vector<Pt2i> &scan)
 int DirectionalScannerO2::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (ltransition)
   {
     lcy --;
@@ -185,7 +201,7 @@ int DirectionalScannerO2::nextOnLeft (vector<Pt2i> &scan)
 int DirectionalScannerO2::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (rtransition)
   {
     rcx ++;
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero2.h b/Code/FBSD/DirectionalScanner/directionalscannero2.h
index 0c7948af288830227a6ff672e8078d4065423257..67c8e4142ec61dab8cbff914398dc3a741340992 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero2.h
+++ b/Code/FBSD/DirectionalScanner/directionalscannero2.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class DirectionalScannerO2 directionalscannero2.h
  * \brief Incremental directional scanner for the 2nd octant.
- * \author {P. Even and B. Kerautret}
  */
 class DirectionalScannerO2 : public DirectionalScanner
 {
@@ -87,13 +86,19 @@ public:
                         int nbs, bool *steps,
                         int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -127,6 +132,15 @@ private:
 
   /** State of the scan */
   bool ltransition, rtransition;
+
+
+  /**
+   * @fn DirectionalScannerO2(DirectionalScannerO2 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  DirectionalScannerO2 (DirectionalScannerO2 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero7.cpp b/Code/FBSD/DirectionalScanner/directionalscannero7.cpp
index c2c40de3bd9438aa494cbda894bd3c5430a25d5c..d0e0ea529b2c6a64501cc794fef9997018f290ef 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero7.cpp
+++ b/Code/FBSD/DirectionalScanner/directionalscannero7.cpp
@@ -114,8 +114,23 @@ ccy = lcy;
 }
 
 
+DirectionalScannerO7::DirectionalScannerO7 (DirectionalScannerO7 *ds)
+                    : DirectionalScanner (ds)
+{
+  lst1 = ds->lst1;
+  rst1 = ds->rst1;
+  ltransition = ds->ltransition;
+  rtransition = ds->rtransition;
+}
+
+
+DirectionalScanner *DirectionalScannerO7::getCopy ()
+{
+  return (new DirectionalScannerO7 (this));
+}
+
 
-int DirectionalScannerO7::first (vector<Pt2i> &scan)
+int DirectionalScannerO7::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -140,7 +155,7 @@ int DirectionalScannerO7::first (vector<Pt2i> &scan)
 int DirectionalScannerO7::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (ltransition)
   {
     lcx --;
@@ -187,7 +202,7 @@ int DirectionalScannerO7::nextOnLeft (vector<Pt2i> &scan)
 int DirectionalScannerO7::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (rtransition)
   {
     rcy --;
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero7.h b/Code/FBSD/DirectionalScanner/directionalscannero7.h
index fb3f7ae9ba9af46fec1a0aa56f0bd5c4df6232d4..fa7dcc8465a463b11d5bf51bce72383272259374 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero7.h
+++ b/Code/FBSD/DirectionalScanner/directionalscannero7.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class DirectionalScannerO7 directionalscannero7.h
  * \brief Incremental directional scanner for the 7th octant.
- * \author {P. Even and B. Kerautret}
  */
 class DirectionalScannerO7 : public DirectionalScanner
 {
@@ -87,13 +86,19 @@ public:
                         int nbs, bool *steps,
                         int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -127,6 +132,15 @@ private:
 
   /** State of the scan */
   bool ltransition, rtransition;
+
+
+  /**
+   * @fn DirectionalScannerO7(DirectionalScannerO7 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  DirectionalScannerO7 (DirectionalScannerO7 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero8.cpp b/Code/FBSD/DirectionalScanner/directionalscannero8.cpp
index 123b84adcf912d9ede1cd1f56c0eb3453d1cdb14..311456f43fdbedacc3bd3813ac648d23dcb5a802 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero8.cpp
+++ b/Code/FBSD/DirectionalScanner/directionalscannero8.cpp
@@ -114,8 +114,23 @@ ccy = lcy;
 }
 
 
+DirectionalScannerO8::DirectionalScannerO8 (DirectionalScannerO8 *ds)
+                    : DirectionalScanner (ds)
+{
+  lst1 = ds->lst1;
+  rst1 = ds->rst1;
+  ltransition = ds->ltransition;
+  rtransition = ds->rtransition;
+}
+
+
+DirectionalScanner *DirectionalScannerO8::getCopy ()
+{
+  return (new DirectionalScannerO8 (this));
+}
+
 
-int DirectionalScannerO8::first (vector<Pt2i> &scan)
+int DirectionalScannerO8::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
   bool *nst = steps;         // Current step in scan direction (jpts)
@@ -140,7 +155,7 @@ int DirectionalScannerO8::first (vector<Pt2i> &scan)
 int DirectionalScannerO8::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (ltransition)
   {
     lcx --;
@@ -186,7 +201,7 @@ int DirectionalScannerO8::nextOnLeft (vector<Pt2i> &scan)
 int DirectionalScannerO8::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   if (rtransition)
   {
     rcy --;
diff --git a/Code/FBSD/DirectionalScanner/directionalscannero8.h b/Code/FBSD/DirectionalScanner/directionalscannero8.h
index b5e5486514dce19e9ef10ffc3e69537cd19b675a..b753ce7a1fb7a0834882f0336e96a5f5431e049b 100755
--- a/Code/FBSD/DirectionalScanner/directionalscannero8.h
+++ b/Code/FBSD/DirectionalScanner/directionalscannero8.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class DirectionalScannerO8 directionalscannero8.h
  * \brief Incremental directional scanner for the 8th octant.
- * \author {P. Even and B. Kerautret}
  */
 class DirectionalScannerO8 : public DirectionalScanner
 {
@@ -87,13 +86,19 @@ public:
                         int nbs, bool *steps,
                         int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -127,6 +132,15 @@ private:
 
   /** State of the scan */
   bool ltransition, rtransition;
+
+
+  /**
+   * @fn DirectionalScannerO8(DirectionalScannerO8 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  DirectionalScannerO8 (DirectionalScannerO8 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/scannerprovider.cpp b/Code/FBSD/DirectionalScanner/scannerprovider.cpp
index dd6f8a057eeb0e769a1c3832235dd91b0984712a..b33cdf02b1b68def614206b303884b15febb2c62 100755
--- a/Code/FBSD/DirectionalScanner/scannerprovider.cpp
+++ b/Code/FBSD/DirectionalScanner/scannerprovider.cpp
@@ -12,16 +12,18 @@
 
 
 
-DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2)
+DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2,
+                                                 bool controlable)
 {
   // Enforces P1 to be lower than P2
-  // or to left of P2 in cas of equality
-  if ((p1.y () > p2.y ())
-      || ((p1.y () == p2.y ()) && (p1.x () > p2.x ())))
+  // or to left of P2 in case of equality
+  lastScanReversed = (p1.y () > p2.y ())
+                     || ((p1.y () == p2.y ()) && (p1.x () > p2.x ()));
+  if (lastScanReversed)
   {
-    Pt2i tmp = p1;
-    p1 = p2;
-    p2 = tmp;
+    Pt2i tmp (p1);
+    p1.set (p2);
+    p2.set (tmp);
   }
 
   // Computes the steps position array
@@ -50,7 +52,11 @@ DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2)
         return (new VHScannerO1 (xmin, ymin, xmax, ymax,
                                  a, b, c2, nbs, steps, repx, repy));
       }
-      else return (new DirectionalScannerO1 (xmin, ymin, xmax, ymax,
+      else return (controlable ?
+                   (DirectionalScanner *)
+                   new AdaptiveScannerO1 (xmin, ymin, xmax, ymax,
+                          a, b, c2, nbs, steps, p1.x (), p1.y ()) :
+                   new DirectionalScannerO1 (xmin, ymin, xmax, ymax,
                           a, b, c2, nbs, steps, p1.x (), p1.y ()));
     }
     else
@@ -63,7 +69,11 @@ DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2)
         return (new VHScannerO2 (xmin, ymin, xmax, ymax,
                                  a, b, c2, nbs, steps, repx, repy));
       }
-      else return (new DirectionalScannerO2 (xmin, ymin, xmax, ymax,
+      else return (controlable ?
+                   (DirectionalScanner *)
+                   new AdaptiveScannerO2 (xmin, ymin, xmax, ymax,
+                          a, b, c2, nbs, steps, p1.x (), p1.y ()) :
+                   new DirectionalScannerO2 (xmin, ymin, xmax, ymax,
                           a, b, c2, nbs, steps, p1.x (), p1.y ()));
     }
   else
@@ -77,7 +87,11 @@ DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2)
         return (new VHScannerO8 (xmin, ymin, xmax, ymax,
                                  a, b, c2, nbs, steps, repx, repy));
       }
-      else return (new DirectionalScannerO8 (xmin, ymin, xmax, ymax,
+      else return (controlable ?
+                   (DirectionalScanner *)
+                   new AdaptiveScannerO8 (xmin, ymin, xmax, ymax,
+                          a, b, c2, nbs, steps, p1.x (), p1.y ()) :
+                   new DirectionalScannerO8 (xmin, ymin, xmax, ymax,
                           a, b, c2, nbs, steps, p1.x (), p1.y ()));
     }
     else
@@ -90,7 +104,11 @@ DirectionalScanner *ScannerProvider::getScanner (Pt2i p1, Pt2i p2)
         return (new VHScannerO7 (xmin, ymin, xmax, ymax,
                                  a, b, c2, nbs, steps, repx, repy));
       }
-      else return (new DirectionalScannerO7 (xmin, ymin, xmax, ymax,
+      else return (controlable ?
+                   (DirectionalScanner *)
+                   new AdaptiveScannerO7 (xmin, ymin, xmax, ymax,
+                          a, b, c2, nbs, steps, p1.x (), p1.y ()) :
+                   new DirectionalScannerO7 (xmin, ymin, xmax, ymax,
                           a, b, c2, nbs, steps, p1.x (), p1.y ()));
     }
 }
diff --git a/Code/FBSD/DirectionalScanner/scannerprovider.h b/Code/FBSD/DirectionalScanner/scannerprovider.h
index ebafab9529b2f04514a0f23d91b238789874adf3..47c227b2dd37943ed50b8f231a9e67f1aa1d0283 100755
--- a/Code/FBSD/DirectionalScanner/scannerprovider.h
+++ b/Code/FBSD/DirectionalScanner/scannerprovider.h
@@ -13,7 +13,6 @@ using namespace std;
  * \brief Directional scanner provider.
  * Provides ad-hoc directional scanners in the relevant octant
  *   and according to static or dynamical needs.
- * \author {P. Even}
  */
 class ScannerProvider
 {
@@ -23,7 +22,7 @@ public:
    * @fn ScannerProvider()
    * \brief Builds a directional scanner provider.
    */
-  ScannerProvider () : isOrtho (false),
+  ScannerProvider () : isOrtho (false), lastScanReversed (false),
                        xmin (0), ymin (0), xmax (100), ymax (100) { }
   
   /**
@@ -54,8 +53,9 @@ public:
    *   defined by control points p1 and p2.
    * @param p1 Start control point.
    * @param p2 End control point.
+   * @param controlable controlability request (true for a dynamical scanner).
    */
-  DirectionalScanner *getScanner (Pt2i p1, Pt2i p2);
+  DirectionalScanner *getScanner (Pt2i p1, Pt2i p2, bool controlable = false);
   
   /**
    * @fn getScanner(Pt2i p1, Pt2i p2, Pt2i v1, Pt2i v2)
@@ -97,6 +97,12 @@ public:
   DirectionalScanner *getScanner (Pt2i centre, Vr2i normal,
                                   int length, bool controlable);
 
+  /**
+   * @fn bool isLastScanReversed () const
+   * \brief Returns whether the currently used scan end points were permutated.
+   */
+  inline bool isLastScanReversed () const { return lastScanReversed; }
+
   /**
    * @fn setOrtho(bool status)
    * \brief Sets the orthogonal scanner modality.
@@ -109,6 +115,8 @@ private:
 
   /** Orthogonal scanner modality. */
   bool isOrtho;
+  /** Last scan end points permutation modality. */
+  bool lastScanReversed;
 
   /** Scan area lowest x coordinate. */
   int xmin;
diff --git a/Code/FBSD/DirectionalScanner/vhscannero1.cpp b/Code/FBSD/DirectionalScanner/vhscannero1.cpp
index c03abfd437f602fd7b058d7bc42f644fd5cb17ea..bf29f152904bcba3ced8644c87a8c29e5148e665 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero1.cpp
+++ b/Code/FBSD/DirectionalScanner/vhscannero1.cpp
@@ -105,7 +105,18 @@ VHScannerO1::VHScannerO1 (int xmin, int ymin, int xmax, int ymax,
 }
 
 
-int VHScannerO1::first (vector<Pt2i> &scan)
+VHScannerO1::VHScannerO1 (VHScannerO1 *ds) : AdaptiveScannerO1 (ds)
+{
+}
+
+
+DirectionalScanner *VHScannerO1::getCopy ()
+{
+  return (new VHScannerO1 (this));
+}
+
+
+int VHScannerO1::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
 
@@ -125,7 +136,7 @@ int VHScannerO1::first (vector<Pt2i> &scan)
 int VHScannerO1::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcx --;
   if (lcx < xmin) return 0;
 
@@ -158,7 +169,7 @@ int VHScannerO1::nextOnLeft (vector<Pt2i> &scan)
 int VHScannerO1::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcx ++;
   if (rcx >= xmax) return 0;
 
diff --git a/Code/FBSD/DirectionalScanner/vhscannero1.h b/Code/FBSD/DirectionalScanner/vhscannero1.h
index e2aac2440ec20e7c710c882e3961f0831b2086c6..6c2b4fbf0302a4b23b120f5b1b98b5937aa52359 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero1.h
+++ b/Code/FBSD/DirectionalScanner/vhscannero1.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class VHScannerO1 vhscannero1.h
  * \brief Vertical/horizontal adaptive DS for the 1st octant.
- * \author {P. Even}
  */
 class VHScannerO1 : public AdaptiveScannerO1
 {
@@ -86,13 +85,19 @@ public:
                int a, int b, int nbs, bool *steps,
                int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -109,6 +114,17 @@ public:
    * @param scan vectors of points to be filled in.
    */
   int nextOnRight (vector<Pt2i> &scan);
+
+
+private:
+
+  /**
+   * @fn VHScannerO1(VHScannerO1 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  VHScannerO1 (VHScannerO1 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/vhscannero2.cpp b/Code/FBSD/DirectionalScanner/vhscannero2.cpp
index 418849733e7ca8eb87b7ae5a78069965cc7712c6..fa518b2d6b4baa59f6486d0af29dcb9fabe29a78 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero2.cpp
+++ b/Code/FBSD/DirectionalScanner/vhscannero2.cpp
@@ -105,7 +105,18 @@ VHScannerO2::VHScannerO2 (int xmin, int ymin, int xmax, int ymax,
 }
 
 
-int VHScannerO2::first (vector<Pt2i> &scan)
+VHScannerO2::VHScannerO2 (VHScannerO2 *ds) : AdaptiveScannerO2 (ds)
+{
+}
+
+
+DirectionalScanner *VHScannerO2::getCopy ()
+{
+  return (new VHScannerO2 (this));
+}
+
+
+int VHScannerO2::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
 
@@ -125,7 +136,7 @@ int VHScannerO2::first (vector<Pt2i> &scan)
 int VHScannerO2::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcy --;
   if (lcy < ymin) return 0;
 
@@ -158,7 +169,7 @@ int VHScannerO2::nextOnLeft (vector<Pt2i> &scan)
 int VHScannerO2::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcy ++;
   if (rcy >= ymax) return 0;
 
diff --git a/Code/FBSD/DirectionalScanner/vhscannero2.h b/Code/FBSD/DirectionalScanner/vhscannero2.h
index 2457c39f77d2db7a68366f679aa9e4433d20a9e9..ef21ca457cebbac162352e3b7c654b8dccbd27ef 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero2.h
+++ b/Code/FBSD/DirectionalScanner/vhscannero2.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class VHScannerO2 vhscannero2.h
  * \brief Vertical / horizontal adpative DS for the 2nd octant.
- * \author {P. Even}
  */
 class VHScannerO2 : public AdaptiveScannerO2
 {
@@ -86,13 +85,19 @@ public:
                int a, int b, int nbs, bool *steps,
                int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -109,6 +114,17 @@ public:
    * @param scan vectors of points to be filled in.
    */
   int nextOnRight (vector<Pt2i> &scan);
+
+
+private :
+
+  /**
+   * @fn VHScannerO2(VHScannerO2 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  VHScannerO2 (VHScannerO2 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/vhscannero7.cpp b/Code/FBSD/DirectionalScanner/vhscannero7.cpp
index 119b7145ed6fcf6f0f156c4172ac0208893ed51d..e24279adbb5e6e36ab1fc0b6acc5014fa95449ab 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero7.cpp
+++ b/Code/FBSD/DirectionalScanner/vhscannero7.cpp
@@ -105,8 +105,18 @@ VHScannerO7::VHScannerO7 (int xmin, int ymin, int xmax, int ymax,
 }
 
 
+VHScannerO7::VHScannerO7 (VHScannerO7 *ds) : AdaptiveScannerO7 (ds)
+{
+}
+
+
+DirectionalScanner *VHScannerO7::getCopy ()
+{
+  return (new VHScannerO7 (this));
+}
+
 
-int VHScannerO7::first (vector<Pt2i> &scan)
+int VHScannerO7::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
 
@@ -126,7 +136,7 @@ int VHScannerO7::first (vector<Pt2i> &scan)
 int VHScannerO7::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcy ++;
   if (lcy >= ymax) return 0;
 
@@ -158,7 +168,7 @@ int VHScannerO7::nextOnLeft (vector<Pt2i> &scan)
 int VHScannerO7::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcy --;
   if (rcy < ymin) return 0;
 
diff --git a/Code/FBSD/DirectionalScanner/vhscannero7.h b/Code/FBSD/DirectionalScanner/vhscannero7.h
index 62c4445b247ce17f9ecaaf9b689e447f5700222d..9d79bc600655884ca43b71715a037ab47eebd6f6 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero7.h
+++ b/Code/FBSD/DirectionalScanner/vhscannero7.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class VHScannerO7 vhscannero7.h
  * \brief Vertical / horizontal adaptive DS for the 7th octant.
- * \author {P. Even}
  */
 class VHScannerO7 : public AdaptiveScannerO7
 {
@@ -86,13 +85,19 @@ public:
                int a, int b, int nbs, bool *steps,
                int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -109,6 +114,17 @@ public:
    * @param scan vectors of points to be filled in.
    */
   int nextOnRight (vector<Pt2i> &scan);
+
+
+private :
+
+  /**
+   * @fn VHScannerO7(VHScannerO7 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  VHScannerO7 (VHScannerO7 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/DirectionalScanner/vhscannero8.cpp b/Code/FBSD/DirectionalScanner/vhscannero8.cpp
index 7efb53d6567105e430b933fe28c0b865c9c6f415..b96c3f2ff42c952c26ac4f77a0b6d0706dd75186 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero8.cpp
+++ b/Code/FBSD/DirectionalScanner/vhscannero8.cpp
@@ -105,8 +105,18 @@ VHScannerO8::VHScannerO8 (int xmin, int ymin, int xmax, int ymax,
 }
 
 
+VHScannerO8::VHScannerO8 (VHScannerO8 *ds) : AdaptiveScannerO8 (ds)
+{
+}
+
+
+DirectionalScanner *VHScannerO8::getCopy ()
+{
+  return (new VHScannerO8 (this));
+}
+
 
-int VHScannerO8::first (vector<Pt2i> &scan)
+int VHScannerO8::first (vector<Pt2i> &scan) const
 {
   int x = lcx, y = lcy;      // Current position coordinates
 
@@ -126,7 +136,7 @@ int VHScannerO8::first (vector<Pt2i> &scan)
 int VHScannerO8::nextOnLeft (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   lcx --;
   if (lcx < xmin) return 0;
 
@@ -158,7 +168,7 @@ int VHScannerO8::nextOnLeft (vector<Pt2i> &scan)
 int VHScannerO8::nextOnRight (vector<Pt2i> &scan)
 {
   // Prepares the next scan
-  scan.clear ();
+  if (clearance) scan.clear ();
   rcx ++;
   if (rcx >= xmax) return 0;
 
diff --git a/Code/FBSD/DirectionalScanner/vhscannero8.h b/Code/FBSD/DirectionalScanner/vhscannero8.h
index 72fcc248938b9502d4e166bb96406795070ff701..f989daf397e4f393ff0639343f9b28792c9d05e3 100755
--- a/Code/FBSD/DirectionalScanner/vhscannero8.h
+++ b/Code/FBSD/DirectionalScanner/vhscannero8.h
@@ -10,7 +10,6 @@ using namespace std;
 /** 
  * @class VHScannerO8 vhscannero8.h
  * \brief Vertical / horizontal adaptive DS for the 8th octant.
- * \author {P. Even}
  */
 class VHScannerO8 : public AdaptiveScannerO8
 {
@@ -86,13 +85,19 @@ public:
                int a, int b, int nbs, bool *steps,
                int cx, int cy, int length);
 
+  /**
+   * @fn DirectionalScanner *getCopy ();
+   * \brief Returns a copy of the directional scanner.
+   */
+  DirectionalScanner *getCopy ();
+
   /**
    * @fn first(vector<Pt2i> &scan)
    * \brief Returns the central scan.
    * Fills in the vector with the central scan and returns its size.
    * @param scan vectors of points to be filled in.
    */
-  int first (vector<Pt2i> &scan);
+  int first (vector<Pt2i> &scan) const;
 
   /**
    * @fn nextOnLeft(vector<Pt2i> &scan)
@@ -109,6 +114,17 @@ public:
    * @param scan vectors of points to be filled in.
    */
   int nextOnRight (vector<Pt2i> &scan);
+
+
+private :
+
+  /**
+   * @fn VHScannerO8(VHScannerO8 *ds)
+   * \brief Creates a copy of given directional scanner.
+   * @param ds Source directional scanner.
+   */
+  VHScannerO8 (VHScannerO8 *ds);
+
 };
 
 #endif
diff --git a/Code/FBSD/FBSD.pro b/Code/FBSD/FBSD.pro
index 38325a9670507fda824e5283e0b9631b148aa1e3..2a03ddf33c75c4bf9c2048c733aa164c52281aca 100644
--- a/Code/FBSD/FBSD.pro
+++ b/Code/FBSD/FBSD.pro
@@ -15,10 +15,10 @@ OBJECTS_DIR = obj
 
 # Input
 HEADERS += BlurredSegment/biptlist.h \
-           BlurredSegment/blurredsegment.h \
            BlurredSegment/blurredsegmentproto.h \
            BlurredSegment/bsdetector.h \
            BlurredSegment/bsfilter.h \
+           BlurredSegment/bsproto.h \
            BlurredSegment/bstracker.h \
            BSTools/bsdetectionwidget.h \
            BSTools/bsgroundtruthitem.h \
@@ -59,9 +59,9 @@ HEADERS += BlurredSegment/biptlist.h \
 SOURCES += main.cpp \
            BlurredSegment/biptlist.cpp \
            BlurredSegment/blurredsegment.cpp \
-           BlurredSegment/blurredsegmentproto.cpp \
            BlurredSegment/bsdetector.cpp \
            BlurredSegment/bsfilter.cpp \
+           BlurredSegment/bsproto.cpp \
            BlurredSegment/bstracker.cpp \
            BSTools/bsdetectionwidget.cpp \
            BSTools/bsgroundtruthitem.cpp \
diff --git a/Code/FBSD/ImageTools/digitalstraightsegment.cpp b/Code/FBSD/ImageTools/digitalstraightsegment.cpp
index 10c9bc885027681119c2d135a512acef85b1c875..22c5dbadb9f2818f02292eb6846673e8235986f7 100755
--- a/Code/FBSD/ImageTools/digitalstraightsegment.cpp
+++ b/Code/FBSD/ImageTools/digitalstraightsegment.cpp
@@ -3,6 +3,14 @@
 
 
 
+DigitalStraightSegment::DigitalStraightSegment ()
+                      : DigitalStraightLine (1, 1, 0, 1)
+{
+  min = 0;
+  max = 1;
+}
+
+
 DigitalStraightSegment::DigitalStraightSegment (Pt2i p1, Pt2i p2, int type,
                                  int xmin, int ymin, int xmax, int ymax)
                       : DigitalStraightLine (p1, p2, type)
@@ -64,6 +72,14 @@ DigitalStraightSegment::DigitalStraightSegment (int va, int vb, int vc,
 }
 
 
+DigitalStraightSegment::DigitalStraightSegment (DigitalStraightSegment *dss)
+                      : DigitalStraightLine (dss->a, dss->b, dss->c, dss->nu)
+{
+  min = dss->min;
+  max = dss->max;
+}
+
+
 Pt2i DigitalStraightSegment::getABoundingPoint (bool upper) const
 {
   int sa = a, sb = b, u1 = 1, v1 = 0, u2 = 0, v2 = 1;
@@ -232,3 +248,24 @@ bool DigitalStraightSegment::contains (Pt2i p, int tol) const
   else
     return (p.y () >= min && p.y () <= max);
 }
+
+
+int DigitalStraightSegment::length2 () const
+{
+  int numin, numax, den;
+  if (a < (b < 0 ? -b : b))
+  {
+    numin = c - a * min;
+    numax = c - a * max;
+    den = b;
+  }
+  else
+  {
+    numin = c - b * min;
+    numax = c - b * max;
+    den = a;
+  }
+  return ((int) (((max - min) * (max - min) * den * den
+                  + (numax - numin) * (numax - numin)
+                  + (den * den) / 2) / (den * den)));
+}
diff --git a/Code/FBSD/ImageTools/digitalstraightsegment.h b/Code/FBSD/ImageTools/digitalstraightsegment.h
index ca8c162518a1eaafa42a3b00b51e755a9edede88..340a30988899b454457230504ed9197b35781fce 100755
--- a/Code/FBSD/ImageTools/digitalstraightsegment.h
+++ b/Code/FBSD/ImageTools/digitalstraightsegment.h
@@ -21,6 +21,11 @@ class DigitalStraightSegment : public DigitalStraightLine
 {
 public:
 
+  /**
+   * Creates a default digital straight segment.
+   */
+  DigitalStraightSegment ();
+
   /**
    * Creates a null-thick segment joining two points.
    * @param p1 First point on the line.
@@ -55,6 +60,20 @@ public:
    */
   DigitalStraightSegment (Pt2i p1, Pt2i p2, int width);
 
+  /**
+   * Sets identical to the gven digital straight segment.
+   * @param dss Original straight segment.
+   */
+  inline void set (const DigitalStraightSegment &dss) {
+    a = dss.a; b = dss.b; c = dss.c; nu = dss.nu;
+    min = dss.min; max = dss.max; }
+
+  /**
+   * Creates a digital straight segment from another one.
+   * @param dss Pointer to the digital straight segment to copy.
+   */
+  DigitalStraightSegment (DigitalStraightSegment *dss);
+
   /**
    * \brief Returns a bounding point of the digital line
    * @param upper true for a upper bounding point, false for a lower one.
@@ -105,6 +124,11 @@ public:
    */
   bool contains (Pt2i p, int tol) const;
 
+  /**
+   * \brief Returns the square length of the digital straight segment.
+   */
+  int length2 () const;
+
 
 protected:
 
diff --git a/Code/FBSD/ImageTools/pt2i.cpp b/Code/FBSD/ImageTools/pt2i.cpp
index 0d60fe436474cd98ad17e9b22520d5abffdb1073..8b2d16133c470fb5dcaf5168f028895512d1bef4 100755
--- a/Code/FBSD/ImageTools/pt2i.cpp
+++ b/Code/FBSD/ImageTools/pt2i.cpp
@@ -1,4 +1,5 @@
 #include "pt2i.h"
+#include <inttypes.h>
 
 
 Pt2i::Pt2i ()
@@ -22,6 +23,91 @@ Pt2i::Pt2i (const Pt2i &p)
 }
 
 
+bool Pt2i::quadCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+{
+  int64_t x1 = (int64_t) xp, y1 = (int64_t) yp;
+  int64_t x2 = (int64_t) (p1.xp), y2 = (int64_t) (p1.yp);
+  int64_t x3 = (int64_t) (p2.xp), y3 = (int64_t) (p2.yp);
+  int64_t x4 = (int64_t) (p3.xp), y4 = (int64_t) (p3.yp);
+
+  return (((x3 - x1) * (y2 - y1) - (y3 - y1) * (x2 - x1))
+           * ((x3 - x1) * (y4 - y1) - (y3 - y1) * (x4 - x1)) <= 0);
+}
+
+
+bool Pt2i::orientCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+{
+  return ((p1.xp - xp) * (p2.xp - p3.xp) + (p1.yp - yp) * (p2.yp - p3.yp) > 0);
+}
+
+
+bool Pt2i::convexCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+{
+  int64_t x1 = (int64_t) (p1.xp - xp), y1 = (int64_t) (p1.yp - yp);
+  int64_t x2 = (int64_t) (p2.xp - p1.xp), y2 = (int64_t) (p2.yp - p1.yp);
+  int64_t x3 = (int64_t) (p3.xp - p2.xp), y3 = (int64_t) (p3.yp - p2.yp);
+  int64_t x4 = (int64_t) (xp - p3.xp), y4 = (int64_t) (yp - p3.yp);
+
+  int64_t o1 = x1 * y2 - x2 * y1;
+  int64_t o2 = x2 * y3 - x3 * y2;
+  int64_t o3 = x3 * y4 - x4 * y3;
+  int64_t o4 = x4 * y1 - x1 * y4;
+
+  int signe = (o1 == 0 ? (o3 > 0 ? 1 : -1) : (o1 > 0 ? 1 : -1));
+  return (signe * o1 >= 0 && signe * o2 >= 0
+          && signe * o3 >= 0 && signe * o4 >= 0);
+}
+
+
+bool Pt2i::inTriangle (Pt2i p1, Pt2i p2, Pt2i p3) const
+{
+  int64_t x0 = (int64_t) xp, y0 = (int64_t) yp;
+  int64_t x1 = (int64_t) (p1.xp), y1 = (int64_t) (p1.yp);
+  int64_t x2 = (int64_t) (p2.xp), y2 = (int64_t) (p2.yp);
+  int64_t x3 = (int64_t) (p3.xp), y3 = (int64_t) (p3.yp);
+
+  // Case of aligned triangle vertices
+  if ((x2 - x1) * (y3 - y1) == (x3 - x1) * (y2 - y1))
+  {
+    int64_t xmin = x1, ymin = y1, xmax = x1, ymax = y1;
+    if (x2 < xmin) xmin = x2;
+    else xmax = x2;
+    if (x3 < xmin) xmin = x3;
+    else if (x3 > xmax) xmax = x3;
+    if (y2 < ymin) ymin = y2;
+    else ymax = y2;
+    if (y3 < ymin) ymin = y3;
+    else if (y3 > ymax) ymax = y3;
+    return ((x2 - x1) * (y0 - y1) == (x0 - x1) * (y2 - y1)
+            && x0 >= xmin && x0 <= xmax && y0 >= ymin && y0 <= ymax);
+  }
+
+  int64_t pv1 = (x0 - x1) * (y2 - y1) - (y0 - y1) * (x2 - x1);
+  int64_t pv2 = (x0 - x2) * (y3 - y2) - (y0 - y2) * (x3 - x2);
+  int64_t pv3 = (x0 - x3) * (y1 - y3) - (y0 - y3) * (x1 - x3);
+  return ((pv1 >= 0 && pv2 >= 0 && pv3 >= 0)
+          || (pv1 <= 0 && pv2 <= 0 && pv3 <= 0));
+}
+
+
+bool Pt2i::inQuad (Pt2i p1, Pt2i p2, Pt2i p3, Pt2i p4) const
+{
+  int64_t x0 = (int64_t) xp, y0 = (int64_t) yp;
+  int64_t x1 = (int64_t) (p1.xp), y1 = (int64_t) (p1.yp);
+  int64_t x2 = (int64_t) (p2.xp), y2 = (int64_t) (p2.yp);
+  int64_t x3 = (int64_t) (p3.xp), y3 = (int64_t) (p3.yp);
+  int64_t x4 = (int64_t) (p4.xp), y4 = (int64_t) (p4.yp);
+
+  int64_t pv1 = (x0 - x1) * (y2 - y1) - (y0 - y1) * (x2 - x1);
+  int64_t pv2 = (x0 - x2) * (y3 - y2) - (y0 - y2) * (x3 - x2);
+  int64_t pv3 = (x0 - x3) * (y4 - y3) - (y0 - y3) * (x4 - x3);
+  int64_t pv4 = (x0 - x4) * (y1 - y4) - (y0 - y4) * (x1 - x4);
+
+  return ((pv1 >= 0 && pv2 >= 0 && pv3 >= 0 && pv4 >= 0)
+          || (pv1 <= 0 && pv2 <= 0 && pv3 <= 0 && pv4 <= 0));
+}
+
+
 Vr2i Pt2i::vectorTo (Pt2i p) const
 {
   return (Vr2i (p.xp - xp, p.yp - yp));
@@ -198,7 +284,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
       dy *= 2;
       while (x1 < x2)
       {
-        pts[i++] = Pt2i (x1, y1);
+        pts[i++].set (x1, y1);
         x1 ++;
         e -= dy;
         if (e < 0)
@@ -207,7 +293,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
           e += dx;
         }
       }
-      pts[i] = Pt2i (x1, y1);
+      pts[i].set (x1, y1);
     }
 
     // Octant 2
@@ -220,7 +306,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
       dy *= 2;
       while (y1 < y2)
       {
-        pts[i++] = Pt2i (x1, y1);
+        pts[i++].set (x1, y1);
         y1 ++;
         e -= dx;
         if (e < 0)
@@ -229,7 +315,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
           e += dy;
         }
       }
-      pts[i] = Pt2i (x1, y1);
+      pts[i].set (x1, y1);
     }
   }
 
@@ -245,7 +331,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
       dy *= 2;
       while (x1 < x2)
       {
-        pts[i++] = Pt2i (x1, y1);
+        pts[i++].set (x1, y1);
         x1 ++;
         e += dy;
         if (e < 0)
@@ -254,7 +340,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
           e += dx;
         }
       }
-      pts[i] = Pt2i (x1, y1);
+      pts[i].set (x1, y1);
     }
 
     // Octant 7
@@ -267,7 +353,7 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
       dy *= 2;
       while (y1 > y2)
       {
-        pts[i++] = Pt2i (x1, y1);
+        pts[i++].set (x1, y1);
         y1 --;
         e -= dx;
         if (e < 0)
@@ -276,9 +362,144 @@ Pt2i *Pt2i::drawing (const Pt2i p, int *n) const
           e -= dy;
         }
       }
-      pts[i] = Pt2i (x1, y1);
+      pts[i].set (x1, y1);
+    }
+  }
+  return (pts);
+}
+
+
+Pt2i *Pt2i::clipLine (const Pt2i p, int left, int low, int right, int up,
+                      int *n) const
+{
+  if (right < left) { int tmp = left; left = right; right = tmp; }
+  if (up < low) { int tmp = low; low = up; up = tmp; }
+  int x1, y1, x2, y2;
+  if (xp > p.xp)
+  {
+    x1 = p.xp;
+    x2 = xp;
+    y1 = p.yp;
+    y2 = yp;
+  }
+  else
+  {
+    x1 = xp;
+    x2 = p.xp;
+    y1 = yp;
+    y2 = p.yp;
+  }
+  int dx = x2 - x1;
+  int dy = y2 - y1;
+  int e, i = 0;
+  Pt2i *pts;
+
+  if (dy > 0)
+  {
+    // Octant 1
+    if (dx >= dy)
+    {
+      *n = dx + 1;
+      pts = new Pt2i[dx + 1];
+      if (x2 >= left && y2 >= low)
+      {
+        e = dx - 1; // car le point limite est trace en dessous de la droite
+        dx *= 2;
+        dy *= 2;
+        while (x1 < x2 && x1 <= right && y1 <= up)
+        {
+          if (x1 >= left && y1 >= low) pts[i++].set (x1, y1);
+          x1 ++;
+          e -= dy;
+          if (e < 0)
+          {
+            y1 ++;
+            e += dx;
+          }
+        }
+        if (x2 <= right && y2 <= up) pts[i++].set (x2, y2);
+      }
+    }
+
+    // Octant 2
+    else
+    {
+      *n = dy + 1;
+      pts = new Pt2i[dy + 1];
+      if (x2 >= left && y2 >= low)
+      {
+        e = dy; // car le point limite est trace a droite de la droite
+        dx *= 2;
+        dy *= 2;
+        while (y1 < y2 && x1 <= right && y1 <= up)
+        {
+          if (x1 >= left && y1 >= low) pts[i++].set (x1, y1);
+          y1 ++;
+          e -= dx;
+          if (e < 0)
+          {
+            x1 ++;
+            e += dy;
+          }
+        }
+        if (x2 <= right && y2 <= up) pts[i++].set (x2, y2);
+      }
+    }
+  }
+
+  else
+  {
+    // Octant 8
+    if (dx >= -dy)
+    {
+      *n = 1 + dx;
+      pts = new Pt2i[dx + 1];
+      if (x2 >= left && y2 <= up)
+      {
+        e = dx - 1; // car le point limite est trace en dessous de la droite
+        dx *= 2;
+        dy *= 2;
+        while (x1 < x2 && x1 <= right && y1 >= low)
+        {
+          if (x1 >= left && y1 <= up) pts[i++].set (x1, y1);
+          x1 ++;
+          e += dy;
+          if (e < 0)
+          {
+            y1 --;
+            e += dx;
+          }
+        }
+        if (x2 <= right && y2 >= low) pts[i++].set (x2, y2);
+      }
+    }
+
+    // Octant 7
+    else
+    {
+      *n = 1 - dy;
+      pts = new Pt2i[1 - dy];
+      if (x2 >= left && y2 <= up)
+      {
+        e = - dy; // car le point limite est trace a gauche de la droite
+        dx *= 2;
+        dy *= 2;
+        while (y1 > y2 && x1 <= right && y1 >= low)
+        {
+          if (x1 >= left && y1 <= up) pts[i++].set (x1, y1);
+          y1 --;
+          e -= dx;
+          if (e < 0)
+          {
+            x1 ++;
+            e -= dy;
+          }
+        }
+        if (x2 <= right && y2 >= low) pts[i].set (x2, y2);
+      }
     }
   }
+  *n = i;
   return (pts);
 }
 
@@ -435,9 +656,9 @@ Pt2i *Pt2i::pathTo (Pt2i p, int *n) const
         {
           y1 ++;
           e += dx;
-          pts[i++] = Pt2i (delta, delta);
+          pts[i++].set (delta, delta);
         }
-        else pts[i++] = Pt2i (delta, 0);
+        else pts[i++].set (delta, 0);
       }
     }
 
@@ -458,9 +679,9 @@ Pt2i *Pt2i::pathTo (Pt2i p, int *n) const
         {
           x1 ++;
           e += dy;
-          pts[i++] = Pt2i (delta, delta);
+          pts[i++].set (delta, delta);
         }
-        else pts[i++] = Pt2i (0, delta);
+        else pts[i++].set (0, delta);
       }
     }
   }
@@ -484,9 +705,9 @@ Pt2i *Pt2i::pathTo (Pt2i p, int *n) const
         {
           y1 --;
           e += dx;
-          pts[i++] = Pt2i (delta, -delta);
+          pts[i++].set (delta, -delta);
         }
-        else pts[i++] = Pt2i (delta, 0);
+        else pts[i++].set (delta, 0);
       }
     }
     // Octant 7
@@ -506,9 +727,9 @@ Pt2i *Pt2i::pathTo (Pt2i p, int *n) const
         {
           x1 ++;
           e -= dy;
-          pts[i++] = Pt2i (delta, -delta);
+          pts[i++].set (delta, -delta);
         }
-        else pts[i++] = Pt2i (0, -delta);
+        else pts[i++].set (0, -delta);
       }
     }
   }
diff --git a/Code/FBSD/ImageTools/pt2i.h b/Code/FBSD/ImageTools/pt2i.h
index 466e45f5543fa4dd6da50bc6802b45eb9a74a902..c08d4828c309e1bb5f91d9502d68fdcf4cb1a13f 100755
--- a/Code/FBSD/ImageTools/pt2i.h
+++ b/Code/FBSD/ImageTools/pt2i.h
@@ -12,7 +12,6 @@ using namespace std;
 /** 
  * @class Pt2i pt2i.h
  * \brief Point in the digital plane.
- * \author {P. Even}
  */
 class Pt2i
 {
@@ -86,7 +85,7 @@ public:
   }
 
   /**
-   * @fn void set (const Pt2i &p y)
+   * @fn void set (const Pt2i &p)
    * \brief Sets the point coordinates.
    * @param point to recopy.
    */
@@ -197,24 +196,84 @@ public:
    */
   AbsRat triangleYRationalHeight (const Pt2i &p1, const Pt2i &p2) const;
 
+  /**
+   * @fn bool quadCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+   * \brief Checks if the quadrilater is well formed.
+   * @param p1 Second quadrilater vertex.
+   * @param p2 Third quadrilater vertex.
+   * @param p3 Fourth quadrilater vertex.
+   */
+  bool quadCheck (Pt2i p1, Pt2i p2, Pt2i p3) const;
+
+  /**
+   * @fn bool orientCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+   * \brief Checks if the vector to p1 is oriented as p2 to p3.
+   * @param p1 Second quadrilater vertex.
+   * @param p2 Third quadrilater vertex.
+   * @param p3 Fourth quadrilater vertex.
+   */
+  bool orientCheck (Pt2i p1, Pt2i p2, Pt2i p3) const;
+
+  /**
+   * @fn bool convexCheck (Pt2i p1, Pt2i p2, Pt2i p3) const
+   * \brief Checks if the quadrilater is convex.
+   * @param p1 Second quadrilater vertex.
+   * @param p2 Third quadrilater vertex.
+   * @param p3 Fourth quadrilater vertex.
+   */
+  bool convexCheck (Pt2i p1, Pt2i p2, Pt2i p3) const;
+
+  /**
+   * @fn bool inTriangle (Pt2i p1, Pt2i p2, Pt2i p3) const
+   * \brief Check if point belongs to the given triangle.
+   * @param p1 First triangle vertex.
+   * @param p2 Second triangle vertex.
+   * @param p3 Third triangle vertex.
+   */
+  bool inTriangle (Pt2i p1, Pt2i p2, Pt2i p3) const;
+
+  /**
+   * @fn bool inQuad (Pt2i p1, Pt2i p2, Pt2i p3, Pt2i p4) const
+   * \brief Check if point belongs to the given convex quadrilater.
+   * @param p1 First quadrilater vertex.
+   * @param p2 Second quadrilater vertex.
+   * @param p3 Third quadrilater vertex.
+   * @param p4 Fourth quadrilater vertex.
+   */
+  bool inQuad (Pt2i p1, Pt2i p2, Pt2i p3, Pt2i p4) const;
+
   /**
    * @fn Vr2i vectorTo (Pt2i p)
-   * \brief Returns the discrete vector to the given point.
-   * @param p the given point.
+   * \brief Returns the vector to given point.
+   * @param p Given point.
    */
   Vr2i vectorTo (Pt2i p) const;
 
   /**
    * @fn Pt2i *drawing (Pt2i p, int *n)
-   * \brief Returns the segment to the given point.
-   * @param p the given point.
-   * @param n size of the returned array.
+   * \brief Returns the straight segment to given point.
+   * @param p Given point.
+   * @param n Size of returned array.
    */
   Pt2i *drawing (const Pt2i p, int *n) const;
 
   /**
-   * @fn Pt2i *drawing (Pt2i p, int *n)
-   * \brief Fills the provided vector with the segment to the given point.
+   * @fn Pt2i *clipLine (Pt2i p, int left, int low, int right, int up, int *n)
+   * \brief Returns the clipped straight segment to given point.
+   *   NB: Always returns a non-null array.
+   * @param p Given point.
+   * @param left Position of left clip bound.
+   * @param low Position of lower clip bound.
+   * @param right Position of right clip bound.
+   * @param up Position of upper clip bound.
+   * @param n Size of returned array.
+   */
+  Pt2i *clipLine (const Pt2i p, int left, int low, int right, int up,
+                  int *n) const;
+
+  /**
+   * @fn void draw (vector<Pt2i> &line, Pt2i p)
+   * \brief Fills provided vector with the straight segment to given point.
    * @param line Vector of points to fill in.
    * @param p Distant point.
    */
@@ -222,26 +281,25 @@ public:
 
   /**
    * @fn Pt2i *pathTo (Pt2i p, int *n)
-   * \brief Returns the succession of local moves to reach the given point.
-   * A local move is the vector from a point to one of his neighbours.
-   * @param p the given point.
-   * @param n size of the returned array.
+   * \brief Returns the position of steps to go straight to given point.
+   * @param p Given point.
+   * @param n Size of returned array.
    */
   Pt2i *pathTo (Pt2i p, int *n) const;
 
   /**
    * @fn bool *stepsTo (Pt2i p, int *n)
-   * \brief Returns the location of the steps to reach the given point.
-   * @param p the given point.
-   * @param n size of the returned array.
+   * \brief Returns the location of steps to go straight to given point.
+   * @param p Given point.
+   * @param n Size of returned array.
    */
   bool *stepsTo (Pt2i p, int *n) const;
 
   /**
    * @fn vector<Pt2i> drawOrtho (const Pt2i p2, int offset) const
-   * \brief Returns the segment orthogonal to the segment to p2.
-   * @param p2 The distant point.
-   * @param offset The orthogonal offset.
+   * \brief Returns an orthogonal segment to straight line joining given point.
+   * @param p2 Given point.
+   * @param offset Orthogonal offset.
    */
   vector<Pt2i> drawOrtho (const Pt2i p2, int offset) const;
 
diff --git a/Code/FBSD/ImageTools/strucel.cpp b/Code/FBSD/ImageTools/strucel.cpp
index 6e11051e98b30f8d128a95268a947e33acb53d95..403ef4981e0d909489c57d8e32afc9635c1f1420 100755
--- a/Code/FBSD/ImageTools/strucel.cpp
+++ b/Code/FBSD/ImageTools/strucel.cpp
@@ -17,11 +17,11 @@ Strucel::Strucel (int type)
     height = 3;
     size = 5;
     pattern = new Vr2i[5];
-    pattern[0] = Vr2i (0, 0);
-    pattern[1] = Vr2i (0, 1);
-    pattern[2] = Vr2i (-1, 0);
-    pattern[3] = Vr2i (1, 0);
-    pattern[4] = Vr2i (0, -1);
+    pattern[0].set (0, 0);
+    pattern[1].set (0, 1);
+    pattern[2].set (-1, 0);
+    pattern[3].set (1, 0);
+    pattern[4].set (0, -1);
   }
   else if (type == TYPE_HOR)
   {
@@ -29,9 +29,9 @@ Strucel::Strucel (int type)
     height = 1;
     size = 3;
     pattern = new Vr2i[3];
-    pattern[0] = Vr2i (0, 0);
-    pattern[1] = Vr2i (0, 1);
-    pattern[2] = Vr2i (0, -1);
+    pattern[0].set (0, 0);
+    pattern[1].set (0, 1);
+    pattern[2].set (0, -1);
   }
   else if (type == TYPE_VER)
   {
@@ -39,9 +39,9 @@ Strucel::Strucel (int type)
     height = 3;
     size = 3;
     pattern = new Vr2i[3];
-    pattern[0] = Vr2i (0, 0);
-    pattern[1] = Vr2i (1, 0);
-    pattern[2] = Vr2i (-1, 0);
+    pattern[0].set (0, 0);
+    pattern[1].set (1, 0);
+    pattern[2].set (-1, 0);
   }
   else
   {
@@ -49,7 +49,7 @@ Strucel::Strucel (int type)
     height = 1;
     size = 1;
     pattern = new Vr2i[1];
-    pattern[0] = Vr2i (1, 1);
+    pattern[0].set (1, 1);
   }
 }
 
@@ -60,9 +60,9 @@ Strucel::~Strucel ()
 }
 
 
-void Strucel::tophatGradient (int *out, int **in, int width, int height)
+void Strucel::tophatGradient (int *out, unsigned char *in,
+                              int width, int height)
 {
-// cout << "TH IN" << endl;
   for (int j = 0; j < height; j++)
     for (int i = 0; i < width; i++)
     {
@@ -72,11 +72,10 @@ void Strucel::tophatGradient (int *out, int **in, int width, int height)
         int x = i - pattern[k].x ();
         int y = j - pattern[k].y ();
         if (x >= 0 && x < width && y >= 0 && y < height)
-          if (min == -1 || in[y][x] < min) min = in[y][x];
+          if (min == -1 || in[y * width + x] < min) min = in[y * width + x];
       }
-      out[j * width + i] = in[j][i] - min;
+      out[j * width + i] = in[j * width + i] - min;
     }
-// cout << "TH OUT" << endl;
 }
 
 
@@ -98,9 +97,27 @@ void Strucel::tophatGradient (int *out, int *in, int width, int height)
 }
 
 
-void Strucel::blackhatGradient (int *out, int **in, int width, int height)
+void Strucel::tophatGradient (int *out, int **in, int width, int height)
+{
+  for (int j = 0; j < height; j++)
+    for (int i = 0; i < width; i++)
+    {
+      int min = -1;
+      for (int k = 0; k < size; k++)
+      {
+        int x = i - pattern[k].x ();
+        int y = j - pattern[k].y ();
+        if (x >= 0 && x < width && y >= 0 && y < height)
+          if (min == -1 || in[y][x] < min) min = in[y][x];
+      }
+      out[j * width + i] = in[j][i] - min;
+    }
+}
+
+
+void Strucel::blackhatGradient (int *out, unsigned char *in,
+                                int width, int height)
 {
-// cout << "BH IN" << endl;
   for (int j = 0; j < height; j++)
     for (int i = 0; i < width; i++)
     {
@@ -110,11 +127,10 @@ void Strucel::blackhatGradient (int *out, int **in, int width, int height)
         int x = i - pattern[k].x ();
         int y = j - pattern[k].y ();
         if (x >= 0 && x < width && y >= 0 && y < height)
-          if (in[y][x] > max) max = in[y][x];
+          if (in[y * width + x] > max) max = in[y * width + x];
       }
-      out[j * width + i] = max - in[j][i];
+      out[j * width + i] = max - in[j * width + i];
     }
-// cout << "BH OUT" << endl;
 }
 
 
@@ -136,7 +152,26 @@ void Strucel::blackhatGradient (int *out, int *in, int width, int height)
 }
 
 
-void Strucel::morphoGradient (int *out, int **in, int width, int height)
+void Strucel::blackhatGradient (int *out, int **in, int width, int height)
+{
+  for (int j = 0; j < height; j++)
+    for (int i = 0; i < width; i++)
+    {
+      int max = 0;
+      for (int k = 0; k < size; k++)
+      {
+        int x = i - pattern[k].x ();
+        int y = j - pattern[k].y ();
+        if (x >= 0 && x < width && y >= 0 && y < height)
+          if (in[y][x] > max) max = in[y][x];
+      }
+      out[j * width + i] = max - in[j][i];
+    }
+}
+
+
+void Strucel::morphoGradient (int *out, unsigned char *in,
+                              int width, int height)
 {
   for (int j = 0; j < height; j++)
     for (int i = 0; i < width; i++)
@@ -149,8 +184,8 @@ void Strucel::morphoGradient (int *out, int **in, int width, int height)
         int y = j - pattern[k].y ();
         if (x >= 0 && x < width && y >= 0 && y < height)
         {
-          if (min == -1 || in[y][x] < min) min = in[y][x];
-          if (in[y][x] > max) max = in[y][x];
+          if (min == -1 || in[y * width + x] < min) min = in[y * width + x];
+          if (in[y * width + x] > max) max = in[y * width + x];
         }
       }
       out[j * width + i] = max - min;
@@ -178,3 +213,25 @@ void Strucel::morphoGradient (int *out, int *in, int width, int height)
       out[j * width + i] = max - min;
     }
 }
+
+
+void Strucel::morphoGradient (int *out, int **in, int width, int height)
+{
+  for (int j = 0; j < height; j++)
+    for (int i = 0; i < width; i++)
+    {
+      int max = 0;
+      int min = -1;
+      for (int k = 0; k < size; k++)
+      {
+        int x = i - pattern[k].x ();
+        int y = j - pattern[k].y ();
+        if (x >= 0 && x < width && y >= 0 && y < height)
+        {
+          if (min == -1 || in[y][x] < min) min = in[y][x];
+          if (in[y][x] > max) max = in[y][x];
+        }
+      }
+      out[j * width + i] = max - min;
+    }
+}
diff --git a/Code/FBSD/ImageTools/strucel.h b/Code/FBSD/ImageTools/strucel.h
index 00844a5b228bf299d2f5ac92953e7acb3db1477c..9482e7cd49a6b8ce44eeaa76193b6bb68503d0e3 100755
--- a/Code/FBSD/ImageTools/strucel.h
+++ b/Code/FBSD/ImageTools/strucel.h
@@ -44,7 +44,6 @@ public:
 
   /** 
    * \brief Returns the pattern height.
-
    */
   inline int getHeight () const
   {
@@ -54,11 +53,11 @@ public:
   /** 
    * \brief Calculates the top hat gradient of an image.
    * @param out Output gradient array (the memory should be allocated before).
-   * @param in Image bi-dimensional array.
+   * @param in Image array.
    * @param width Image width.
    * @param height Image height.
    */
-  void tophatGradient (int *out, int **in, int width, int height);
+  void tophatGradient (int *out, unsigned char *in, int width, int height);
 
   /** 
    * \brief Calculates the top hat gradient of an image.
@@ -70,13 +69,22 @@ public:
   void tophatGradient (int *out, int *in, int width, int height);
 
   /** 
-   * \brief Calculates the black hat gradient of an image.
+   * \brief Calculates the top hat gradient of an image.
    * @param out Output gradient array (the memory should be allocated before).
    * @param in Image bi-dimensional array.
    * @param width Image width.
    * @param height Image height.
    */
-  void blackhatGradient (int *out, int **in, int width, int height);
+  void tophatGradient (int *out, int **in, int width, int height);
+
+  /** 
+   * \brief Calculates the black hat gradient of an image.
+   * @param out Output gradient array (the memory should be allocated before).
+   * @param in Image array.
+   * @param width Image width.
+   * @param height Image height.
+   */
+  void blackhatGradient (int *out, unsigned char *in, int width, int height);
 
   /** 
    * \brief Calculates the black hat gradient of an image.
@@ -88,13 +96,22 @@ public:
   void blackhatGradient (int *out, int *in, int width, int height);
 
   /** 
-   * \brief Calculates the morphological gradient of an image.
+   * \brief Calculates the black hat gradient of an image.
    * @param out Output gradient array (the memory should be allocated before).
    * @param in Image bi-dimensional array.
    * @param width Image width.
    * @param height Image height.
    */
-  void morphoGradient (int *out, int **in, int width, int height);
+  void blackhatGradient (int *out, int **in, int width, int height);
+
+  /** 
+   * \brief Calculates the morphological gradient of an image.
+   * @param out Output gradient array (the memory should be allocated before).
+   * @param in Image array.
+   * @param width Image width.
+   * @param height Image height.
+   */
+  void morphoGradient (int *out, unsigned char *in, int width, int height);
 
   /** 
    * \brief Calculates the morphological gradient of an image.
@@ -105,6 +122,15 @@ public:
    */
   void morphoGradient (int *out, int *in, int width, int height);
 
+  /** 
+   * \brief Calculates the morphological gradient of an image.
+   * @param out Output gradient array (the memory should be allocated before).
+   * @param in Image bi-dimensional array.
+   * @param width Image width.
+   * @param height Image height.
+   */
+  void morphoGradient (int *out, int **in, int width, int height);
+
 
 private:
 
diff --git a/Code/FBSD/ImageTools/vmap.cpp b/Code/FBSD/ImageTools/vmap.cpp
index 7f5afed6d04e08c2389aa748ddd676007fd45308..42cb4a1155dad6c5b73d1101784d5cb2a7fd4bd2 100755
--- a/Code/FBSD/ImageTools/vmap.cpp
+++ b/Code/FBSD/ImageTools/vmap.cpp
@@ -1,10 +1,12 @@
 // #include <iostream>
 #include "vmap.h"
 #include "math.h"
+#include <inttypes.h>
 
 using namespace std;
 
 
+const int VMap::TYPE_UNKNOWN = -1;
 const int VMap::TYPE_SOBEL_3X3 = 0;
 const int VMap::TYPE_SOBEL_5X5 = 1;
 const int VMap::TYPE_TOP_HAT = 2;
@@ -23,6 +25,102 @@ const int VMap::NB_DILATIONS = 5;
 const int VMap::DEFAULT_DILATION = 4;
 
 
+VMap::VMap (int width, int height, unsigned char *data, int type)
+{
+  this->width = width;
+  this->height = height;
+  this->gtype = type;
+  init ();
+  imap = new int[width * height];
+  if (type == TYPE_TOP_HAT)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.tophatGradient (imap, data, width, height);
+    buildSobel5x5Map (data);
+  }
+  else if (type == TYPE_FULL_TOP_HAT)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.tophatGradient (imap, data, width, height);
+    int *jmap = new int[width * height];
+    Strucel seh (Strucel::TYPE_HOR);
+    seh.tophatGradient (jmap, data, width, height);
+    int *kmap = new int[width * height];
+    Strucel sev (Strucel::TYPE_VER);
+    sev.tophatGradient (kmap, data, width, height);
+    map = new Vr2i[width * height];
+    Vr2i *tmpmap = map;
+    for (int i = 0; i < width * height; i ++)
+    {
+      tmpmap->set (jmap[i], kmap[i]);
+      tmpmap ++;
+    }
+  }
+  else if (type == TYPE_BLACK_HAT)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.blackhatGradient (imap, data, width, height);
+    buildSobel5x5Map (data);
+  }
+  else if (type == TYPE_FULL_BLACK_HAT)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.blackhatGradient (imap, data, width, height);
+    int *jmap = new int[width * height];
+    Strucel seh (Strucel::TYPE_HOR);
+    seh.blackhatGradient (jmap, data, width, height);
+    int *kmap = new int[width * height];
+    Strucel sev (Strucel::TYPE_VER);
+    sev.blackhatGradient (kmap, data, width, height);
+    map = new Vr2i[width * height];
+    Vr2i *tmpmap = map;
+    for (int i = 0; i < width * height; i ++)
+    {
+      tmpmap->set (jmap[i], kmap[i]);
+      tmpmap ++;
+    }
+  }
+  else if (type == TYPE_MORPHO)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.morphoGradient (imap, data, width, height);
+    buildSobel5x5Map (data);
+  }
+  else if (type == TYPE_FULL_MORPHO)
+  {
+    Strucel se (Strucel::TYPE_PLUS_3X3);
+    se.morphoGradient (imap, data, width, height);
+    int *jmap = new int[width * height];
+    Strucel seh (Strucel::TYPE_HOR);
+    seh.morphoGradient (jmap, data, width, height);
+    int *kmap = new int[width * height];
+    Strucel sev (Strucel::TYPE_VER);
+    sev.morphoGradient (kmap, data, width, height);
+    map = new Vr2i[width * height];
+    Vr2i *tmpmap = map;
+    for (int i = 0; i < width * height; i ++)
+    {
+      tmpmap->set (jmap[i], kmap[i]);
+      tmpmap ++;
+    }
+  }
+  else if (type == TYPE_SOBEL_5X5)
+  {
+    buildSobel5x5Map (data);
+    for (int i = 0; i < width * height; i++)
+      imap[i] = (int) sqrt (map[i].norm2 ());
+    gmagThreshold *= gradientThreshold;
+  }
+  else if (type == TYPE_SOBEL_3X3)
+  {
+    buildGradientMap (data);
+    for (int i = 0; i < width * height; i++)
+      imap[i] = (int) sqrt (map[i].norm2 ());
+    gmagThreshold *= gradientThreshold;
+  }
+}
+
+
 VMap::VMap (int width, int height, int *data, int type)
 {
   this->width = width;
@@ -215,6 +313,20 @@ VMap::VMap (int width, int height, int **data, int type)
 }
 
 
+VMap::VMap (int width, int height, Vr2i *map)
+{
+  this->width = width;
+  this->height = height;
+  this->gtype = TYPE_UNKNOWN;
+  this->map = map;
+  init ();
+  imap = new int[width * height];
+  for (int i = 0; i < width * height; i++)
+    imap[i] = (int) sqrt (map[i].norm2 ());
+  gmagThreshold *= gradientThreshold;
+}
+
+
 VMap::~VMap ()
 {
   delete [] map;
@@ -236,26 +348,26 @@ void VMap::init ()
   angleThreshold = NEAR_SQ_ANGLE;
   orientedGradient = true;
   bowl = new Vr2i[MAX_BOWL];
-  bowl[0] = Vr2i (1, 0);
-  bowl[1] = Vr2i (0, 1);
-  bowl[2] = Vr2i (-1, 0);
-  bowl[3] = Vr2i (0, -1);
-  bowl[4] = Vr2i (1, 1);
-  bowl[5] = Vr2i (1, -1);
-  bowl[6] = Vr2i (-1, -1);
-  bowl[7] = Vr2i (-1, 1);
-  bowl[8] = Vr2i (2, 0);
-  bowl[9] = Vr2i (0, 2);
-  bowl[10] = Vr2i (-2, 0);
-  bowl[11] = Vr2i (0, -2);
-  bowl[12] = Vr2i (2, 1);
-  bowl[13] = Vr2i (1, 2);
-  bowl[14] = Vr2i (-1, 2);
-  bowl[15] = Vr2i (-2, 1);
-  bowl[16] = Vr2i (-2, -1);
-  bowl[17] = Vr2i (-1, -2);
-  bowl[18] = Vr2i (1, -2);
-  bowl[19] = Vr2i (2, -1);
+  bowl[0].set (1, 0);
+  bowl[1].set (0, 1);
+  bowl[2].set (-1, 0);
+  bowl[3].set (0, -1);
+  bowl[4].set (1, 1);
+  bowl[5].set (1, -1);
+  bowl[6].set (-1, -1);
+  bowl[7].set (-1, 1);
+  bowl[8].set (2, 0);
+  bowl[9].set (0, 2);
+  bowl[10].set (-2, 0);
+  bowl[11].set (0, -2);
+  bowl[12].set (2, 1);
+  bowl[13].set (1, 2);
+  bowl[14].set (-1, 2);
+  bowl[15].set (-2, 1);
+  bowl[16].set (-2, -1);
+  bowl[17].set (-1, -2);
+  bowl[18].set (1, -2);
+  bowl[19].set (2, -1);
   dilations = new int[NB_DILATIONS];
   dilations[0] = 0;
   dilations[1] = 4;
@@ -266,6 +378,47 @@ void VMap::init ()
 }
 
 
+void VMap::buildGradientMap (unsigned char *data)
+{
+  map = new Vr2i[width * height];
+  Vr2i *gm = map;
+
+  for (int j = 0; j < width; j++)
+  {
+    gm->set (0, 0);
+    gm++;
+  }
+  for (int i = 1; i < height - 1; i++)
+  {
+    gm->set (0, 0);
+    gm++;
+    for (int j = 1; j < width - 1; j++)
+    {
+      gm->set (data[(i - 1) * height + j + 1]
+               + 2 * data[i * height + j + 1]
+               + data[(i + 1) * height + j + 1]
+               - data[(i - 1) * height + j - 1]
+               - 2 * data[i * height + j - 1]
+               - data[(i + 1) * height + j - 1],
+               data[(i + 1) * height + j - 1]
+               + 2 * data[(i + 1) * height + j]
+               + data[(i + 1) * height + j + 1]
+               - data[(i - 1) * height + j - 1]
+               - 2 * data[(i - 1) * height + j]
+               - data[(i - 1) * height + j + 1]);
+      gm++;
+    }
+    gm->set (0, 0);
+    gm++;
+  }
+  for (int j = 0; j < width; j++)
+  {
+    gm->set (0, 0);
+    gm++;
+  }
+}
+
+
 void VMap::buildGradientMap (int *data)
 {
   map = new Vr2i[width * height];
@@ -340,6 +493,79 @@ void VMap::buildGradientMap (int **data)
 }
 
 
+void VMap::buildSobel5x5Map (unsigned char *data)
+{
+  map = new Vr2i[width * height];
+  Vr2i *gm = map;
+
+  for (int j = 0; j < 2 * width; j++)
+  {
+    gm->set (0, 0);
+    gm++;
+  }
+  for (int i = 2; i < height - 2; i++)
+  {
+    gm->set (0, 0);
+    gm++;
+    gm->set (0, 0);
+    gm++;
+    for (int j = 2; j < width - 2; j++)
+    {
+      gm->set (5 * data[(i - 2) * height + j + 2]
+                 + 8 * data[(i - 1) * height + j + 2]
+                 + 10 * data[i * height + j + 2]
+                 + 8 * data[(i + 1) * height + j + 2]
+                 + 5 * data[(i + 2) * height + j + 2]
+               + 4 * data[(i - 2) * height + j + 1]
+                 + 10 * data[(i - 1) * height + j + 1]
+                 + 20 * data[i * height + j + 1]
+                 + 10 * data[(i + 1) * height + j + 1]
+                 + 4 * data[(i + 2) * height + j + 1]
+               - 4 * data[(i - 2) * height + j - 1]
+                 - 10 * data[(i - 1) * height + j - 1]
+                 - 20 * data[i * height + j - 1]
+                 - 10 * data[(i + 1) * height + j - 1]
+                 - 4 * data[(i + 2) * height + j - 1] 
+               - 5 * data[(i - 2) * height + j - 2]
+                 - 8 * data[(i - 1) * height + j - 2]
+                 - 10 * data[i * height + j - 2]
+                 - 8 * data[(i + 1) * height + j - 2]
+                 - 5 * data[(i + 2) * height + j - 2],
+               5 * data[(i + 2) * height + j - 2]
+                 + 8 * data[(i + 2) * height + j - 1]
+                 + 10 * data[(i + 2) * height + j]
+                 + 8 * data[(i + 2) * height + j + 1]
+                 + 5 * data[(i + 2) * height + j + 2]
+               + 4 * data[(i + 1) * height + j - 2]
+                 + 10 * data[(i + 1) * height + j - 1]
+                 + 20 * data[(i + 1) * height + j]
+                 + 10 * data[(i + 1) * height + j + 1]
+                 + 4 * data[(i + 1) * height + j + 2]
+               - 4 * data[(i - 1) * height + j - 2]
+                 - 10 * data[(i - 1) * height + j - 1]
+                 - 20 * data[(i - 1) * height + j]
+                 - 10 * data[(i - 1) * height + j + 1]
+                 - 4 * data[(i - 1) * height + j + 2]
+               - 5 * data[(i - 2) * height + j - 2]
+                 - 8 * data[(i - 2) * height + j - 1]
+                 - 10 * data[(i - 2) * height + j]
+                 - 8 * data[(i - 2) * height + j + 1]
+                 - 5 * data[(i - 2) * height + j + 2]);
+      gm++;
+    }
+    gm->set (0, 0);
+    gm++;
+    gm->set (0, 0);
+    gm++;
+  }
+  for (int j = 0; j < 2 * width; j++)
+  {
+    gm->set (0, 0);
+    gm++;
+  }
+}
+
+
 void VMap::buildSobel5x5Map (int *data)
 {
   map = new Vr2i[width * height];
@@ -514,8 +740,9 @@ int VMap::keepFreeElementsIn (const vector<Pt2i> &pix, int n, int *ind) const
 
 int VMap::keepContrastedMax (int *lmax, int n, int *in) const
 {
-  int min[n-1];
-  bool fired[n];
+  if (n == 0) return 0;
+  int *min = new int[n-1];
+  bool *fired = new bool[n];
   int nbfired = 0;
   int sleft = 0;
 
@@ -584,17 +811,17 @@ int VMap::keepDirectedElementsAs (const vector<Pt2i> &pix,
 int VMap::keepOrientedElementsAs (const vector<Pt2i> &pix,
                                   int n, int *ind, const Vr2i &ref) const
 {
-  int vx = ref.x ();
-  int vy = ref.y ();
-  long vn2 = vx * vx + vy * vy;
+  int64_t vx = (int64_t) ref.x ();
+  int64_t vy = (int64_t) ref.y ();
+  int64_t vn2 = (int64_t) vx * (int64_t) vx + (int64_t) vy * (int64_t) vy;
 
   int i = 0;
   while (i < n)
   {
     Pt2i p = pix.at (ind[i]);
     Vr2i gr = map[p.y () * width + p.x ()];
-    long gx = gr.x ();
-    long gy = gr.y ();
+    int64_t gx = (int64_t) gr.x ();
+    int64_t gy = (int64_t) gr.y ();
     if ((vx * vx * gx * gx + vy * vy * gy * gy + 2 * vx * vy * gx * gy) * 100
         < vn2 * (gx * gx + gy * gy) * angleThreshold)
       ind[i] = ind[--n];
diff --git a/Code/FBSD/ImageTools/vmap.h b/Code/FBSD/ImageTools/vmap.h
index aa733e773834f53b491317bfa09c9927b7c27554..34555300776e42df293d862c5e47985d4a638c57 100755
--- a/Code/FBSD/ImageTools/vmap.h
+++ b/Code/FBSD/ImageTools/vmap.h
@@ -16,6 +16,8 @@ class VMap
 {
 public:
 
+  /** Gradient extraction method : Undeterminated. */
+  static const int TYPE_UNKNOWN;
   /** Gradient extraction method : Sobel with 3x3 kernel. */
   static const int TYPE_SOBEL_3X3;
   /** Gradient extraction method : Sobel with 5x5 kernel. */
@@ -39,7 +41,16 @@ public:
    * @param width Map width.
    * @param height Map height.
    * @param data Scalar data array.
-   * @param type Gradient extraction method (default is Soble with 3x3 kernel).
+   * @param type Gradient extraction method (default is Sobel with 3x3 kernel).
+   */
+  VMap (int width, int height, unsigned char *data, int type = 0);
+
+  /** 
+   * \brief Creates a gradient map from scalar data.
+   * @param width Map width.
+   * @param height Map height.
+   * @param data Scalar data array.
+   * @param type Gradient extraction method (default is Sobel with 3x3 kernel).
    */
   VMap (int width, int height, int *data, int type = 0);
 
@@ -51,6 +62,14 @@ public:
    */
   VMap (int width, int height, int **data, int type = 0);
 
+  /** 
+   * \brief Creates a gradient map from given vector map.
+   * @param width Map width.
+   * @param height Map height.
+   * @param map Vector map.
+   */
+  VMap (int width, int height, Vr2i *map);
+
   /** 
    * \brief Deletes the vector map.
    */
@@ -95,19 +114,20 @@ public:
    * @param i Column number.
    * @param j Line number.
    */
-  Vr2i getValue (int i, int j) const
-  {
-    return (map[j * width + i]);
-  }
+  inline Vr2i *getVectorMap () const { return (map); }
+
+  /**
+   * \brief Returns the vector at point (i,j).
+   * @param i Column number.
+   * @param j Line number.
+   */
+  inline Vr2i getValue (int i, int j) const { return (map[j * width + i]); }
 
   /**
    * \brief Returns the vector at given point.
    * @param p Point position.
    */
-  Vr2i getValue (Pt2i p) const
-  {
-    return (map[p.y () * width + p.x ()]);
-  }
+  inline Vr2i getValue (Pt2i p) const { return (map[p.y () * width + p.x ()]); }
 
   /**
    * \brief Returns the squared norm of the vector magnitude at point (i,j).
@@ -339,6 +359,13 @@ private:
    */
   void init ();
 
+  /** 
+   * \brief Builds the vector map as a gradient map from provided data.
+   * Uses a Sobel 3x3 kernel.
+   * @param data Initial scalar data.
+   */
+  void buildGradientMap (unsigned char *data);
+
   /** 
    * \brief Builds the vector map as a gradient map from provided data.
    * Uses a Sobel 3x3 kernel.
@@ -353,6 +380,13 @@ private:
    */
   void buildGradientMap (int **data);
 
+  /** 
+   * \brief Builds the vector map as a gradient map from provided data.
+   * Uses a Sobel 5x5 kernel.
+   * @param data Initial scalar data.
+   */
+  void buildSobel5x5Map (unsigned char *data);
+
   /** 
    * \brief Builds the vector map as a gradient map from provided data.
    * Uses a Sobel 5x5 kernel.
diff --git a/Code/FBSD/main.cpp b/Code/FBSD/main.cpp
index 5bececd1c2be3ed70116f2b59903a0a1422afd5c..37f8df491ed9d8757edabd853a7b2fb91102b8f8 100755
--- a/Code/FBSD/main.cpp
+++ b/Code/FBSD/main.cpp
@@ -3,7 +3,6 @@
 #include <cmath>
 #include "bswindow.h"
 #include "bsrandomtester.h"
-// #include "scanwindow.h"  // DEV
 
 
 int main (int argc, char *argv[])
@@ -16,18 +15,6 @@ int main (int argc, char *argv[])
   bool out = false;
   QApplication app (argc, argv);
 
-// DEV IN
-/*
-  if (argc == 2 && string (argv[1]) == string ("scans"))
-  {
-    ScanWindow sw (&val);
-    sw.setWindowTitle ("Scans");
-    sw.show ();
-    return app.exec ();
-  }
-*/
-// DEV OUT
-
   BSWindow window (&val);   // val : necessary argument !
   for (int i = 1; i < argc; i++)
   {