From f29e84dcc0d313a690d1c78a568cc3659360bd80 Mon Sep 17 00:00:00 2001
From: Jean-Francois Remacle <jean-francois.remacle@uclouvain.be>
Date: Mon, 4 Jul 2005 15:07:41 +0000
Subject: [PATCH] *** empty log message ***

---
 Fltk/Callbacks.cpp          | 116 +++++++++++++++++++++++++++++++++++-
 Fltk/Callbacks.h            |   5 +-
 Fltk/GUI.cpp                | 113 +++++++++++++++++++++++++++++------
 Fltk/GUI.h                  |   8 ++-
 Graphics/Geom.cpp           |  24 ++++++--
 Mesh/BDS.cpp                | 106 +++++++++++++++++++++++++-------
 Mesh/BDS.h                  |  27 +++++----
 Mesh/DiscreteSurface.cpp    |  70 +++++++++++++++++++++-
 Parser/OpenFile.cpp         |  21 ++++---
 benchmarks/stl/ktoolcav.stl | Bin 205083 -> 205307 bytes
 10 files changed, 416 insertions(+), 74 deletions(-)

diff --git a/Fltk/Callbacks.cpp b/Fltk/Callbacks.cpp
index 311f5b4d32..9655ebbd5c 100644
--- a/Fltk/Callbacks.cpp
+++ b/Fltk/Callbacks.cpp
@@ -1,4 +1,4 @@
-// $Id: Callbacks.cpp,v 1.359 2005-07-03 08:02:23 geuzaine Exp $
+// $Id: Callbacks.cpp,v 1.360 2005-07-04 15:07:40 remacle Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -836,13 +836,125 @@ void wizard_update_edges_cb(CALLBACK_ARGS)
     if (THEM && THEM->bds && WID)
     {
 	const double angle = WID->swiz_value[0]->value() * M_PI / 180;
-	printf ("value = %g\n",angle);
 	THEM->bds->classify (angle);
 	BDS_To_Mesh (THEM); 
+	char a[25];
+	sprintf(a,"%d",Tree_Nbr(THEM->Points));
+	WID->swiz_output[1]->value(a);
+	sprintf(a,"%d",Tree_Nbr(THEM->Curves));
+	WID->swiz_output[2]->value(a);
+	sprintf(a,"%d",Tree_Nbr(THEM->Surfaces));
+	WID->swiz_output[3]->value(a);
 	Draw();
     }
 }
 
