From 0aea97d82c70dcc54b581197c9624833b3a4d091 Mon Sep 17 00:00:00 2001
From: gdd <gdd>
Date: Thu, 31 Mar 2011 14:47:35 +0000
Subject: [PATCH] rnc: The size of the the grid in the discretized parametric
 space of the surface to be meshed is now calculated in a different way. The
 number of grid points is now constant.

---
 src/BLSURFPlugin/BLSURFPlugin_Attractor.cxx  | 77 +++++++++++++++-----
 src/BLSURFPlugin/BLSURFPlugin_Hypothesis.cxx |  2 +-
 2 files changed, 59 insertions(+), 20 deletions(-)

diff --git a/src/BLSURFPlugin/BLSURFPlugin_Attractor.cxx b/src/BLSURFPlugin/BLSURFPlugin_Attractor.cxx
index ad23179..0f9478f 100644
--- a/src/BLSURFPlugin/BLSURFPlugin_Attractor.cxx
+++ b/src/BLSURFPlugin/BLSURFPlugin_Attractor.cxx
@@ -1,3 +1,34 @@
+//  Copyright (C) 2007-2010  CEA/DEN, EDF R&D
+//
+//  This library is free software; you can redistribute it and/or
+//  modify it under the terms of the GNU Lesser General Public
+//  License as published by the Free Software Foundation; either
+//  version 2.1 of the License.
+//
+//  This library is distributed in the hope that it will be useful,
+//  but WITHOUT ANY WARRANTY; without even the implied warranty of
+//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+//  Lesser General Public License for more details.
+//
+//  You should have received a copy of the GNU Lesser General Public
+//  License along with this library; if not, write to the Free Software
+//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+//
+//  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+
+// ---
+// File    : BLSURFPlugin_Attractor.cxx
+// Authors : Renaud Nédélec (OCC)
+// ---
+// 
+// The idea of the algorithm used to calculate the distance on a 
+// non-euclidian parametric surface has been found in the ref. below:
+//
+// Ref:"Accurate Anisotropic Fast Marching for Diffusion-Based Geodesic Tractography"
+// S. Jbabdi, P. Bellec, R. Toro, Daunizeau, M. Pélégrini-Issac, and H. Benali1
+//
+
 #include "BLSURFPlugin_Attractor.hxx"
 #include <utilities.h>
 #include <algorithm>
@@ -32,11 +63,11 @@ BLSURFPlugin_Attractor::BLSURFPlugin_Attractor ()
   _isMapBuilt(false),
   _isEmpty(true){ MESSAGE("construction of a void attractor"); }
 
-BLSURFPlugin_Attractor::BLSURFPlugin_Attractor (const TopoDS_Face& Face, const TopoDS_Shape& Attractor, const std::string& attEntry, double Step)
+BLSURFPlugin_Attractor::BLSURFPlugin_Attractor (const TopoDS_Face& Face, const TopoDS_Shape& Attractor, const std::string& attEntry, double Step) // TODO Step is now unused -> remove it if testing is OK
   : _face(),
   _attractorShape(),
   _attEntry(attEntry),
-  _step(Step),
+  _step(),
   _gridU(),
   _gridV(),
   _vectU(),
