diff --git a/Mesh/2D_Mesh.cpp b/Mesh/2D_Mesh.cpp
index b47cc4483d5fd119fc1d837d4b0a5c016d87ccbd..0000cf8d88d27534a8051f13c65ef53853ec60d2 100644
--- a/Mesh/2D_Mesh.cpp
+++ b/Mesh/2D_Mesh.cpp
@@ -1,4 +1,4 @@
-// $Id: 2D_Mesh.cpp,v 1.65 2004-06-23 18:52:45 geuzaine Exp $
+// $Id: 2D_Mesh.cpp,v 1.66 2004-06-23 22:25:01 geuzaine Exp $
//
// Copyright (C) 1997-2004 C. Geuzaine, J.-F. Remacle
//
@@ -715,12 +715,10 @@ void filldel(Delaunay * deladd, int aa, int bb, int cc,
deladd->v.voisin2 = NULL;
deladd->v.voisin3 = NULL;
- CircumCircle(points[aa].where.h,
- points[aa].where.v,
- points[bb].where.h,
- points[bb].where.v,
- points[cc].where.h,
- points[cc].where.v, &deladd->t.xc, &deladd->t.yc);
+ CircumCircle(points[aa].where.h, points[aa].where.v,
+ points[bb].where.h, points[bb].where.v,
+ points[cc].where.h, points[cc].where.v,
+ &deladd->t.xc, &deladd->t.yc);
pt2.h = deladd->t.xc;
pt2.v = deladd->t.yc;
@@ -789,85 +787,24 @@ void ActionInvertQua(void *a, void *b)
q->V[3] = tmp;
}
-int isPointOnPlanarSurface(Surface * S, double X, double Y, double Z, double n[3]);
-
-int isMiddlePointOnPlanarSurface(Surface *s, Vertex *v1, Vertex *v2)
-{
- double n[3] = {s->a, s->b, s->c};
- norme(n);
- double x = 0.5 * (v1->Pos.X + v2->Pos.X);
- double y = 0.5 * (v1->Pos.Y + v2->Pos.Y);
- double z = 0.5 * (v1->Pos.Z + v2->Pos.Z);
- return isPointOnPlanarSurface(s, x, y, z, n);
-}
-
-void Get_SurfaceNormal(Surface *s, double n[3])
-{
- double t1[3], t2[3];
-
- if(s->Typ == MSH_SURF_PLAN){
- // don't use s->plan and co.: these are computed in MeanPlane,
- // which borks the orientation... we cannot use
- // InterpolateSurface euther, since we use Calcule_Z_plan in
- // there, which relies in the MeanPlane stuff too.
- List_T *points = List_Create(10, 10, sizeof(Vertex *));
- for(int i = 0; i < List_Nbr(s->Generatrices); i++) {
- Curve *c;
- List_Read(s->Generatrices, i, &c);
- // no need to play with beg/end: negative curves are already inverted
- List_Add(points, &c->beg);
- }
- if(List_Nbr(points) > 2){
- Vertex *v1 = *(Vertex**)List_Pointer(points, 0);
- Vertex *v2 = *(Vertex**)List_Pointer(points, 1);
- t1[0] = v2->Pos.X - v1->Pos.X;
- t1[1] = v2->Pos.Y - v1->Pos.Y;
- t1[2] = v2->Pos.Z - v1->Pos.Z;
- for(int i = 2; i < List_Nbr(points); i++){
- Vertex *v3 = *(Vertex**)List_Pointer(points, i);
- t2[0] = v3->Pos.X - v1->Pos.X;
- t2[1] = v3->Pos.Y - v1->Pos.Y;
- t2[2] = v3->Pos.Z - v1->Pos.Z;
- prodve(t1, t2, n);
- if(norme(n) && isMiddlePointOnPlanarSurface(s, v1, v3))
- break;
- }
- }
- if(List_Nbr(points) <= 2 || !norme(n)){
- Msg(WARNING, "Couldn't compute normal to Surface %d using control points: "
- "reverting to mean plane", s->Num);
- n[0] = s->a;
- n[1] = s->b;
- n[2] = s->c;
- norme(n);
- }
- List_Delete(points);
- }
- else{
- Vertex v1 = InterpolateSurface(s, 0.5, 0.5, 0, 0);
- Vertex v2 = InterpolateSurface(s, 0.6, 0.5, 0, 0);
- Vertex v3 = InterpolateSurface(s, 0.5, 0.6, 0, 0);
- t1[0] = v2.Pos.X - v1.Pos.X;
- t1[1] = v2.Pos.Y - v1.Pos.Y;
- t1[2] = v2.Pos.Z - v1.Pos.Z;
- t2[0] = v3.Pos.X - v1.Pos.X;
- t2[1] = v3.Pos.Y - v1.Pos.Y;
- t2[2] = v3.Pos.Z - v1.Pos.Z;
- prodve(t1, t2, n);
- norme(n);
- }
-}
-
-// A little routine to re-orient the surface meshe to match the
-// geometrical orientation. It would be better to actually generate
-// the meshes correctly in the first place (!), but this is handy
-// since we use many different algorithms since the mean plane
-// computation doesn't take the original orientation into account.
+// A little routine that re-orients the surface meshes so that they
+// match the geometrical orientations. It would be better to generate
+// the meshes correctly in the first place (!), but it's actually
+// pretty nice to do it here since we use many algorithms (some
+// closed) and since the mean plane computation doesn't take the
+// original orientation into account.
void ReOrientSurfaceMesh(Surface *s)
{
+ if(!List_Nbr(s->Generatrices))
+ return;
+
Curve *c;
List_Read(s->Generatrices, 0, &c);
+
+ if(List_Nbr(c->Vertices) < 2)
+ return;
+
Vertex *v1, *v2;
List_Read(c->Vertices, 0, &v1);
List_Read(c->Vertices, 1, &v2);
@@ -879,15 +816,22 @@ void ReOrientSurfaceMesh(Surface *s)
EdgesContainer edges;
edges.AddSimplexTree(s->Simplexes);
Edge *e = (Edge*)Tree_PQuery(edges.AllEdges, &edge);
- if(e && e->V[0]->Num == v2->Num && e->V[1]->Num == v1->Num)
+ if(e && e->V[0]->Num == v2->Num && e->V[1]->Num == v1->Num){
Tree_Action(s->Simplexes, ActionInvertTri);
+ if(Tree_Nbr(s->Quadrangles)) // for partially recombined meshes
+ Tree_Action(s->Quadrangles, ActionInvertQua);
+ }
}
+
if(Tree_Nbr(s->Quadrangles)){
EdgesContainer edges;
edges.AddQuadrangleTree(s->Quadrangles);
Edge *e = (Edge*)Tree_PQuery(edges.AllEdges, &edge);
- if(e && e->V[0]->Num == v2->Num && e->V[1]->Num == v1->Num)
+ if(e && e->V[0]->Num == v2->Num && e->V[1]->Num == v1->Num){
Tree_Action(s->Quadrangles, ActionInvertQua);
+ if(Tree_Nbr(s->Simplexes)) // for partially recombined meshes
+ Tree_Action(s->Simplexes, ActionInvertTri);
+ }
}
}
@@ -996,6 +940,4 @@ void Maillage_Surface(void *data, void *dum)
ReOrientSurfaceMesh(s);
}
-
}
-