+void wizard_update_more_edges_cb(CALLBACK_ARGS)
+{
+    Vertex *v;
+    Curve *c;
+    Surface *s;
+    int n,p[100];
+    extern  void BDS_To_Mesh(Mesh *m);
+    if (THEM && THEM->bds && WID)
+    {
+	const double angle = WID->swiz_value[2]->value() * M_PI / 180;
+	THEM->bds->classify (angle);
+	BDS_To_Mesh (THEM); 
+	Draw();
+	n=0;
+	while(1) {
+	    Msg(STATUS3N, "Adding new Model Edges");
+	    if(n == 0)
+		Msg(ONSCREEN, "Select Model Edges\n"
+		    "[Press 'q' to abort or 'e' end]");
+	    if(n == 1)
+		Msg(ONSCREEN, "Select Model Edge\n"
+		    "[Press 'u' to undo last selection, 'q' to abort, 'e' end]");
+	    char ib = SelectEntity(ENT_LINE, &v, &c, &s);
+	    printf("ib = %c\n",ib);
+	    if(ib == 'l') {
+		p[n++] = c->Num;
+		printf("line %d has been selected\n",c->Num);
+	    }
+	    if(ib == 'u') {
+		if(n > 0){
+		    ZeroHighlightEntityNum(p[n-1], 0, 0);
+		    Draw(); // inefficient, but hard to do otherwise
+		    n--;
+		}
+	    }
+	    if(ib == 'q') {
+		ZeroHighlight(THEM);
+		Msg(ONSCREEN, "");
+		Draw(); 
+		Msg(STATUS3N, "Ready");
+		const double angle = WID->swiz_value[0]->value() * M_PI / 180;
+		THEM->bds->classify (angle);
+		BDS_To_Mesh (THEM); 
+		Draw();
+		break;
+	    }
+	    if(ib == 'e') {
+		for (int i=0;i<n;i++)
+		{
+		    BDS_GeomEntity *g = THEM->bds->get_geom(p[i],1);
+		    std::list<BDS_Edge*>::iterator it  = g->e.begin();
+		    std::list<BDS_Edge*>::iterator ite = g->e.end();
+		    while (it!=ite){			
+			BDS_Edge *e = (*it);
+			e->status = 1;
+			++it;
+		    }
+		}
+
+		ZeroHighlight(THEM);
+		Draw();
+		n = 0;
+	    }
+	}
+    }
+}
+
+void wizard_update_tolerance_cb(CALLBACK_ARGS)
+{
+    extern  void BDS_To_Mesh(Mesh *m);    
+    
+    if (THEM && THEM->bds && WID)
+    {
+	const double tol = WID->swiz_value[1]->value();
+	if (THEM->bds)delete THEM->bds;
+	THEM->bds = new BDS_Mesh;
+	printf("reading file %s\n",WID->surfmesh_filename.c_str());
+	THEM->bds->read_stl ( WID->surfmesh_filename.c_str(), tol );
+	BDS_To_Mesh (THEM); 
+	SetBoundingBox(); 
+	char a[25];
+	sprintf(a,"%d",THEM->bds->points.size());
+	WID->swiz_output[0]->value(a);
+	WID->swiz_output[1]->value("0");
+	WID->swiz_output[2]->value("0");
+	WID->swiz_output[3]->value("0");
+	Draw();
+    }
+}
+
+void wizard_update_next_cb(CALLBACK_ARGS)
+{
+    if (WID)
+    {
+	WID->swiz_wiz->next();
+    }
+}
+
+void wizard_update_prev_cb(CALLBACK_ARGS)
+{
+    if (WID)
+    {
+	WID->swiz_wiz->prev();
+    }
+}
+
 void options_ok_cb(CALLBACK_ARGS)
 {
   general_options_ok_cb(w, data);
diff --git a/Fltk/Callbacks.h b/Fltk/Callbacks.h
index 6d4b6fe31d..cb77b1b4b1 100644
--- a/Fltk/Callbacks.h
+++ b/Fltk/Callbacks.h
@@ -293,7 +293,10 @@ void solver_ok_cb(CALLBACK_ARGS);
 
 // SURFACE MESH WIZARD
 
+void wizard_update_more_edges_cb(CALLBACK_ARGS);
 void wizard_update_edges_cb(CALLBACK_ARGS);
-
+void wizard_update_tolerance_cb(CALLBACK_ARGS);
+void wizard_update_prev_cb(CALLBACK_ARGS);
+void wizard_update_next_cb(CALLBACK_ARGS);
 #endif
 
diff --git a/Fltk/GUI.cpp b/Fltk/GUI.cpp
index 9f2726458a..fb511d3eed 100644
--- a/Fltk/GUI.cpp
+++ b/Fltk/GUI.cpp
@@ -1,4 +1,4 @@
-// $Id: GUI.cpp,v 1.446 2005-07-03 08:02:23 geuzaine Exp $
+// $Id: GUI.cpp,v 1.447 2005-07-04 15:07:40 remacle Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -1612,62 +1612,137 @@ void GUI::reset_external_view_list()
   }
 }
 