@@ -74,8 +105,13 @@ bool BLSURFPlugin_Attractor::init(){
   ShapeAnalysis::GetFaceUVBounds(_face,_u1,_u2,_v1,_v2);
   MESSAGE("u1 = "<<_u1<<" ,u2  = "<<_u2);
   MESSAGE("v1 = "<<_v1<<" ,v2  = "<<_v2);
-  _gridU = floor (_u2 - _u1) / _step;
-  _gridV = floor (_v2 - _v1) / _step;
+//   _gridU = floor (_u2 - _u1) / _step;
+//   _gridV = floor (_v2 - _v1) / _step;
+  // TEST
+  _gridU = 300;
+  _gridV = 300;
+  _step = std::min((_u2-_u1)/_gridU,(_v2-_v1)/_gridV);
+
   for (i=0; i<=_gridU; i++){
     _vectU.push_back(_u1+i*(_u2-_u1)/_gridU) ;
   }
@@ -116,7 +152,10 @@ bool BLSURFPlugin_Attractor::init(){
     ShapeConstruct_ProjectCurveOnSurface curveProjector;
     curveProjector.Init(aSurf, Precision::Confusion());
     curveProjector.PerformAdvanced (aCurve3d, first, last, aCurve2d);
-    int N = 20 * (last - first) / _step;  // If the edge is a circle : 4>Pi so the number of points on the edge should be good -> 5 for ellipses
+    // TEST
+    //int N = 20 * (last - first) / _step;  // If the edge is a circle : 4>Pi so the number of points on the edge should be good -> 5 for ellipses
+    int N = 1200;
+    MESSAGE("Initialisation des points de départ")
     for (i=0; i<=N; i++){
       P2 = aCurve2d->Value(first + i * (last-first) / N);
       i0 = floor( (P2.X() - _u1) * _gridU / (_u2 - _u1) + 0.5 );
@@ -161,23 +200,23 @@ double BLSURFPlugin_Attractor::_distance(double u, double v){
 //   
 //   double u_inf = _vectU[i_inf];
 //   double v_inf = _vectV[j_inf];
-//   
-// //   MESSAGE("i_inf , i_sup, j_inf, j_sup = "<<i_inf<<" , "<<i_sup<<" , "<<j_inf<<" , "<<j_sup)
-// //   MESSAGE("u = "<<u<<", _u1 ="<<_u1)
-// //   MESSAGE("(u - _u1) / stepU = "<< (u - _u1) / stepU)
+// //   
+// // //   MESSAGE("i_inf , i_sup, j_inf, j_sup = "<<i_inf<<" , "<<i_sup<<" , "<<j_inf<<" , "<<j_sup)
+// // //   MESSAGE("u = "<<u<<", _u1 ="<<_u1)
+// // //   MESSAGE("(u - _u1) / stepU = "<< (u - _u1) / stepU)
 //   double d1 = _DMap[i_sup][j_sup]; 
 //   double d2 = _DMap[i_sup][j_inf];
 //   double d3 = _DMap[i_inf][j_sup];
 //   double d4 = _DMap[i_inf][j_inf];
-//   
-//   //double d = 0.25 * (d1 + d2 + d3 + d4);
-//   //double d = d1;
-//  //
-// //   if (fabs(v_inf-v_sup) < 1e-14 || fabs(u_inf-u_sup) < 1e-14 ){
-// //     MESSAGE("division par zero v_inf-v_sup = "<< fabs(v_inf-v_sup)<<" , u_inf-u_sup"<<fabs(u_inf-u_sup))
-// //     MESSAGE("v_inf = "<< v_inf<<" , v_sup"<<v_sup)
-// //     MESSAGE("u_inf = "<< u_inf<<" , u_sup"<<u_sup)
-// //   }
+// //   
+// //   double d = 0.25 * (d1 + d2 + d3 + d4);
+// //   //double d = d1;
+// //  //
+// // //   if (fabs(v_inf-v_sup) < 1e-14 || fabs(u_inf-u_sup) < 1e-14 ){
+// // //     MESSAGE("division par zero v_inf-v_sup = "<< fabs(v_inf-v_sup)<<" , u_inf-u_sup"<<fabs(u_inf-u_sup))
+// // //     MESSAGE("v_inf = "<< v_inf<<" , v_sup"<<v_sup)
+// // //     MESSAGE("u_inf = "<< u_inf<<" , u_sup"<<u_sup)
+// // //   }
 //   double P_inf =   d4 * ( 1 + (u_inf - u) / stepU + (v_inf - v) / stepV )
 //                  + d3 * ( (v - v_inf) / stepV )
 //                  + d2 * ( (u - u_inf) / stepU ) ;
@@ -198,7 +237,7 @@ double BLSURFPlugin_Attractor::_distance(double u, double v){
 //   
 //   return d; 
   
-  // BLSURF seems to perform a linear interpolation so it's sufficient to give it a non-continuous distance map
+  //   BLSURF seems to perform a linear interpolation so it's sufficient to give it a non-continuous distance map
   int i = floor ( (u - _u1) * _gridU / (_u2 - _u1) + 0.5 );
   int j = floor ( (v - _v1) * _gridV / (_v2 - _v1) + 0.5 );
   
diff --git a/src/BLSURFPlugin/BLSURFPlugin_Hypothesis.cxx b/src/BLSURFPlugin/BLSURFPlugin_Hypothesis.cxx
index 4e51e65..90dc81d 100644
--- a/src/BLSURFPlugin/BLSURFPlugin_Hypothesis.cxx
+++ b/src/BLSURFPlugin/BLSURFPlugin_Hypothesis.cxx
@@ -413,7 +413,7 @@ void BLSURFPlugin_Hypothesis::SetClassAttractorEntry(const std::string& entry, c
   if ( !attExists || (attExists && _classAttractors[entry]->GetAttractorEntry().compare(attEntry) != 0)){ 
     ShapeAnalysis::GetFaceUVBounds(FaceShape,u1,u2,v1,v2);
     diag = sqrt((u2 - u1) * (u2 - u1) + (v2 - v1) * (v2 - v1));  
-    BLSURFPlugin_Attractor* myAttractor = new BLSURFPlugin_Attractor(FaceShape, AttractorShape, attEntry, 0.002 * diag); 
+    BLSURFPlugin_Attractor* myAttractor = new BLSURFPlugin_Attractor(FaceShape, AttractorShape, attEntry, 0.1 ); // test 0.002 * diag); 
     myAttractor->BuildMap();
     myAttractor->SetParameters(StartSize, EndSize, ActionRadius, ConstantRadius);
     _classAttractors[entry] = myAttractor;
-- 
2.39.2