-void GUI::create_surface_mesh_wizard()
+void GUI::create_surface_mesh_wizard(const char *name)
 {
+    if (name)
+	surfmesh_filename = name; 
+
     int width = 42 * fontsize;
-    int height = 12 * BH + 5 * WB;
+    int height = 5 * BH + 5 * WB;
 
-    if(swiz_window) {
+    if(swiz_window) {	
 	swiz_window->show();
 	return;
     }
-    swiz_window = new Dialog_Window(width, height);
+    swiz_window = new Dialog_Window(width, height, "Surface Mesh Wizard");
     swiz_window->box(GMSH_WINDOW_BOX);
+//    swiz_window->set_modal();
 
-    swiz_wiz = new Fl_Wizard(0, 0, width, height, "Surface Mesh Wizard");
+    swiz_wiz = new Fl_Wizard(0, 0, width, height, "Surface Mesh Wizard");    
     {
-	Fl_Tabs *o = new Fl_Tabs(WB, WB, width - 2 * WB, height - 2 * WB);
+//	Fl_Tabs *o = new Fl_Tabs(WB, WB, width - 2 * WB, height - 2 * WB);
+	{
+	    Fl_Group *o = new Fl_Group(WB, WB + BH, width - 2 * WB, height - 2 * WB - BH, "STL File");
+	    swiz_value[1] = new Fl_Value_Input(2 * WB, 2 * WB + 0 * BH, IW, BH, "Distance Tolerance");
+	    swiz_value[1]->value(5.e-7);
+	    swiz_value[1]->minimum(0);
+	    swiz_value[1]->maximum(1.e-4);
+	    swiz_value[1]->step(1.e-7);
+	    swiz_value[1]->align(FL_ALIGN_RIGHT);
+	    swiz_output[0] = new Fl_Output(2 * WB, 2 * WB + 1 * BH, IW, BH, "Number of Nodes after Merge");
+	    swiz_output[0]->align(FL_ALIGN_RIGHT);
+	    swiz_output[0]->value("0");
+	    {
+		Fl_Return_Button *b = new Fl_Return_Button(width - 4 * BB - 4 * WB, height - BH - WB, BB, BH, "Apply");
+		b->callback(wizard_update_tolerance_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Next");
+		b->callback(wizard_update_next_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Prev");
+		b->callback(wizard_update_prev_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - BB - WB, height - BH - WB, BB, BH, "Cancel");
+		b->callback(cancel_cb, (void *)swiz_window);
+	    }
+	    
+	    o->end();
+	}	
 	{
 	    Fl_Group *o = new Fl_Group(WB, WB + BH, width - 2 * WB, height - 2 * WB - BH, "Model Edge Detection");
-	    swiz_value[0] = new Fl_Value_Input(2 * WB, 2 * WB + 1 * BH, IW, BH, "Treshold Dihedral Angle");
+	    swiz_value[0] = new Fl_Value_Input(2 * WB, 2 * WB + 0 * BH, IW, BH, "Treshold Dihedral Angle");
 	    swiz_value[0]->value(180/8);
 	    swiz_value[0]->minimum(0);
 	    swiz_value[0]->maximum(90);
 	    swiz_value[0]->step(1);
 	    swiz_value[0]->align(FL_ALIGN_RIGHT);
+	    swiz_output[1] = new Fl_Output(2 * WB, 2 * WB + 1 * BH, IW, BH, "Number of Model Vertices");
+	    swiz_output[1]->align(FL_ALIGN_RIGHT);
+	    swiz_output[1]->value("0");
+	    swiz_output[2] = new Fl_Output(2 * WB, 2 * WB + 2 * BH, IW, BH, "Number of Model Edges");
+	    swiz_output[2]->align(FL_ALIGN_RIGHT);
+	    swiz_output[2]->value("0");
+	    swiz_output[3] = new Fl_Output(2 * WB, 2 * WB + 3 * BH, IW, BH, "Number of Model Faces");
+	    swiz_output[3]->align(FL_ALIGN_RIGHT);
+	    swiz_output[3]->value("0");
 	    {
-		Fl_Return_Button *b = new Fl_Return_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Apply");
+		Fl_Return_Button *b = new Fl_Return_Button(width - 4 * BB - 4 * WB, height - BH - WB, BB, BH, "Apply");
 		b->callback(wizard_update_edges_cb);
 	    }
 	    {
-		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Next");
-//		b->callback(options_save_cb);
+		Fl_Button *b = new Fl_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Next");
+		b->callback(wizard_update_next_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Prev");
+		b->callback(wizard_update_prev_cb);
 	    }
 	    {
 		Fl_Button *b = new Fl_Button(width - BB - WB, height - BH - WB, BB, BH, "Cancel");
 		b->callback(cancel_cb, (void *)swiz_window);
+	    }	    
+	    o->end();
+	}	
+	{ 
+	    Fl_Group *o = new Fl_Group(WB, WB + BH, width - 2 * WB, height - 2 * WB - BH, "Select More Edges");
+	    swiz_value[2] = new Fl_Value_Input(2 * WB, 2 * WB + 0 * BH, IW, BH, "Treshold Dihedral Angle For Extra Edges");
+	    swiz_value[2]->value(180/8);
+	    swiz_value[2]->minimum(0);
+	    swiz_value[2]->maximum(90);
+	    swiz_value[2]->step(1);
+	    swiz_value[2]->align(FL_ALIGN_RIGHT);
+	    {
+		Fl_Return_Button *b = new Fl_Return_Button(width - 4 * BB - 4 * WB, height - BH - WB, BB, BH, "Apply");
+		b->callback(wizard_update_more_edges_cb);
 	    }
-
+	    {
+		Fl_Button *b = new Fl_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Next");
+		b->callback(wizard_update_next_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Prev");
+		b->callback(wizard_update_prev_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - BB - WB, height - BH - WB, BB, BH, "Cancel");
+		b->callback(cancel_cb, (void *)swiz_window);
+	    }	    
 	    o->end();
 	}	
 	{
 	    Fl_Group *o = new Fl_Group(WB, WB + BH, width - 2 * WB, height - 2 * WB - BH, "Reverse Engineering the CAD");
 	    {
-		Fl_Return_Button *b = new Fl_Return_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Apply");
-//		b->callback(wizard_update_edges_cb);
+		Fl_Return_Button *b = new Fl_Return_Button(width - 4 * BB - 4 * WB, height - BH - WB, BB, BH, "Apply");
+		b->callback(wizard_update_edges_cb);
+	    }
+	    {
+		Fl_Button *b = new Fl_Button(width - 3 * BB - 3 * WB, height - BH - WB, BB, BH, "Next");
+		b->callback(wizard_update_next_cb);
 	    }
 	    {
-		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Next");
-//		b->callback(options_save_cb);
+		Fl_Button *b = new Fl_Button(width - 2 * BB - 2 * WB, height - BH - WB, BB, BH, "Prev");
+		b->callback(wizard_update_prev_cb);
 	    }
 	    {
 		Fl_Button *b = new Fl_Button(width - BB - WB, height - BH - WB, BB, BH, "Cancel");
 		b->callback(cancel_cb, (void *)swiz_window);
-	    }
-
+	    }	    
 	    o->end();
 	}	
-	o->end();
+//	o->end();
     } 
     swiz_wiz->end();
 }
diff --git a/Fltk/GUI.h b/Fltk/GUI.h
index 0aab262135..6f82ff9b62 100644
--- a/Fltk/GUI.h
+++ b/Fltk/GUI.h
@@ -53,6 +53,7 @@
 #endif
 
 #include <vector>
+#include <string>
 
 #include "Opengl_Window.h"
 #include "Colorbar_Window.h"
@@ -162,8 +163,9 @@ public:
   Fl_Button        *swiz_push_butt[20] ;
   Fl_Value_Input   *swiz_value[50] ;
   Fl_Input         *swiz_input[20] ;
-  Fl_Choice        *swiz_choice[20] ;
-
+  Fl_Output        *swiz_output[20] ;
+  Fl_Choice        *swiz_choice[20] ;  
+  std::string surfmesh_filename;
   // Option window
   Fl_Window        *opt_window ;
   Fl_Hold_Browser  *opt_browser ;
@@ -269,7 +271,7 @@ public:
   void create_menu_window();
   void create_graphic_window();
   void create_option_window();
-  void create_surface_mesh_wizard();
+  void create_surface_mesh_wizard(const char * = 0);
   void hide_all_option_groups();
   void create_general_options_window();
   void create_geometry_options_window();
diff --git a/Graphics/Geom.cpp b/Graphics/Geom.cpp
index 278b8df514..e116352874 100644
--- a/Graphics/Geom.cpp
+++ b/Graphics/Geom.cpp
@@ -1,4 +1,4 @@
-// $Id: Geom.cpp,v 1.88 2005-07-03 07:55:32 geuzaine Exp $
+// $Id: Geom.cpp,v 1.89 2005-07-04 15:07:40 remacle Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -156,10 +156,23 @@ void Draw_Curve(void *a, void *b)
       std::list<BDS_Edge*>::iterator ite = g->e.end();
       while (it!=ite){
 	BDS_Edge *e = (*it);
-	glBegin(GL_LINES);
-	glVertex3d(e->p1->X,e->p1->Y,e->p1->Z);
-	glVertex3d(e->p2->X,e->p2->Y,e->p2->Z);
-	glEnd();
+	if(CTX.geom.line_type < 1) {
+	    glBegin(GL_LINES);
+	    glVertex3d(e->p1->X,e->p1->Y,e->p1->Z);
+	    glVertex3d(e->p2->X,e->p2->Y,e->p2->Z);
+	    glEnd();
+	}
+	else
+	{
+	    x[0] = e->p1->X;
+	    y[0] = e->p1->Y;
+	    z[0] = e->p1->Z;
+	    x[1] = e->p2->X;
+	    y[1] = e->p2->Y;
+	    z[1] = e->p2->Z;
+	    Draw_Cylinder(c->ipar[3] ? CTX.geom.line_sel_width : CTX.geom.line_width,
+			  x, y, z, CTX.geom.light);
+	}
 	++it;
       }
     }
@@ -337,6 +350,7 @@ void Draw_Polygonal_Surface(Surface * s)
     glEnable(GL_POLYGON_OFFSET_FILL); // always!
 
     BDS_GeomEntity *g = s->bds->get_geom ( s->Num,2);	
+//    if (g->surf) glColor4ubv((GLubyte *) & CTX.color.geom.line);
     std::list<BDS_Triangle*>::iterator it  = g->t.begin();
     std::list<BDS_Triangle*>::iterator ite = g->t.end();
     while (it!=ite) {
diff --git a/Mesh/BDS.cpp b/Mesh/BDS.cpp
index 1c41d3bfab..7270cddc97 100644
--- a/Mesh/BDS.cpp
+++ b/Mesh/BDS.cpp
@@ -13,6 +13,60 @@
 
 */
 
+void BDS_Quadric::projection ( double xa, double ya, double za,
+			       double &x, double &y, double &z) const 
+{
+    x = xa;
+    y = ya;
+    z = za;
+    const double residual = signedDistanceTo (xa,ya,za);
+
+    int ITER = 0;
+    while (1)
+    {
+	// the equation to solve is 
+	// X = XA + z * N (XA)
+	// N = Gradient F (XA)
+	// F(X) is the quadric :: X^t A X + X^t B + C
+	// z^2 ( N^t A N ) + z ( XA^t A N + N^t A XA + N^t B) + F ( XA )  
+	BDS_Vector grad = Gradient (x,y,z);
+	double coef1 = 
+	    a * grad.x * grad.x +  
+	    b * grad.y * grad.y +  
+	    c * grad.z * grad.z +
+	    2 * d *  grad.x * grad.y +    
+	    2 * e *  grad.x * grad.z +  
+	    2 * f *  grad.y * grad.z ;
+	double coef2 = 
+	    2* a * xa * grad.x +  
+	    2* b * ya * grad.y +  
+	    2* c * za * grad.z +  
+	    2 * d *  (xa * grad.y + ya *grad.x)+   
+	    2 * e *  (xa * grad.z + za *grad.x)+   
+	    2 * f *  (za * grad.y + ya *grad.z)+
+	    grad.x * g +grad.y * h +grad.z * i;
+	double delta = coef2*coef2 - 4.*coef1*residual;
+	if (delta < 0){
+	    printf ("aaargh error projection");
+	}
+	else
+	{
+	    double alpha;
+	    double alpha1 = (-coef2 + sqrt (delta) ) / (2*coef1);
+	    double alpha2 = (-coef2 - sqrt (delta) ) / (2*coef1);
+	    if (fabs(alpha1) > fabs(alpha2))alpha = alpha2;
+	    else alpha = alpha1;
+	    x = xa + alpha * grad.x;
+	    y = ya + alpha * grad.y;
+	    z = za + alpha * grad.z;
+	}
+//	printf("point (%g,%g,%g) ==> (%g,%g,%g) dist %g\n",
+//	       xa,ya,za,x,y,z,signedDistanceTo (x,y,z));;
+	if (ITER++ == 1)	
+	    break;
+   }
+}
+
 BDS_Vector::BDS_Vector (const BDS_Point &p2,const BDS_Point &p1)
     : x(p2.X-p1.X),y(p2.Y-p1.Y),z(p2.Z-p1.Z)
 {    
@@ -393,7 +447,7 @@ void BDS_Mesh :: reverseEngineerCAD ( )
 		    for (int i=0;i<pts.size();i++)
 		    {
 			const double dist = PLANE  ( i , 0 ) * RSLT(0)+PLANE  ( i , 1 ) * RSLT(1)+PLANE  ( i , 2 ) * RSLT(2)+1;
-			if (fabs(dist) > 1.e-6 * LC)plane = false;
+			if (fabs(dist) > 1.e-5 * LC)plane = false;
 		    }
 		    
 		    if ( plane ) 
@@ -402,8 +456,9 @@ void BDS_Mesh :: reverseEngineerCAD ( )
 			(*it)->surf = new BDS_Plane (RSLT(0),RSLT(1),RSLT(2));
 		    }
 		}
-		if (!(*it)->surf && pts.size() > 20)
+		if (!(*it)->surf && pts.size() > 10)
 		{
+		    printf("coucou quadrique\n");
 		    Double_Matrix QUADRIC  ( pts.size() , 9 );
 		    Double_Vector ONES  ( pts.size() );
 		    Double_Vector RSLT  ( 9 );
@@ -438,8 +493,9 @@ void BDS_Mesh :: reverseEngineerCAD ( )
 			    QUADRIC  ( i , 5 ) * RSLT(5) +
 			    QUADRIC  ( i , 6 ) * RSLT(6) +
 			    QUADRIC  ( i , 7 ) * RSLT(7) +
-			    QUADRIC  ( i , 8 ) * RSLT(8) - 1;			    
-			if (fabs(dist) > 1.e-4 * LC )quad = false;
+			    QUADRIC  ( i , 8 ) * RSLT(8) - 1;
+//			printf("dist = %g LC %g\n",dist,LC);
+			if (fabs(dist) > 1.e-3 * LC )quad = false;
 		    }		    
 		    if ( quad ) 
 		    {
@@ -592,13 +648,22 @@ void BDS_Mesh :: classify ( double angle )
 	while (it!=ite)
 	{
 	    BDS_Edge &e = *((BDS_Edge *) *it);
-	    if ( e.numfaces() == 1) 
+	    if ( e.status == 1)
+	    {
+		e.g = &EDGE_CLAS;
+	    }
+
+	    else if ( e.status == -1)
+	    {
+	    }
+
+	    else if ( e.numfaces() == 1) 
 	    {
 		e.g = &EDGE_CLAS;		
 	    }
 	    
 	    else if (e.numfaces() == 2 && 
-		     e.faces(0)->g != e.faces(1)->g )
+		     e.faces(0)->status != e.faces(1)->status )
 	    {
 		e.g = &EDGE_CLAS;
 	    }
@@ -768,6 +833,7 @@ void BDS_Mesh :: classify ( double angle )
 	    else {
 		add_geom (vertextag, 0);
 		(*it)->g = get_geom(vertextag++,0);
+		(*it)->g->p = (*it);
 	    }	    
 	    if (!(*it)->g)printf("argh\n");
 	    ++it;
@@ -796,6 +862,7 @@ void BDS_Mesh :: classify ( double angle )
 		{
 		    add_geom (vertextag, 0);
 		  (*it)->g = get_geom(vertextag++,0);
+		  (*it)->g->p = (*it);
 		}
 	    }
 	    ++it;
@@ -1086,7 +1153,7 @@ bool BDS_Mesh :: read_stl ( const char *filename , const double tolerance)
 	delete [] DATA;
     }
     fclose (f);    
-
+    classify ( M_PI );
     return true;
 }
 
@@ -1152,12 +1219,7 @@ bool BDS_Mesh :: read_mesh ( const char *filename )
 		      fgets (buffer,255,f);
 		      sscanf (buffer,"%d %d %d %d",&n1,&n2,&n3,&cl);
 		      BDS_Triangle *t = add_triangle(n1,n2,n3);
-		      t->g = get_geom  (cl,2);
-		      if (!t->g)
-			{
-			  add_geom (cl,2);
-			  t->g = get_geom  (cl,2);
-			}
+		      t->status = cl;
 		    }
 		}
 	      else if (!strcmp (name,"Quadrilaterals"))
@@ -1175,13 +1237,8 @@ bool BDS_Mesh :: read_mesh ( const char *filename )
 		      // pas 12
 		      // pas 13
 		      BDS_Triangle *t2 = add_triangle(n1,n3,n4);
-		      t1->g = get_geom  (cl,2);
-		      if (!t1->g)
-			{
-			  add_geom (cl,2);
-			  t1->g = get_geom  (cl,2);
-			}
-		      t2->g = t1->g;
+		      t1->status = cl;
+		      t2->status = cl;
 		    }
 		}
 	    }
@@ -1447,6 +1504,11 @@ bool BDS_Mesh ::split_edge ( BDS_Edge *e, double coord)
 
     mid->g = ge;
 
+    if (mid->g->surf)
+    {
+	mid->g->surf->projection ( mid->X,mid->Y,mid->Z,mid->X,mid->Y,mid->Z);
+    }
+
     triangles.push_back(t1); 
     triangles.push_back (t2); 
     triangles.push_back (t3); 
@@ -1742,7 +1804,7 @@ bool BDS_Mesh ::smooth_point_b ( BDS_Point *p  )
 
     if (p->g->surf)
     {
-//	p->g->surf->projection ( X,Y,Z,p->X,p->Y,p->Z);
+	p->g->surf->projection ( p->X,p->Y,p->Z,p->X,p->Y,p->Z);
 //	printf("coucou\n");
     }
     else
@@ -1781,7 +1843,7 @@ bool BDS_Mesh ::smooth_point ( BDS_Point *p , BDS_Mesh *geom_mesh )
     }
     else
     {
-//	return false;
+	return false;
 	p->X = X;
 	p->Y = Y;
 	p->Z = Z;
diff --git a/Mesh/BDS.h b/Mesh/BDS.h
index 7eac5c726d..2186b7ab89 100644
--- a/Mesh/BDS.h
+++ b/Mesh/BDS.h
@@ -34,7 +34,7 @@ public:
 
     std::list<BDS_Triangle *> t;
     std::list<BDS_Edge *>     e;
-
+    BDS_Point   *p;
     BDS_Surface *surf;
 
     inline bool operator <  ( const BDS_GeomEntity & other ) const
@@ -45,7 +45,7 @@ public:
 	    return false;
 	}
     BDS_GeomEntity (int a, int b)
-	: classif_tag (a),classif_degree(b),surf(0)
+	: classif_tag (a),classif_degree(b),p(0),surf(0)
 	{
 	}
 };
@@ -185,6 +185,7 @@ class BDS_Edge
     std::vector <BDS_Triangle *> _faces;
 public:
     bool deleted;
+    int status;
     BDS_Point *p1,*p2;
     BDS_GeomEntity *g;
     inline BDS_Triangle* faces(int i) const
@@ -241,7 +242,7 @@ public:
     inline void oppositeof (BDS_Point * oface[2]) const; 
 
     BDS_Edge ( BDS_Point *A, BDS_Point *B )
-	: deleted(false), g(0)
+	: deleted(false), status(0),g(0)
 	{	    
 	    if (*A < *B) 
 	    {
@@ -262,6 +263,7 @@ class BDS_Triangle
 {
 public:
     bool deleted;
+    int status;
     BDS_Edge *e1,*e2,*e3;
     BDS_Vector N() const ;
     BDS_GeomEntity *g;
@@ -282,7 +284,7 @@ public:
 	  n[2] = e2->commonvertex (e3);
 	}
     BDS_Triangle ( BDS_Edge *A, BDS_Edge *B, BDS_Edge *C)
-	: deleted (false) , e1(A),e2(B),e3(C),g(0)
+	: deleted (false) , status(0), e1(A),e2(B),e3(C),g(0)
 	{	
 	    e1->addface(this);
 	    e2->addface(this);
@@ -320,6 +322,14 @@ public :
 	: a(A),b(B),c(C),d(D),e(E),f(F),g(G),h(H),i(I)
 	{
 	}
+
+    virtual BDS_Vector Gradient ( double x, double y, double z ) const 
+	{
+	    return BDS_Vector ( 2* ( a * x + d * y + e * z ) + g ,
+				2* ( d * x + b * y + f * z ) + h ,
+				2* ( e * x + f * y + c * z ) + i );
+	}
+
     virtual double signedDistanceTo (  double x, double y, double z ) const {
 	const double q = 
 	    a * x * x +  
@@ -334,14 +344,7 @@ public :
 	return q;
     }
     virtual void projection ( double xa, double ya, double za,
-			      double &x, double &y, double &z) const 
-	{
-	    // not done yet
-	    // this is not as simple as for the plane
-	    // you shoud have min (signedDistance), this can
-	    // be done using the GSL... 2 BE DONE !!!!!
-	    throw;
-	}
+			      double &x, double &y, double &z) const ;
     virtual std::string nameOf () const {return std::string("Quadric");}
 };
 
diff --git a/Mesh/DiscreteSurface.cpp b/Mesh/DiscreteSurface.cpp
index 82d320df74..2ed3a62235 100644
--- a/Mesh/DiscreteSurface.cpp
+++ b/Mesh/DiscreteSurface.cpp
@@ -1,4 +1,4 @@
-// $Id: DiscreteSurface.cpp,v 1.16 2005-06-27 15:03:46 remacle Exp $
+// $Id: DiscreteSurface.cpp,v 1.17 2005-07-04 15:07:41 remacle Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -61,12 +61,69 @@ void Mesh_To_BDS(Surface *s, BDS_Mesh *m)
 
 }
 
+void BDS_To_Mesh_2(Mesh *m)
+{
+    Tree_Action(m->Vertices, Free_Vertex);  
+    Tree_Delete(m->Vertices);
+    m->Vertices = Tree_Create(sizeof(Vertex *), compareVertex);
+
+    {
+	std::set<BDS_Point*, PointLessThan>::iterator it   = m->bds_mesh->points.begin();
+	std::set<BDS_Point*, PointLessThan>::iterator ite  = m->bds_mesh->points.end();
+	while (it != ite)
+	{
+	    Vertex *vert = Create_Vertex((*it)->iD, (*it)->X,(*it)->Y,(*it)->Z, 1.0, 0.0);
+	    Tree_Add (m->Vertices,&vert);
+	    ++it;
+	}
+    }
+    {
+	std::set<BDS_Edge*, EdgeLessThan>::iterator it  = m->bds_mesh->edges.begin();
+	std::set<BDS_Edge*, EdgeLessThan>::iterator ite = m->bds_mesh->edges.end();
+	while(it!=ite)
+	{
+	    BDS_GeomEntity *g = (*it)->g;
+	    if (g && g->classif_degree == 1)
+	    {
+		Vertex *v1 = FindVertex((*it)->p1->iD, m);
+		Vertex *v2 = FindVertex((*it)->p2->iD, m);
+		SimplexBase *simp = Create_SimplexBase(v1,v2,NULL, NULL);
+		Curve *c = FindCurve (g->classif_tag,m);
+		simp->iEnt = g->classif_tag;
+		Tree_Insert(c->SimplexesBase, &simp);
+	    }
+	    ++it;
+	}
+    }
+    {
+	std::list<BDS_Triangle*>::iterator it  = m->bds_mesh->triangles.begin();
+	std::list<BDS_Triangle*>::iterator ite = m->bds_mesh->triangles.end();
+	while (it!=ite){
+	    BDS_Point *nod[3];
+	    (*it)->getNodes (nod);
+	    Vertex *v1 = FindVertex(nod[0]->iD, m);
+	    Vertex *v2 = FindVertex(nod[1]->iD, m);
+	    Vertex *v3 = FindVertex(nod[2]->iD, m);
+	    SimplexBase *simp = Create_SimplexBase(v1,v2,v3, NULL);
+	    BDS_GeomEntity *g = (*it)->g;
+	    Surface *s = FindSurface (g->classif_tag,m);
+	    simp->iEnt = g->classif_tag;
+	    Tree_Add(s->Simplexes, &simp);
+	    ++it;
+	}
+    }
+
+
+}
 void BDS_To_Mesh(Mesh *m)
 {
+    Tree_Action(m->Points, Free_Vertex);  
+    Tree_Delete(m->Points);
     Tree_Action(m->Surfaces, Free_Surface);
     Tree_Action(m->Curves, Free_Curve);
     Tree_Delete(m->Surfaces);
     Tree_Delete(m->Curves);
+    m->Points = Tree_Create(sizeof(Vertex *), compareVertex);
     m->Curves = Tree_Create(sizeof(Curve *), compareCurve);
     m->Surfaces = Tree_Create(sizeof(Surface *), compareSurface);
 
@@ -84,6 +141,9 @@ void BDS_To_Mesh(Mesh *m)
 	    _Surf->bds = m->bds;
 	    End_Surface(_Surf);
 	    Tree_Add(m->Surfaces, &_Surf);
+	    
+	   
+
 	}
 	else if ((*it)->classif_degree ==1  )
 	{
@@ -95,6 +155,12 @@ void BDS_To_Mesh(Mesh *m)
 	    End_Curve(_c);
 	    Tree_Add(m->Curves, &_c);
 	}
+	else if ((*it)->classif_degree == 0 )
+	{
+	    BDS_Point *p = (*it)->p;
+	    Vertex *_v = Create_Vertex(p->iD, p->X,p->Y,p->Z,1,0);	
+	    Tree_Add(m->Points, &_v);
+	}
 	++it;
     }
 
@@ -118,6 +184,8 @@ int MeshDiscreteSurface(Surface *s)
 	printf("iter %d done\n",iter);
 	iter ++;
       }
+      BDS_To_Mesh_2(THEM);
+      printf("mesh has %d vertices (%d)\n",Tree_Nbr(THEM->Vertices),THEM->bds->points.size());
       THEM->bds_mesh->save_gmsh_format ( "3.msh" );
     }
     return 1;
diff --git a/Parser/OpenFile.cpp b/Parser/OpenFile.cpp
index 9eced090e5..e1915ab7e6 100644
--- a/Parser/OpenFile.cpp
+++ b/Parser/OpenFile.cpp
@@ -1,4 +1,4 @@
-// $Id: OpenFile.cpp,v 1.80 2005-06-27 15:03:46 remacle Exp $
+// $Id: OpenFile.cpp,v 1.81 2005-07-04 15:07:41 remacle Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -285,17 +285,20 @@ int MergeProblem(char *name, int warn_if_missing)
     if(!strcmp(ext, ".mesh"))
     {
       THEM->bds->read_mesh ( name );
+#if defined(HAVE_FLTK)
+      WID->create_surface_mesh_wizard(name);
+#endif
     }
-    else       
-      {
-	// STEREO LITOGRAPHY (STL) FILE
-	THEM->bds->read_stl ( name , 5.e-7 );
-      }
-    THEM->bds->save_gmsh_format ( "1.msh" );
-    THEM->bds->classify ( M_PI / 8 );
+    else
+    {
 #if defined(HAVE_FLTK)
-    WID->create_surface_mesh_wizard();
+	WID->create_surface_mesh_wizard(name);
 #endif
+    }
+
+
+    THEM->bds->save_gmsh_format ( "1.msh" );
+    THEM->bds->classify ( M_PI / 8 );
     THEM->bds->save_gmsh_format ( "2.msh" );
     BDS_To_Mesh (THEM);
     SetBoundingBox();
diff --git a/benchmarks/stl/ktoolcav.stl b/benchmarks/stl/ktoolcav.stl
index 357cddaef8d1106a893478e95e201d116775d3f5..65ffc1f12f178c9533f6e8b4152be2665a602cf3 100644
GIT binary patch
delta 184
zcmbPzi0AiVo`x327N#xC(u(yisX3`7sS4Fx3JN}%d8yS#mKLV9)?Bq*SR@Q^NSNc&
hVUA0O8BB*5y6t9QEz|Agne`Y=r{~Kvi?bSngaETtHWL5<

delta 15
Wcmex;m}mAOo`x327N#xC(ux2%6$P>Y

-- 
GitLab