diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/algprim.cpp netgen-4.5_new/libsrc/csg/algprim.cpp --- netgen-4.5_orig/libsrc/csg/algprim.cpp 2006-01-25 16:30:28.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/algprim.cpp 2010-06-23 12:35:22.000000000 +0400 @@ -108,7 +108,7 @@ void Plane :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "plane"; + classname = (char*)"plane"; coeffs.SetSize (6); coeffs.Elem(1) = p(0); coeffs.Elem(2) = p(1); @@ -355,7 +355,7 @@ void Sphere :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "sphere"; + classname = (char*)"sphere"; coeffs.SetSize (4); coeffs.Elem(1) = c(0); coeffs.Elem(2) = c(1); @@ -760,7 +760,7 @@ void Cylinder :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "cylinder"; + classname = (char*)"cylinder"; coeffs.SetSize (7); coeffs.Elem(1) = a(0); coeffs.Elem(2) = a(1); @@ -1243,7 +1243,7 @@ void Cone :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "cone"; + classname = (char*)"cone"; coeffs.SetSize (8); coeffs.Elem(1) = a(0); coeffs.Elem(2) = a(1); @@ -1446,7 +1446,7 @@ void Torus :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "torus"; + classname = (char*)"torus"; coeffs.SetSize (8); coeffs.Elem(1) = c(0); coeffs.Elem(2) = c(1); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/brick.cpp netgen-4.5_new/libsrc/csg/brick.cpp --- netgen-4.5_orig/libsrc/csg/brick.cpp 2006-02-08 15:23:15.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/brick.cpp 2010-06-23 12:35:35.000000000 +0400 @@ -345,7 +345,7 @@ void Brick :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "brick"; + classname = (char*)"brick"; coeffs.SetSize(12); coeffs.Elem(1) = p1(0); coeffs.Elem(2) = p1(1); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/meshsurf.cpp netgen-4.5_new/libsrc/csg/meshsurf.cpp --- netgen-4.5_orig/libsrc/csg/meshsurf.cpp 2006-02-14 11:54:35.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/meshsurf.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -77,11 +77,12 @@ } -void MeshOptimize2dSurfaces :: ProjectPoint (INDEX surfind, Point3d & p) const +bool MeshOptimize2dSurfaces :: ProjectPoint (INDEX surfind, Point3d & p) const { Point<3> hp = p; geometry.GetSurface(surfind)->Project (hp); p = hp; + return true; } void MeshOptimize2dSurfaces :: ProjectPoint2 (INDEX surfind, INDEX surfind2, diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/meshsurf.hpp netgen-4.5_new/libsrc/csg/meshsurf.hpp --- netgen-4.5_orig/libsrc/csg/meshsurf.hpp 2004-01-20 14:49:44.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/meshsurf.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -45,7 +45,7 @@ MeshOptimize2dSurfaces (const CSGeometry & ageometry); /// - virtual void ProjectPoint (INDEX surfind, Point3d & p) const; + virtual bool ProjectPoint (INDEX surfind, Point3d & p) const; /// virtual void ProjectPoint2 (INDEX surfind, INDEX surfind2, Point3d & p) const; /// diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/polyhedra.cpp netgen-4.5_new/libsrc/csg/polyhedra.cpp --- netgen-4.5_orig/libsrc/csg/polyhedra.cpp 2006-02-09 13:33:11.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/polyhedra.cpp 2010-06-23 12:36:00.000000000 +0400 @@ -287,7 +287,7 @@ void Polyhedra :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "Polyhedra"; + classname = (char*)"Polyhedra"; coeffs.SetSize(0); coeffs.Append (points.Size()); coeffs.Append (faces.Size()); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/csg/surface.cpp netgen-4.5_new/libsrc/csg/surface.cpp --- netgen-4.5_orig/libsrc/csg/surface.cpp 2006-02-08 15:23:16.000000000 +0300 +++ netgen-4.5_new/libsrc/csg/surface.cpp 2010-06-23 12:35:47.000000000 +0400 @@ -215,7 +215,7 @@ void Primitive :: GetPrimitiveData (char *& classname, ARRAY & coeffs) const { - classname = "undef"; + classname = (char*)"undef"; coeffs.SetSize (0); } diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/general/profiler.cpp netgen-4.5_new/libsrc/general/profiler.cpp --- netgen-4.5_orig/libsrc/general/profiler.cpp 2006-01-11 13:05:59.000000000 +0300 +++ netgen-4.5_new/libsrc/general/profiler.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -34,8 +34,14 @@ { StopTimer (total_timer); - ofstream prof ("netgen.prof"); - Print (prof); + char* env; + if ((env = getenv("NETGEN_PROF")) && !strcmp(env, "1")) { + ofstream prof ("netgen.prof"); + Print (prof); + } + else if ((env = getenv("NETGEN_PROF")) && !strcmp(env, "0")) { + Print (std::cout); + } } diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/geom2d/genmesh2d.cpp netgen-4.5_new/libsrc/geom2d/genmesh2d.cpp --- netgen-4.5_orig/libsrc/geom2d/genmesh2d.cpp 2006-02-16 19:17:47.000000000 +0300 +++ netgen-4.5_new/libsrc/geom2d/genmesh2d.cpp 2010-06-23 12:36:59.000000000 +0400 @@ -121,11 +121,11 @@ int hsteps = mp.optsteps2d; - mp.optimize2d = "smcm"; + mp.optimize2d = (char*)"smcm"; mp.optsteps2d = hsteps/2; Optimize2d (*mesh, mp); - mp.optimize2d = "Smcm"; + mp.optimize2d = (char*)"Smcm"; mp.optsteps2d = (hsteps+1)/2; Optimize2d (*mesh, mp); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/gprim/geom2d.hpp netgen-4.5_new/libsrc/gprim/geom2d.hpp --- netgen-4.5_orig/libsrc/gprim/geom2d.hpp 2004-01-20 14:49:44.000000000 +0300 +++ netgen-4.5_new/libsrc/gprim/geom2d.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -53,7 +53,7 @@ int IsOnLongLine (const Line2d & l, const Point2d & p); int Hit (const Line2d & l1, const Line2d & l2, double heps = EPSGEOM); ostream & operator<<(ostream & s, const Line2d & l); -Point2d CrossPoint (const PLine2d & l1, const PLine2d & l2); +Point2d CrossPoint (const Line2d & l1, const Line2d & l2); int Parallel (const PLine2d & l1, const PLine2d & l2, double peps = EPSGEOM); int IsOnLine (const PLine2d & l, const Point2d & p, double heps = EPSGEOM); int IsOnLongLine (const PLine2d & l, const Point2d & p); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/gprim/geom3d.hpp netgen-4.5_new/libsrc/gprim/geom3d.hpp --- netgen-4.5_orig/libsrc/gprim/geom3d.hpp 2004-08-30 16:04:04.000000000 +0400 +++ netgen-4.5_new/libsrc/gprim/geom3d.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -25,6 +25,7 @@ inline Point3d Center (const Point3d & p1, const Point3d & p2, const Point3d & p3); inline Point3d Center (const Point3d & p1, const Point3d & p2, const Point3d & p3, const Point3d & p4); +inline double Dist2 (const Point3d & p1, const Point3d & p2); ostream & operator<<(ostream & s, const Point3d & p); inline Vec3d operator- (const Vec3d & p1, const Vec3d & v); inline Vec3d operator+ (const Vec3d & p1, const Vec3d & v); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/interface/Makefile netgen-4.5_new/libsrc/interface/Makefile --- netgen-4.5_orig/libsrc/interface/Makefile 2005-08-09 18:14:59.000000000 +0400 +++ netgen-4.5_new/libsrc/interface/Makefile 2010-06-23 13:19:48.000000000 +0400 @@ -1,4 +1,5 @@ -src = nginterface.cpp writeuser.cpp writediffpack.cpp writeabaqus.cpp writefluent.cpp writepermas.cpp writetochnog.cpp writetecplot.cpp wuchemnitz.cpp writetochnog.cpp writefeap.cpp writeelmer.cpp writegmsh.cpp writejcm.cpp readuser.cpp importsolution.cpp +#src = nginterface.cpp writeuser.cpp writediffpack.cpp writeabaqus.cpp writefluent.cpp writepermas.cpp writetochnog.cpp writetecplot.cpp wuchemnitz.cpp writetochnog.cpp writefeap.cpp writeelmer.cpp writegmsh.cpp writejcm.cpp readuser.cpp importsolution.cpp +src = writeuser.cpp writediffpack.cpp writeabaqus.cpp writefluent.cpp writepermas.cpp writetochnog.cpp writetecplot.cpp wuchemnitz.cpp writetochnog.cpp writefeap.cpp writeelmer.cpp writegmsh.cpp writejcm.cpp readuser.cpp nglib.cpp ngnewdelete.cpp # lib = nginterface libpath = libsrc/interface diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/interface/nglib.cpp netgen-4.5_new/libsrc/interface/nglib.cpp --- netgen-4.5_orig/libsrc/interface/nglib.cpp 2005-10-18 17:53:18.000000000 +0400 +++ netgen-4.5_new/libsrc/interface/nglib.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -56,7 +56,8 @@ void Ng_Exit () { - ; + delete testout; + testout = NULL; } diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/interface/writeuser.cpp netgen-4.5_new/libsrc/interface/writeuser.cpp --- netgen-4.5_orig/libsrc/interface/writeuser.cpp 2005-08-09 18:14:59.000000000 +0400 +++ netgen-4.5_new/libsrc/interface/writeuser.cpp 2010-06-23 12:37:42.000000000 +0400 @@ -17,7 +17,7 @@ void RegisterUserFormats (ARRAY & names) { - char *types[] = + const char *types[] = { "Neutral Format", "Surface Mesh Format" , diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/linalg/densemat.hpp netgen-4.5_new/libsrc/linalg/densemat.hpp --- netgen-4.5_orig/libsrc/linalg/densemat.hpp 2005-12-09 15:26:19.000000000 +0300 +++ netgen-4.5_new/libsrc/linalg/densemat.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -14,6 +14,8 @@ #include +class DenseMatrix; +void CalcAtA (const DenseMatrix & a, DenseMatrix & m2); class DenseMatrix { diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/makefile.inc netgen-4.5_new/libsrc/makefile.inc --- netgen-4.5_orig/libsrc/makefile.inc 2005-09-02 17:17:51.000000000 +0400 +++ netgen-4.5_new/libsrc/makefile.inc 2010-06-23 13:19:48.000000000 +0400 @@ -8,17 +8,14 @@ LIBSRC_DIR=$(CPP_DIR)/libsrc LIB_DIR=$(CPP_DIR)/lib/$(MACHINE) -#OCC_DIR=../../occ -#OCCINC_DIR=$(OCC_DIR)/inc -#OCCLIB_DIR=$(OCC_DIR)/lib -# OCC_DIR=/opt/OpenCASCADE5.2/ros -# OCC_DIR=/home/joachim/download/occ/Linux -# OCCINC_DIR=$(OCC_DIR)/inc -I$(OCC_DIR)/ros/inc -# OCCLIB_DIR=$(OCC_DIR)/Linux/lib +OCC_DIR=$(CASROOT) +OCCINC_DIR=$(OCC_DIR)/inc +OCCLIB_DIR=$(OCC_DIR)/Linux/lib # include $(LIBSRC_DIR)/makefile.mach.$(MACHINE) # -CPLUSPLUSFLAGS1 = -c -I$(LIBSRC_DIR)/include -I$(OCCINC_DIR) +CPLUSPLUSFLAGS1 = -c -fPIC -I$(LIBSRC_DIR)/include -I$(OCCINC_DIR) \ + -DOCCGEOMETRY -DOCC52 -DHAVE_IOSTREAM -DHAVE_LIMITS_H # ARFLAGS = r # diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/makefile.mach.LINUX netgen-4.5_new/libsrc/makefile.mach.LINUX --- netgen-4.5_orig/libsrc/makefile.mach.LINUX 2004-10-11 23:49:26.000000000 +0400 +++ netgen-4.5_new/libsrc/makefile.mach.LINUX 2010-06-23 13:19:48.000000000 +0400 @@ -16,7 +16,7 @@ # CFLAGS2 = -CPLUSPLUSFLAGS2 = -O2 -I/usr/include/GL3.5 -DLINUX -DOPENGL \ +CPLUSPLUSFLAGS2 = -O2 -I/usr/include/GL3.5 -DLINUX \ -ftemplate-depth-99 -finline-limit=10000 \ -Wdisabled-optimization -funroll-loops -DnoNGSOLVE diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/meshing/improve2.cpp netgen-4.5_new/libsrc/meshing/improve2.cpp --- netgen-4.5_orig/libsrc/meshing/improve2.cpp 2006-01-11 19:08:19.000000000 +0300 +++ netgen-4.5_new/libsrc/meshing/improve2.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -4,7 +4,7 @@ #include #ifndef SMALLLIB -#include +//#include #endif namespace netgen diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/meshing/improve2.hpp netgen-4.5_new/libsrc/meshing/improve2.hpp --- netgen-4.5_orig/libsrc/meshing/improve2.hpp 2004-10-12 23:22:55.000000000 +0400 +++ netgen-4.5_new/libsrc/meshing/improve2.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -32,17 +32,16 @@ /// virtual void SelectSurfaceOfPoint (const Point3d & p, const PointGeomInfo & gi); - /// - virtual void ProjectPoint (INDEX /* surfind */, Point3d & /* p */) const { }; + + /// project point on surface, returns true if success + virtual bool ProjectPoint (INDEX /* surfind */, Point3d & /* p */) const { return false; } + /// fast project point on surface using point geom info of a neighboring point + /// if gi.trignum != 0, + /// returns true if success, gi is updated + virtual bool ProjectPoint (INDEX surfind, Point3d & p, PointGeomInfo& gi) const + { gi.trignum = 1; return ProjectPoint (surfind, p); } /// virtual void ProjectPoint2 (INDEX /* surfind */, INDEX /* surfind2 */, Point3d & /* p */) const { }; - /// liefert zu einem 3d-Punkt die geominfo (Dreieck) und liefert 1, wenn erfolgreich, - /// 0, wenn nicht (Punkt ausserhalb von chart) - virtual int CalcPointGeomInfo(PointGeomInfo& gi, const Point3d& /*p3*/) const - { gi.trignum = 1; return 1;}; - - virtual int CalcPointGeomInfo(int /* surfind */, PointGeomInfo& gi, const Point3d& p3) const - { return CalcPointGeomInfo (gi, p3); } /// virtual void GetNormalVector(INDEX surfind, const Point3d & p, PointGeomInfo & gi, Vec3d & n) const; diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/meshing/meshtype.cpp netgen-4.5_new/libsrc/meshing/meshtype.cpp --- netgen-4.5_orig/libsrc/meshing/meshtype.cpp 2006-02-10 13:11:08.000000000 +0300 +++ netgen-4.5_new/libsrc/meshing/meshtype.cpp 2010-06-23 12:39:02.000000000 +0400 @@ -1,4 +1,5 @@ #include +#include #include "meshing.hpp" @@ -774,7 +775,7 @@ frob /= 2; double det = trans.Det(); - if (det <= 0) + if (det <= DBL_MIN) err += 1e12; else err += frob * frob / det; @@ -2222,9 +2223,9 @@ MeshingParameters :: MeshingParameters () { - optimize3d = "cmdmstm"; + optimize3d = (char*)"cmdmstm"; optsteps3d = 3; - optimize2d = "smsmsmSmSmSm"; + optimize2d = (char*)"smsmsmSmSmSm"; optsteps2d = 3; opterrpow = 2; blockfill = 1; diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/meshing/meshtype.hpp netgen-4.5_new/libsrc/meshing/meshtype.hpp --- netgen-4.5_orig/libsrc/meshing/meshtype.hpp 2006-02-10 13:11:08.000000000 +0300 +++ netgen-4.5_new/libsrc/meshing/meshtype.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -13,7 +13,7 @@ Classes for NETGEN */ - +class Mesh; enum ELEMENT_TYPE { SEGMENT = 1, SEGMENT3 = 2, TRIG = 10, QUAD=11, TRIG6 = 12, QUAD6 = 13, QUAD8 = 14, diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/meshing/smoothing2.cpp netgen-4.5_new/libsrc/meshing/smoothing2.cpp --- netgen-4.5_orig/libsrc/meshing/smoothing2.cpp 2006-01-11 19:08:20.000000000 +0300 +++ netgen-4.5_new/libsrc/meshing/smoothing2.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -300,7 +300,7 @@ double Opti2SurfaceMinFunction :: FuncGrad (const Vector & x, Vector & grad) const { - Vec3d n, vgrad; + Vec3d vgrad; Point3d pp1; double g1x, g1y; double badness, hbadness; @@ -308,8 +308,6 @@ vgrad = 0; badness = 0; - meshthis -> GetNormalVector (surfi, sp1, gi1, n); - pp1 = sp1; pp1.Add2 (x.Get(1), t1, x.Get(2), t2); @@ -360,7 +358,7 @@ double Opti2SurfaceMinFunction :: FuncDeriv (const Vector & x, const Vector & dir, double & deriv) const { - Vec3d n, vgrad; + Vec3d vgrad; Point3d pp1; double g1x, g1y; double badness, hbadness; @@ -368,8 +366,6 @@ vgrad = 0; badness = 0; - meshthis -> GetNormalVector (surfi, sp1, gi1, n); - pp1 = sp1; pp1.Add2 (x.Get(1), t1, x.Get(2), t2); @@ -520,7 +516,7 @@ // from 2d: int j, k, lpi, gpi; - Vec3d n, vgrad; + Vec3d vgrad; Point3d pp1; Vec2d g1, vdir; double badness, hbadness, hbad, hderiv; @@ -528,8 +524,6 @@ vgrad = 0; badness = 0; - meshthis -> GetNormalVector (surfi, sp1, gi1, n); - pp1 = sp1; pp1.Add2 (x.Get(1), t1, x.Get(2), t2); @@ -593,7 +587,7 @@ // from 2d: int j, k, lpi, gpi; - Vec3d n, vgrad; + Vec3d vgrad; Point3d pp1; Vec2d g1, vdir; double badness, hbadness, hbad, hderiv; @@ -601,8 +595,6 @@ vgrad = 0; badness = 0; - meshthis -> GetNormalVector (surfi, sp1, gi1, n); - pp1 = sp1; pp1.Add2 (x.Get(1), t1, x.Get(2), t2); @@ -859,19 +851,21 @@ locelements.SetSize(0); locrots.SetSize (0); lochs.SetSize (0); + ngi.trignum = 0; for (j = 0; j < elementsonpoint[pi].Size(); j++) { sei = elementsonpoint[pi][j]; const Element2d & bel = mesh[sei]; surfi = mesh.GetFaceDescriptor(bel.GetIndex()).SurfNr(); - + locelements.Append (sei); for (k = 1; k <= bel.GetNP(); k++) if (bel.PNum(k) == pi) { locrots.Append (k); + ngi = bel.GeomInfoPi(k); break; } @@ -942,7 +936,7 @@ } //optimizer loop (if not whole distance is not possible, move only a bit!!!!) - while (loci <= 5 && !moveisok) + while (loci <= 5 && !moveisok) { loci ++; mesh[pi].X() = origp.X() + (x.Get(1) * t1.X() + x.Get(2) * t2.X())*fact; @@ -951,11 +945,9 @@ fact = fact/2.; - ProjectPoint (surfi, mesh[pi]); + moveisok = ProjectPoint (surfi, mesh[pi], ngi); - moveisok = CalcPointGeomInfo(surfi, ngi, mesh[pi]); - // point lies on same chart in stlsurface - + // point lies on same chart in stlsurface if (moveisok) { for (j = 0; j < locelements.Size(); j++) diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Inter2d.cxx netgen-4.5_new/libsrc/occ/Partition_Inter2d.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Inter2d.cxx 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Inter2d.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -29,10 +29,10 @@ // $Header$ //using namespace std; -#include "Partition_Inter2d.ixx" - #include "utilities.h" +#include "Partition_Inter2d.ixx" + #include #include #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Inter2d.cxx.orig netgen-4.5_new/libsrc/occ/Partition_Inter2d.cxx.orig --- netgen-4.5_orig/libsrc/occ/Partition_Inter2d.cxx.orig 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/Partition_Inter2d.cxx.orig 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,677 @@ +#ifdef OCCGEOMETRY + +// GEOM PARTITION : partition algorithm +// +// Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R& D, LEG, PRINCIPIA R& D, BUREAU VERITAS +// +// 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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org +// +// +// +// File : Partition_Inter2d.cxx +// Author : Benedicte MARTIN +// Module : GEOM +// $Header$ + +//using namespace std; +#include "Partition_Inter2d.ixx" + +#include "utilities.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#ifdef DEB +static Standard_Boolean TestEdges = 0; +static Standard_Integer NbF2d = 0; +static Standard_Integer NbE2d = 0; +#endif + +//======================================================================= +//function : getOtherShape +//purpose : +//======================================================================= + +static TopoDS_Shape getOtherShape(const TopoDS_Shape& theS, + const TopTools_ListOfShape& theSList) +{ + TopTools_ListIteratorOfListOfShape anIt( theSList ); + for ( ; anIt.More(); anIt.Next() ) + if (!theS.IsSame( anIt.Value() )) + return anIt.Value(); + + return TopoDS_Shape(); +} + +//======================================================================= +//function : findVOnE +//purpose : on theE, find a vertex close to theV, such that an edge +// passing through it is an itersection of theF1 and theF2. +// theE intersects theE2 at theV +//======================================================================= + +static Standard_Boolean findVOnE(const TopoDS_Vertex & theV, + const TopoDS_Edge& theE, + const TopoDS_Edge& theE2, + const TopoDS_Shape& theF1, + const TopoDS_Shape& theF2, + const Handle(BRepAlgo_AsDes)& theAsDes, + TopoDS_Vertex & theFoundV) +{ + Standard_Real MinDist2 = ::RealLast(); + gp_Pnt P; + + // check all vertices on theE + const TopTools_ListOfShape& aVList = theAsDes->Descendant( theE ); + TopTools_ListIteratorOfListOfShape anIt( aVList ); + if (anIt.More()) + P = BRep_Tool::Pnt( theV ); + for ( ; anIt.More(); anIt.Next() ) + { + // check by distance + TopoDS_Vertex & V = TopoDS::Vertex( anIt.Value() ); + Standard_Real dist2 = P.SquareDistance( BRep_Tool::Pnt( V )); + if (dist2 < MinDist2) + MinDist2 = dist2; + else + continue; + + // V is a candidate if among edges passing through V there is one + // which is an intersection of theF1 and theF2 + TopTools_ListIteratorOfListOfShape anEIt( theAsDes->Ascendant( V )); + Standard_Boolean isOk = Standard_False; + for ( ; !isOk && anEIt.More(); anEIt.Next() ) + { + const TopoDS_Shape & E2 = anEIt.Value(); + if ( theE2.IsSame( E2 )) + continue; + const TopTools_ListOfShape & aFList = theAsDes->Ascendant( E2 ); + if (aFList.IsEmpty()) + continue; + if ( theF1.IsSame( aFList.First() )) + isOk = theF2.IsSame( aFList.Last() ); + else + isOk = theF2.IsSame( aFList.First() ) && theF1.IsSame( aFList.Last() ); + } + if (isOk) + theFoundV = V; + } + + if (theFoundV.IsNull()) + return Standard_False; + + // check that MinDist2 is not too large + Standard_Real f, l; + TopLoc_Location L; + Handle(Geom_Curve) aCurve = BRep_Tool::Curve( theE, L, f, l ); + gp_Pnt P1 = aCurve->Value( f ); + gp_Pnt P2 = aCurve->Value( 0.3 * f + 0.7 * l ); + //gp_Pnt P2 = aCurve->Value( 0.5 * ( f + l )); + if (MinDist2 > P1.SquareDistance( P2 )) + return Standard_False; + +#ifdef DEB + MESSAGE("findVOnE: found MinDist = " << sqrt (MinDist2)); +#endif + + return Standard_True; +} + +//======================================================================= +//function : AddVonE +//purpose : Put V in AsDes as intersection of E1 and E2. +// Check that vertex equal to V already exists on one +// of edges, in such a case, V is not added but +// existing vertex is updated to be on E1 and E2 and +// is returned insead of V. +//======================================================================= + +TopoDS_Vertex Partition_Inter2d::AddVonE(const TopoDS_Vertex& theV, + const TopoDS_Edge& E1, + const TopoDS_Edge& E2, + const Handle(BRepAlgo_AsDes)& AsDes, + const TopoDS_Face& theF) + +{ + //------------------------------------------------------------- + // test if the points of intersection already exist. If not, + // add as descendants of the edges. + // nb: theses points are only vertices of intersection. + //------------------------------------------------------------- + const TopTools_ListOfShape& VOnE1 = AsDes->Descendant(E1); + const TopTools_ListOfShape& VOnE2 = AsDes->Descendant(E2); + gp_Pnt P1,P2; + TopoDS_Vertex V1,V2; + TopTools_ListIteratorOfListOfShape it; + BRep_Builder B; + TopAbs_Orientation O1,O2; + Standard_Real U1,U2; + Standard_Real Tol,Tol1,Tol2; + Standard_Boolean OnE1,OnE2; + + TopoDS_Vertex V = theV; + + U1 = BRep_Tool::Parameter(V,E1); + U2 = BRep_Tool::Parameter(V,E2); + O1 = V.Orientation(); + O2 = O1; + P1 = BRep_Tool::Pnt(V); + Tol = BRep_Tool::Tolerance( V ); + OnE1 = OnE2 = Standard_False; + + //----------------------------------------------------------------- + // Search if the point of intersection is a vertex of E1. + //----------------------------------------------------------------- + for (it.Initialize(VOnE1); it.More(); it.Next()) { + const TopoDS_Vertex& CV = TopoDS::Vertex( it.Value() ); + if (V.IsSame( CV )) { + V1 = V; + OnE1 = Standard_True; + break; + } + P2 = BRep_Tool::Pnt( CV ); + Tol1 = 1.1*(Tol + BRep_Tool::Tolerance( CV )); + if (P1.SquareDistance(P2) <= Tol1*Tol1) { + V = CV; + V1 = V; + OnE1 = Standard_True; + break; + } + } + if (OnE1) { + //----------------------------------------------------------------- + // Search if the vertex found is still on E2. + //----------------------------------------------------------------- + for (it.Initialize(VOnE2); it.More(); it.Next()) { + if (V.IsSame( it.Value() )) { + OnE2 = Standard_True; + V2 = V; + break; + } + } + } + if (!OnE2) { + for (it.Initialize(VOnE2); it.More(); it.Next()) { + //----------------------------------------------------------------- + // Search if the point of intersection is a vertex of E2. + //----------------------------------------------------------------- + const TopoDS_Vertex& CV = TopoDS::Vertex( it.Value() ); + P2 = BRep_Tool::Pnt( CV ); + Tol2 = 1.1*(Tol + BRep_Tool::Tolerance( CV )); + if (P1.SquareDistance(P2) <= Tol2*Tol2) { + V = CV; + V2 = V; + OnE2 = Standard_True; + break; + } + } + } + + + if (!OnE1 && !OnE2 && !theF.IsNull()) + { + // if 3 faces intersects each others, 3 new edges on them must pass + // through one vertex but real intersection points of each + // pair of edges are sometimes more far than a tolerance. + // Try to analitically find vertices that E1 and E2 must pass trough + + TopoDS_Shape F1 = getOtherShape( theF, AsDes->Ascendant( E1 )); + TopoDS_Shape F2 = getOtherShape( theF, AsDes->Ascendant( E2 )); + if (!F1.IsNull() && !F2.IsNull() && !F1.IsSame( F2 )) + { + OnE1 = findVOnE ( theV, E1, E2, F1, F2, AsDes, V1 ); + OnE2 = findVOnE ( theV, E2, E1, F1, F2, AsDes, V2 ); + if (OnE2) V = V2; + if (OnE1) V = V1; + } + } + + if (OnE1 && OnE2) { + if (!V1.IsSame(V2)) { + // replace V1 with V2 on all edges V1 is on + Standard_Real UV1; + TopoDS_Edge EWE1; + TopoDS_Vertex VI; + const TopTools_ListOfShape& EdgeWithV1 = AsDes->Ascendant(V1); + + for (it.Initialize(EdgeWithV1); it.More(); it.Next()) { + EWE1 = TopoDS::Edge(it.Value()); + VI = V1; + VI.Orientation(TopAbs_INTERNAL); + UV1 = BRep_Tool::Parameter(VI,EWE1); + VI = V2; + VI.Orientation(TopAbs_INTERNAL); + B.UpdateVertex( VI, UV1, EWE1, GetTolerance( VI, UV1, EWE1, AsDes)); + } + AsDes->Replace(V1,V2); + V = V2; + } + } + + // add existing vertices instead of new ones + if (!OnE1) { + if (OnE2) { + V.Orientation(TopAbs_INTERNAL); + B.UpdateVertex (V, U1, E1, GetTolerance( V, U1, E1, AsDes)); + } + V.Orientation(O1); + AsDes->Add(E1,V); + } + if (!OnE2) { + if (OnE1) { + V.Orientation(TopAbs_INTERNAL); + B.UpdateVertex (V, U2, E2, GetTolerance( V, U2, E2, AsDes )); + } + V.Orientation(O2); + AsDes->Add(E2,V); + } + + return V; +} + +//======================================================================= +//function : FindEndVertex +//purpose : Returns a vertex from having parameter on +// closest to or . is True if +// found vertex is closer to . returns parameter +// difference. +//======================================================================= + +TopoDS_Vertex Partition_Inter2d::FindEndVertex(const TopTools_ListOfShape& LV, + const Standard_Real f, + const Standard_Real l, + const TopoDS_Edge& E, + Standard_Boolean& isFirst, + Standard_Real& minDU) +{ + TopoDS_Vertex endV; + Standard_Real U, endU, min; + minDU = 1.e10; + + TopTools_ListIteratorOfListOfShape it; + it.Initialize(LV); + for (; it.More(); it.Next()) { + const TopoDS_Vertex& v = TopoDS::Vertex(it.Value()); + U = BRep_Tool::Parameter(v, E); + min = Min( Abs(U-f), Abs(U-l) ); + if (min < minDU) { + endV = v; + endU = U; + minDU = min; + } + } + if (Abs(endU-f) < Abs(endU-l)) + isFirst = Standard_True; + else + isFirst = Standard_False; + + return endV; +} + +//======================================================================= +//function : treatClosed +//purpose : add second vertex to closed edge. Vertex is one of +//======================================================================= + +static void treatClosed (const TopoDS_Edge& E1, + const Standard_Real f, + const Standard_Real l, + TopTools_ListOfShape& LV1, + TopTools_ListOfShape& /*LV2*/) +{ + Standard_Boolean isFirst=0; + Standard_Real minDU = 1.e10; + TopoDS_Vertex endV; + endV = Partition_Inter2d::FindEndVertex(LV1, f,l, E1, isFirst,minDU); + + if (minDU > Precision::PConfusion()) + return; // not end point + + Standard_Real newU; + if (isFirst) + newU = f + (l - f); + else + newU = l - (l - f); + + // update end parameter + BRep_Builder B; + endV.Orientation(TopAbs_INTERNAL); + B.UpdateVertex(endV,newU,E1,BRep_Tool::Tolerance(endV)); +} + +//======================================================================= +//function : EdgesPartition +//purpose : +//======================================================================= + +static void EdgesPartition(const TopoDS_Face& F, + const TopoDS_Edge& E1, + const TopoDS_Edge& E2, + const Handle(BRepAlgo_AsDes)& AsDes, + const TopTools_MapOfShape& NewEdges, + const Standard_Boolean WithOri) +{ + + Standard_Real f[3],l[3]; + Standard_Real MilTol2; + Standard_Real Tol = Max (BRep_Tool::Tolerance(E1), + BRep_Tool::Tolerance(E2)); + MilTol2 = Tol * Tol * 10; + + BRep_Tool::Range(E1, f[1], l[1]); + BRep_Tool::Range(E2, f[2], l[2]); + + BRepAdaptor_Curve CE1(E1,F); + BRepAdaptor_Curve CE2(E2,F); + + TopoDS_Edge EI[3]; EI[1] = E1; EI[2] = E2; + TopTools_ListOfShape LV1; // new vertices at intersections on E1 + TopTools_ListOfShape LV2; // ... on E2 + BRep_Builder B; + + // if E1 and E2 are results of intersection of F and two connex faces then + // no need to intersect edges, they can contact by vertices only + // (encounted an exception in TopOpeBRep_EdgesIntersector in such a case) + Standard_Boolean intersect = Standard_True; + TopTools_IndexedMapOfShape ME; + TopExp::MapShapes(F, TopAbs_EDGE, ME); + if (!ME.Contains(E1) && ! ME.Contains(E2)) { // if E1 and E2 are new on F + TopoDS_Shape F1, F2; + const TopTools_ListOfShape& LF1 = AsDes->Ascendant( E1 ); + F1 = F.IsSame( LF1.First() ) ? LF1.Last() : LF1.First(); + const TopTools_ListOfShape& LF2 = AsDes->Ascendant( E2 ); + F2 = F.IsSame( LF2.First() ) ? LF2.Last() : LF2.First(); + if (!F.IsSame(F2) && !F.IsSame(F1) ) { + TopExp_Explorer exp(F2, TopAbs_EDGE); + TopExp::MapShapes(F1, TopAbs_EDGE, ME); + for (; exp.More(); exp.Next()) { + if (ME.Contains( exp.Current())) { + intersect = Standard_False; + break; + } + } + } + } + + if (intersect) { + //------------------------------------------------------ + // compute the points of Intersection in 2D + //----------------------------------------------------- + // i.e. fill LV1 and LV2 + TopOpeBRep_EdgesIntersector EInter; + EInter.SetFaces(F,F); + Standard_Real TolDub = 1.e-7; + EInter.ForceTolerances(TolDub,TolDub); + Standard_Boolean reducesegments = Standard_False; + EInter.Perform (E1,E2,reducesegments); + + Standard_Boolean rejectreducedsegmentpoints = Standard_False; + EInter.InitPoint(rejectreducedsegmentpoints); + for ( ; EInter.MorePoint(); EInter.NextPoint() ) + { + const TopOpeBRep_Point2d& P2D = EInter.Point(); + const gp_Pnt& P = P2D.Value(); + TopoDS_Vertex V = BRepLib_MakeVertex(P); + + //------------------------- + // control the point found. + //------------------------- + gp_Pnt P1 = CE1.Value(P2D.Parameter(1)); + gp_Pnt P2 = CE2.Value(P2D.Parameter(2)); + Standard_Real sqd1 = P1.SquareDistance(P); + Standard_Real sqd2 = P2.SquareDistance(P); + if (sqd1 > MilTol2 || sqd2 > MilTol2 ) + continue; + + // add a new vertex to the both edges + Standard_Real toler = Max( Tol, sqrt( Max( sqd1, sqd2 ))); + Standard_Integer i; + for (i = 1; i <= 2; i++) { + Standard_Real U = P2D.Parameter(i); + V.Orientation(TopAbs_INTERNAL); + B.UpdateVertex( V,U,EI[i], toler); + TopAbs_Orientation OO = TopAbs_REVERSED; + if (WithOri) { + if (P2D.IsVertex(i)) + OO = P2D.Vertex(i).Orientation(); + else if (P2D.Transition(i).Before() == TopAbs_OUT) { + OO = TopAbs_FORWARD; + } + V.Orientation(OO); + if (i == 1) LV1.Append(V); + else LV2.Append(V); + } + } + } + } // if (intersect) + + //---------------------------------- + // Test the extremities of the edges. + //---------------------------------- + // add to LV* vertices for vertex-vertex closeness + Standard_Real U1,U2; + Standard_Real TolConf2, TolConf; + TopoDS_Vertex V1[2],V2[2]; + TopExp::Vertices(E1,V1[0],V1[1]); + TopExp::Vertices(E2,V2[0],V2[1]); + + Standard_Integer i,j,k; + for (j = 0; j < 2; j++) { + if (V1[j].IsNull()) continue; + for ( k = 0; k < 2; k++) { + if (V2[k].IsNull()) continue; + gp_Pnt P1 = BRep_Tool::Pnt(V1[j]); + gp_Pnt P2 = BRep_Tool::Pnt(V2[k]); + TolConf = BRep_Tool::Tolerance(V1[j]) + BRep_Tool::Tolerance(V2[k]); + TolConf = Max (Tol, TolConf); + TolConf2 = TolConf * TolConf; + if (!intersect) + TolConf2 *= 100; + Standard_Real SqDist = P1.SquareDistance(P2); + + if (SqDist <= TolConf2) { + TopoDS_Vertex V = BRepLib_MakeVertex(P1); + V.Orientation(TopAbs_INTERNAL); + U1 = (j == 0) ? f[1] : l[1]; + U2 = (k == 0) ? f[2] : l[2]; + B.UpdateVertex(V,U1,E1,TolConf); + B.UpdateVertex(V,U2,E2,TolConf); + LV1.Prepend(V.Oriented(V1[j].Orientation())); + LV2.Prepend(V.Oriented(V2[k].Orientation())); + } + } + } + + Standard_Boolean AffichPurge = Standard_False; + + if ( LV1.IsEmpty()) return; + + //---------------------------------- + // Purge of all the vertices. + //---------------------------------- + // remove one of close vertices + TopTools_ListIteratorOfListOfShape it1LV1,it1LV2,it2LV1; + gp_Pnt P1,P2; + Standard_Boolean Purge = Standard_True; + + while (Purge) { + i = 1; + Purge = Standard_False; + for (it1LV1.Initialize(LV1),it1LV2.Initialize(LV2); + it1LV1.More(); + it1LV1.Next(),it1LV2.Next()) { + j = 1; + it2LV1.Initialize(LV1); + while (j < i) { + const TopoDS_Vertex& VE1 = TopoDS::Vertex(it1LV1.Value()); + const TopoDS_Vertex& VE2 = TopoDS::Vertex(it2LV1.Value()); + Standard_Real Tol1 = BRep_Tool::Tolerance( VE1 ); + Standard_Real Tol2 = BRep_Tool::Tolerance( VE2 ); + P1 = BRep_Tool::Pnt( VE1 ); + P2 = BRep_Tool::Pnt( VE2 ); + if (P1.IsEqual(P2, Tol1 + Tol2)) { + LV1.Remove(it1LV1); + LV2.Remove(it1LV2); + Purge = Standard_True; + break; + } + j++; + it2LV1.Next(); + } + if (Purge) break; + i++; + } + } + + // care of new closed edges, they always intersect with seam at end + if (V1[0].IsSame( V1[1] ) && NewEdges.Contains(E1) ) + treatClosed (E1, f[1], l[1], LV1, LV2); + if (V2[0].IsSame( V2[1] ) && NewEdges.Contains(E2) ) + treatClosed (E2, f[2], l[2], LV2, LV1); + + //---------------- + // Stocking vertex + //---------------- + + for ( it1LV1.Initialize( LV1 ); it1LV1.More(); it1LV1.Next()) + Partition_Inter2d::AddVonE (TopoDS::Vertex( it1LV1.Value()), + E1, E2, AsDes, F); +} + +//======================================================================= +//function : CompletPart2d +//purpose : Computes the intersections between the edges stored +// is AsDes as descendants of . Intersections is computed +// between two edges if one of them is bound in NewEdges. +//======================================================================= + +void Partition_Inter2d::CompletPart2d (const Handle(BRepAlgo_AsDes)& AsDes, + const TopoDS_Face& F, + const TopTools_MapOfShape& NewEdges) +{ + +#ifdef DEB + NbF2d++; + NbE2d = 0; +#endif + + //Do not intersect the edges of a face + TopTools_IndexedMapOfShape EdgesOfFace; + TopExp::MapShapes( F, TopAbs_EDGE , EdgesOfFace); + + //------------------------------------------------------------------- + // compute the intersection2D on the faces touched by the intersection3D + //------------------------------------------------------------------- + TopTools_ListIteratorOfListOfShape it1LE ; + TopTools_ListIteratorOfListOfShape it2LE ; + + //----------------------------------------------- + // Intersection edge-edge. + //----------------------------------------------- + const TopTools_ListOfShape& LE = AsDes->Descendant(F); + TopoDS_Vertex V1,V2; + Standard_Integer j, i = 1; + + TopoDS_Face FF = F; + FF.Orientation(TopAbs_FORWARD); + + for ( it1LE.Initialize(LE) ; it1LE.More(); it1LE.Next()) { + const TopoDS_Edge& E1 = TopoDS::Edge(it1LE.Value()); + j = 1; + it2LE.Initialize(LE); + + while (j < i && it2LE.More()) { + const TopoDS_Edge& E2 = TopoDS::Edge(it2LE.Value()); + //---------------------------------------------------------- + // Intersections of the new edges obtained by intersection + // between them and with the restrictions edges + //---------------------------------------------------------- + if ( (!EdgesOfFace.Contains(E1) || !EdgesOfFace.Contains(E2)) && + (NewEdges.Contains(E1) || NewEdges.Contains(E2)) ) { + EdgesPartition(FF,E1,E2,AsDes,NewEdges,Standard_True); + } + it2LE.Next(); + j++; + } + i++; + } +} + +//======================================================================= +//function : GetTolerance +//purpose : Returns tolerance theV must have atfer its +// addition to theE with theU parameter. theAsDes is +// used to find pcurves of theE +//======================================================================= + +Standard_Real Partition_Inter2d::GetTolerance + (const TopoDS_Vertex & theV, + const Standard_Real theU, + const TopoDS_Edge & theE, + const Handle(BRepAlgo_AsDes)& theAsDes) +{ + Standard_Real aTol = BRep_Tool::Tolerance( theV ); + gp_Pnt aPnt = BRep_Tool::Pnt( theV ); + + // check point on 3D curve + Standard_Real f,l; + Handle(Geom_Curve) C = BRep_Tool::Curve( theE, f, l ); + if (!C.IsNull()) + aTol = Max ( aTol, aPnt.Distance( C->Value( theU ))); + + // check points on pcurves + const TopTools_ListOfShape& aFList = theAsDes->Ascendant( theE ); + TopTools_ListIteratorOfListOfShape aFIt( aFList ); + for ( ; aFIt.More(); aFIt.Next() ) + { + const TopoDS_Face& F = TopoDS::Face( aFIt.Value() ); + Handle(Geom2d_Curve) pcurve = BRep_Tool::CurveOnSurface( theE, F, f, l ); + if (!pcurve.IsNull()) + { + gp_Pnt2d aPnt2d = pcurve->Value( theU ); + TopLoc_Location L; + Handle(Geom_Surface) S = BRep_Tool::Surface( F, L ); + gp_Pnt aPntOnS = S->Value( aPnt2d.X(), aPnt2d.Y() ); + if (!L.IsIdentity()) + aPntOnS.Transform( L.Transformation() ); + aTol = Max ( aTol, aPnt.Distance( aPntOnS )); + } + } + + return aTol; +} + +#endif diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Inter3d.cxx netgen-4.5_new/libsrc/occ/Partition_Inter3d.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Inter3d.cxx 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Inter3d.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -29,13 +29,17 @@ // $Header$ //using namespace std; + +#include "utilities.h" + #include "Partition_Inter2d.hxx" #include "Partition_Inter3d.ixx" -#include "utilities.h" #include #include #include +//using namespace std; + #include #include #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Inter3d.cxx.orig netgen-4.5_new/libsrc/occ/Partition_Inter3d.cxx.orig --- netgen-4.5_orig/libsrc/occ/Partition_Inter3d.cxx.orig 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/Partition_Inter3d.cxx.orig 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,943 @@ +#ifdef OCCGEOMETRY + +// GEOM PARTITION : partition algorithm +// +// Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// 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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org +// +// +// +// File : Partition_Inter3d.cxx +// Author : Benedicte MARTIN +// Module : GEOM +// $Header$ + +//using namespace std; +#include "Partition_Inter2d.hxx" +#include "Partition_Inter3d.ixx" +#include "utilities.h" + +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#ifdef DEB +#include +#endif + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +//======================================================================= +//function : Partition_Inter3d +//purpose : +//======================================================================= + +Partition_Inter3d::Partition_Inter3d() +{ +} +//======================================================================= +//function : Partition_Inter3d +//purpose : +//======================================================================= + +Partition_Inter3d::Partition_Inter3d(const Handle(BRepAlgo_AsDes)& AsDes) + :myAsDes(AsDes) +{ + mySectionEdgesAD = new BRepAlgo_AsDes; +} + +//======================================================================= +//function : CompletPart3d +//purpose : FaceShapeMap is just to know the shape a face belongs to +//======================================================================= + +void Partition_Inter3d::CompletPart3d(const TopTools_ListOfShape& SetOfFaces1, + const TopTools_DataMapOfShapeShape& FaceShapeMap) +{ + if (myAsDes.IsNull()) + myAsDes = new BRepAlgo_AsDes; + + TopTools_ListIteratorOfListOfShape it; + + //--------------------------------------------------------------- + // Construction of bounding boxes. + //--------------------------------------------------------------- + + BRep_Builder B; + TopoDS_Compound CompOS; + B.MakeCompound(CompOS); + for (it.Initialize(SetOfFaces1); it.More(); it.Next()) + B.Add(CompOS, it.Value()); + + TopOpeBRepTool_BoxSort BOS; + BOS.AddBoxesMakeCOB(CompOS,TopAbs_FACE); + + for (it.Initialize(SetOfFaces1); it.More(); it.Next()) { + TopoDS_Face F1 = TopoDS::Face(it.Value()); + + // avoid intersecting faces of one shape + TopoDS_Shape S1; + if (FaceShapeMap.IsBound(F1)) S1 = FaceShapeMap.Find(F1); + + // to filter faces sharing an edge + TopTools_IndexedMapOfShape EM; + TopExp::MapShapes( F1, TopAbs_EDGE, EM); + + TColStd_ListIteratorOfListOfInteger itLI = BOS.Compare(F1); + for (; itLI.More(); itLI.Next()) { + TopoDS_Face F2 = TopoDS::Face(BOS.TouchedShape(itLI)); + if (F1.IsSame(F2) || IsDone(F1,F2)) + continue; + + TopoDS_Shape S2; + if (FaceShapeMap.IsBound(F2)) S2 = FaceShapeMap.Find(F2); + if (!S1.IsNull() && S1.IsSame(S2)) + continue; // descendants of one shape + + TopExp_Explorer expE (F2, TopAbs_EDGE); + for ( ; expE.More(); expE.Next()) + if (EM.Contains( expE.Current() )) + break; + if (expE.More()) + { + // faces have a common edge, check if they are a tool and a face + // generated by the tool in another shape; in that case they are + // to be intersected + TopLoc_Location L1, L2; + Handle(Geom_Surface) S1 = BRep_Tool::Surface( F1, L1 ); + Handle(Geom_Surface) S2 = BRep_Tool::Surface( F2, L2 ); + if ( S1 != S2 || L1 != L2 ) + continue; + } + + F1.Orientation(TopAbs_FORWARD); + F2.Orientation(TopAbs_FORWARD); + FacesPartition(F1,F2); + } + + // mark as modified a face which has at least one new edge + if (!myAsDes->HasDescendant( F1 )) + continue; + TopTools_ListIteratorOfListOfShape itE (myAsDes->Descendant( F1 )); + for ( ; itE.More(); itE.Next()) { + if (myNewEdges.Contains( itE.Value())) { + myTouched.Add( F1 ); + break; + } + } + } +} + +//======================================================================= +//function : PutInBounds +//purpose : +//======================================================================= + +static void PutInBounds (const TopoDS_Face& F, + const TopoDS_Edge& E, + Handle(Geom2d_Curve)& C2d) +{ + Standard_Real umin,umax,vmin,vmax; + Standard_Real f,l; + BRep_Tool::Range(E,f,l); + + TopLoc_Location L; // Recup S avec la location pour eviter la copie. + Handle (Geom_Surface) S = BRep_Tool::Surface(F,L); + + if (S->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) { + S = (*(Handle_Geom_RectangularTrimmedSurface*)&S)->BasisSurface(); + } + if (!S->IsUPeriodic() && !S->IsVPeriodic()) + return; + + BRepTools::UVBounds(F,umin,umax,vmin,vmax); + + gp_Pnt2d Pf = C2d->Value(f); + gp_Pnt2d Pl = C2d->Value(l); + const Standard_Real Um = 0.34*f + 0.66*l; + gp_Pnt2d Pm = C2d->Value( Um ); + + // sometimes on shpere, pcurve is out of domain by V though S is + // UPeriodic, sometimes it is in domain but nontheless it has + // wrong position. + // Check pcurve position by 3D point + if (S->IsKind(STANDARD_TYPE( Geom_SphericalSurface ))) + { + // get point on the surface + gp_Pnt Ps = S->Value( Pm.X(), Pm.Y() ); + // get point on the edge + Handle(Geom_Curve) C = BRep_Tool::Curve( E, f, l ); + gp_Pnt Pc = C->Value( Um ); + // compare points + Standard_Real TolE = BRep_Tool::Tolerance( E ); + if ( Pc.SquareDistance( Ps ) * 0.95 < TolE * TolE ) + return; // OK + + // find good UV for Pc: project Pc on S + GeomAdaptor_Surface SA (S); + Extrema_ExtPS anExtPS (Pc, SA, + SA.UResolution( TolE ), SA.VResolution( TolE )); + if (anExtPS.IsDone()) + { + Standard_Integer i, nbExt = anExtPS.NbExt(); + Extrema_POnSurf aPOnSurf; + for (i = 1; i <= nbExt; ++i ) + if (anExtPS.Value( i ) <= TolE) { + aPOnSurf = anExtPS.Point( i ); + break; + } + if (i <= nbExt) { + // a point found + Standard_Real u, v; + aPOnSurf.Parameter( u, v ); + gp_Pnt2d aGoodPm ( u, v ); + C2d->Translate( Pm , aGoodPm ); + } + } + } + + //--------------- + // Recadre en U. + //--------------- + if (S->IsUPeriodic()) { + Standard_Real period = S->UPeriod(); + Standard_Real eps = period*1.e-6; + Standard_Real minC = Min(Pf.X(),Pl.X()); minC = Min(minC,Pm.X()); + Standard_Real maxC = Max(Pf.X(),Pl.X()); maxC = Max(maxC,Pm.X()); + Standard_Real du = 0.; + if (minC< umin - eps) { + du = (int((umin - minC)/period) + 1)*period; + } + if (minC > umax + eps) { + du = -(int((minC - umax)/period) + 1)*period; + } + if (du != 0) { + gp_Vec2d T1(du,0.); + C2d->Translate(T1); + minC += du; maxC += du; + } + // Ajuste au mieux la courbe dans le domaine. + if (maxC > umax +100*eps) { + Standard_Real d1 = maxC - umax; + Standard_Real d2 = umin - minC + period; + if (d2 < d1) du =-period; + if ( du != 0.) { + gp_Vec2d T2(du,0.); + C2d->Translate(T2); + } + } + } + //------------------ + // Recadre en V. + //------------------ + if (S->IsVPeriodic()) { + Standard_Real period = S->VPeriod(); + Standard_Real eps = period*1.e-6; + Standard_Real minC = Min(Pf.Y(),Pl.Y()); minC = Min(minC,Pm.Y()); + Standard_Real maxC = Max(Pf.Y(),Pl.Y()); maxC = Max(maxC,Pm.Y()); + Standard_Real dv = 0.; + if (minC< vmin - eps) { + dv = (int((vmin - minC)/period) + 1)*period; + } + if (minC > vmax + eps) { + dv = -(int((minC - vmax)/period) + 1)*period; + } + if (dv != 0) { + gp_Vec2d T1(0.,dv); + C2d->Translate(T1); + minC += dv; maxC += dv; + } + // Ajuste au mieux la courbe dans le domaine. + if (maxC > vmax +100*eps) { + Standard_Real d1 = maxC - vmax; + Standard_Real d2 = vmin - minC + period; + if (d2 < d1) dv =-period; + if ( dv != 0.) { + gp_Vec2d T2(0.,dv); + C2d->Translate(T2); + } + } + } +} + +//======================================================================= +//function : Inter3D +//purpose : +//======================================================================= + +void Partition_Inter3d::Inter3D(const TopoDS_Face& F1, + const TopoDS_Face& F2, + TopTools_ListOfShape& L) +{ + BRep_Builder B; + + // fill the data Structure + Handle(TopOpeBRepDS_HDataStructure) DatStr = new TopOpeBRepDS_HDataStructure(); + TopOpeBRep_DSFiller DSFiller; + DSFiller.Insert(F1,F2,DatStr); + + // define the GeomTool used by the DSFiller : + // compute BSpline of degree 1 on intersection curves. + Standard_Real tol3dAPPROX = 1e-7; + Standard_Real tol2dAPPROX = 1e-7; + TopOpeBRepTool_GeomTool GT2 (TopOpeBRepTool_APPROX); + GT2.SetTolerances(tol3dAPPROX,tol2dAPPROX); + TopOpeBRepDS_BuildTool BT(GT2); + + // Perform Section + TopOpeBRepBuild_Builder TopB(BT); + TopB.Perform(DatStr); + + // =============== + // Store new edges + // =============== + + L.Clear(); + TopOpeBRepDS_CurveExplorer cex(DatStr->DS()); + for (; cex.More(); cex.Next()) { + const TopOpeBRepDS_Curve& CDS = cex.Curve(); + Standard_Integer ic = cex.Index(); + Handle(Geom2d_Curve) pc1 = CDS.Curve1(); + Handle(Geom2d_Curve) pc2 = CDS.Curve2(); + + TopTools_ListIteratorOfListOfShape itLE = TopB.NewEdges(ic); + while (itLE.More()) { + TopoDS_Edge E = TopoDS::Edge(itLE.Value()); + + PutInBounds (F1,E,pc1); + PutInBounds (F2,E,pc2); + + B.UpdateEdge (E,pc1,F1,0.); + B.UpdateEdge (E,pc2,F2,0.); + + L.Append (E); + + itLE.Next(); + if (itLE.More()) { + pc1 = Handle(Geom2d_Curve)::DownCast(pc1->Copy()); + pc2 = Handle(Geom2d_Curve)::DownCast(pc2->Copy()); + } + } + } + + // ======================== + // store same domain faces + // ======================== + + + if ( DatStr->HasSameDomain( F1 )) + { + TopTools_ListOfShape emptyList; + if (!mySameDomainFM.IsBound(F1)) + mySameDomainFM.Bind(F1,emptyList); + if (!mySameDomainFM.IsBound(F2)) + mySameDomainFM.Bind(F2,emptyList); + mySameDomainFM(F1).Append(F2); + mySameDomainFM(F2).Append(F1); + } + + // ==================== + // Store section edges + // ==================== + + const TopOpeBRepDS_DataStructure& DS = DatStr->DS(); + Standard_Integer j,i,nse = DS.NbSectionEdges(); + if (nse == 0) return; + + + TopoDS_Vertex V, sdeV1, sdeV2; + TopTools_MapOfShape MV; + TopTools_ListOfShape LSE; // list of section edges + TopoDS_Face dummyF; + + for (i = 1; i <= nse; i++) + { + const TopoDS_Edge & se = DS.SectionEdge(i); + if (! TopB.IsSplit(se,TopAbs_ON)) + continue; + LSE.Append( se ); + + // add vertices where section edges interferes with other + // edges as its descendant in myAsDes + + TopoDS_Edge sde, oe; // same domain, other edge + if (DatStr->HasSameDomain(se)) { + sde = TopoDS::Edge( DatStr->SameDomain(se).Value() ); + TopExp::Vertices( sde, sdeV1, sdeV2); + } + TColStd_MapOfInteger MIV; // indices of added edges + TopOpeBRepDS_PointIterator itP (DS.ShapeInterferences( se )); + itP.SupportKind( TopOpeBRepDS_EDGE ); + // loop on intersections of se + for (; itP.More(); itP.Next()) { + oe = TopoDS::Edge( DS.Shape( itP.Support())); + if (itP.IsVertex()) { + // there is a vertex at intersection + if ( !MIV.Add( itP.Current() )) + continue; + V = TopoDS::Vertex( DS.Shape( itP.Current())); + if ( !sde.IsNull() && (V.IsSame(sdeV1) || V.IsSame(sdeV2)) ) + oe = sde; + V = ReplaceSameDomainV( V , oe ); + V.Orientation( TopAbs_INTERNAL); + B.UpdateVertex( V, itP.Parameter(), se, 0.); // AddVonE() sets real U + } + else { + // create a new vertex at the intersection point + const TopOpeBRepDS_Point& DSP = DS.Point( itP.Current()); + V = BRepLib_MakeVertex( DSP.Point() ); + V.Orientation( TopAbs_INTERNAL); + B.UpdateVertex( V, itP.Parameter(), se, DSP.Tolerance()); + // make V be on the other edge + TopOpeBRepDS_PointIterator itOP (DS.ShapeInterferences( oe )); + for (; itOP.More(); itOP.Next()) { + const TopOpeBRepDS_Point& ODSP = DS.Point( itOP.Current()); + if ( DSP.IsEqual (ODSP)) { + B.UpdateVertex( V, itOP.Parameter(), TopoDS::Edge(oe), ODSP.Tolerance()); + break; + } + } + } + // add V on the both intersecting edges + TopoDS_Vertex addedV = Partition_Inter2d::AddVonE( V,se,oe,myAsDes,dummyF); + if (!addedV.IsSame( V )) + mySameDomainVM.Bind (V, addedV); // equal vertex is already there + + MV.Add( addedV ); // to ease storage of vertices of ON splits + } + } + + // add section edge to the face it intersects and find + // splits ON that do not have same domain pair + + TopB.SplitSectionEdges(); // let TopB find ON splits + + TopTools_MapOfShape SPM; // map of ON splits + TopTools_IndexedMapOfShape ME[2]; + TopExp::MapShapes( F1, TopAbs_EDGE, ME[1]); + TopExp::MapShapes( F2, TopAbs_EDGE, ME[0]); + + TopTools_ListIteratorOfListOfShape itSP, itLSE (LSE); + while ( itLSE.More() ) { + + TopoDS_Edge se = TopoDS::Edge( itLSE.Value() ); + + // move itLSE to the next se + Standard_Integer ancRank = DS.AncestorRank(se); + if (ME[ancRank-1].Contains( se )) + { + LSE.Remove( itLSE ); // se is an edge of face it intersects + continue; + } + else + { + itLSE.Next(); + } + + const TopoDS_Face& F = (ancRank == 1) ? F2 : F1; + + // add se to face but dont add twice + TopTools_ListIteratorOfListOfShape itE( myAsDes->Descendant( F )); + if (myAsDes->HasDescendant( F )) { + for ( ; itE.More(); itE.Next()) + if (se.IsSame( itE.Value() )) + break; + } + if (!itE.More()) + { + myAsDes->Add( F, se ); + + // check se pcurve on F + Standard_Real tol, f,l, umin=1e100, umax=-1e100; + Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( se, F, f,l); + if (pc.IsNull()) { + itSP.Initialize( TopB.Splits(se,TopAbs_ON) ); + for ( ; itSP.More(); itSP.Next()) { + const TopoDS_Edge& E = TopoDS::Edge ( itSP.Value()); + BRep_Tool::Range(E, f, l); + umin = Min( umin, f); + umax = Max( umax, l); + } + Handle(Geom_Curve) C3d = BRep_Tool::Curve( se, f, l); + if (umin < umax) // sometimes umin == umax for closed edge + C3d = new Geom_TrimmedCurve( C3d, umin, umax); + pc = TopOpeBRepTool_CurveTool::MakePCurveOnFace (F,C3d,tol); + if (pc.IsNull()) { + MESSAGE (" CANT BUILD PCURVE "); + } + B.UpdateEdge( se, pc, F, tol); + } + } + + // to detect splits that do not have same domain pair + // ie which split a face into parts and not pass by its boundary + itSP.Initialize( TopB.Splits(se,TopAbs_ON) ); + for ( ; itSP.More(); itSP.Next()) { + const TopoDS_Shape& SP = itSP.Value(); + if (!SPM.Add( SP )) + SPM.Remove( SP ); + } + } + + // store vertices of ON splits and bind section edges to faces + + for (itLSE.Initialize (LSE); itLSE.More(); itLSE.Next()) + { + const TopoDS_Shape& se = itLSE.Value(); + + Standard_Integer ancRank = DS.AncestorRank(se); + TopoDS_Face F = (ancRank == 1) ? F2 : F1; + + // add vertices of ON splits which have no same domain pair + Standard_Boolean added = Standard_False; + itSP.Initialize( TopB.Splits(se,TopAbs_ON) ); + for ( ; itSP.More(); itSP.Next()) + { + if (!SPM.Contains( itSP.Value() )) + continue; + + const TopoDS_Edge& S = TopoDS::Edge ( itSP.Value()); + + added = Standard_True; + mySectionEdgesAD->Add( F, se ); + + TopoDS_Vertex VS[2]; + TopExp::Vertices (S, VS[0], VS[1]); + for (j=0; j<2; ++j) + { + if (mySameDomainVM.IsBound( VS[j] )) + VS[j] = TopoDS::Vertex( mySameDomainVM( VS[j] )); + if ( !MV.Contains( VS[j] )) { + // find equal vertex on se - point interference + gp_Pnt P1 = BRep_Tool::Pnt( VS[j] ); + TopTools_ListIteratorOfListOfShape itV( myAsDes->Descendant(se) ); + for (; itV.More(); itV.Next()) { + V = TopoDS::Vertex( itV.Value() ); + if ( V.IsSame( VS[j] )) + break; + gp_Pnt P2 = BRep_Tool::Pnt( V ); + if (P1.IsEqual( P2, Precision::Confusion())) { + mySameDomainVM.Bind (VS[j], V); + VS[j] = V; + break; + } + } + if (!itV.More()) // no interferences with edges + myAsDes->Add( se, VS[j]); + } + + // add ends of ON splits to F in order to detect later + // if a split is on face in IsSplitOn() + mySectionEdgesAD->Add( F, VS[j]); + } + // in the descendants of F, first go ends of an ON split and + // then a split itself + mySectionEdgesAD->Add( F, S ); + } + if (!added) + mySectionEdgesAD->Add( F, se ); + + myNewEdges.Add( se ); + } +} + +//======================================================================= +//function : FacesPartition +//purpose : +//======================================================================= + +void Partition_Inter3d::FacesPartition(const TopoDS_Face& F1, + const TopoDS_Face& F2) + //(const TopTools_DataMapOfShapeListOfShape& /*SetOfFaces2*/) +{ + TopTools_ListOfShape LInt; + + Inter3D (F1,F2,LInt); + + StorePart3d (F1,F2,LInt); +} + +//======================================================================= +//function : SetDone +//purpose : +//======================================================================= + +void Partition_Inter3d::SetDone(const TopoDS_Face& F1, + const TopoDS_Face& F2) +{ + if (!myDone.IsBound(F1)) { + TopTools_ListOfShape emptyList; + myDone.Bind(F1,emptyList); + } + myDone(F1).Append(F2); + if (!myDone.IsBound(F2)) { + TopTools_ListOfShape emptyList; + myDone.Bind(F2,emptyList); + } + myDone(F2).Append(F1); +} + +//======================================================================= +//function : IsDone +//purpose : +//======================================================================= + +Standard_Boolean Partition_Inter3d::IsDone(const TopoDS_Face& F1, + const TopoDS_Face& F2) + + const +{ + if (myDone.IsBound(F1)) { + TopTools_ListIteratorOfListOfShape it (myDone(F1)); + for (; it.More(); it.Next()) { + if (it.Value().IsSame(F2)) return Standard_True; + } + } + return Standard_False; +} + +//======================================================================= +//function : StorePart3d +//purpose : +//======================================================================= + +void Partition_Inter3d::StorePart3d(const TopoDS_Face& F1, + const TopoDS_Face& F2, + const TopTools_ListOfShape& LInt) +{ + if (!LInt.IsEmpty()) { + myAsDes->Add( F1,LInt); + myAsDes->Add( F2,LInt); + + TopTools_ListIteratorOfListOfShape it(LInt); + for (; it.More(); it.Next()) { + + TopoDS_Edge E = TopoDS::Edge(it.Value()); + + BRep_Builder B; + B.SameParameter(E,Standard_False); + BRepLib::SameParameter(E,1.0e-7); + + myNewEdges.Add(E); + } + } + SetDone(F1,F2); +} + +//======================================================================= +//function : TouchedFaces +//purpose : +//======================================================================= + +TopTools_MapOfShape& Partition_Inter3d::TouchedFaces() +{ + return myTouched; +} + +//======================================================================= +//function : AsDes +//purpose : +//======================================================================= + +Handle(BRepAlgo_AsDes) Partition_Inter3d::AsDes() const +{ + return myAsDes; +} + +//======================================================================= +//function : NewEdges +//purpose : +//======================================================================= + +TopTools_MapOfShape& Partition_Inter3d::NewEdges() +{ + return myNewEdges; +} + +//======================================================================= +//function : Affiche +//purpose : +//======================================================================= + +void Partition_Inter3d::Affiche(const TopTools_ListOfShape& SetOfFaces) const +{ +#ifdef DEB + char PSection[1024]; + char *section=PSection; + Standard_Integer i = 0; + Standard_Real j=1; + TopTools_ListOfShape aList; + TopTools_ListIteratorOfListOfShape it; + for (it.Initialize(SetOfFaces); it.More(); it.Next()) { + const TopoDS_Shape& OS = it.Value(); + aList=myAsDes->Descendant(OS); + MESSAGE ( " the number of items stored in the list " << j << " : " << aList.Extent() ) + j++; + TopTools_ListIteratorOfListOfShape itaList; + for (itaList.Initialize(aList); itaList.More(); itaList.Next()) { + const TopoDS_Shape& SS = itaList.Value(); + i++; + sprintf(PSection,"section_%d",i); + DBRep::Set(section,SS); + } + } +#endif +} + +//======================================================================= +//function : SameDomain +//purpose : +//======================================================================= + +const TopTools_ListOfShape& Partition_Inter3d::SameDomain(const TopoDS_Face& F) const +{ + if (mySameDomainFM.IsBound( F )) + return mySameDomainFM (F); + + static TopTools_ListOfShape emptyList; + return emptyList; +} + +//======================================================================= +//function : HasSameDomainF +//purpose : Return true if F has same domain faces +//======================================================================= + +Standard_Boolean Partition_Inter3d::HasSameDomainF(const TopoDS_Shape& F) const +{ + return mySameDomainFM.IsBound( F ); +} + +//======================================================================= +//function : IsSameDomain +//purpose : Return true if F1 and F2 are same domain faces +//======================================================================= + +Standard_Boolean Partition_Inter3d::IsSameDomainF(const TopoDS_Shape& F1, + const TopoDS_Shape& F2) const +{ + if (mySameDomainFM.IsBound( F1 )) { + TopTools_ListIteratorOfListOfShape it (mySameDomainFM( F1 )); + for (; it.More(); it.Next()) + if (F2.IsSame( it.Value())) + return Standard_True; + } + return F1.IsSame( F2 ); +} + +//======================================================================= +//function : ReplaceSameDomainV +//purpose : return same domain vertex of V if it was replaced +// and make this vertex to be on E too, else return V +//======================================================================= + +TopoDS_Vertex Partition_Inter3d::ReplaceSameDomainV(const TopoDS_Vertex& V, + const TopoDS_Edge& E) const +{ + TopoDS_Vertex SDV = V; + if (mySameDomainVM.IsBound( V )) { + + TopoDS_Vertex V1,V2; + TopExp::Vertices(E,V1,V2); + Standard_Boolean isClosed = V1.IsSame( V2 ) && V.IsSame(V1); + + SDV = TopoDS::Vertex( mySameDomainVM(V) ); + Standard_Real tol = BRep_Tool::Tolerance( V ); + BRep_Builder B; + SDV.Orientation( V.Orientation()); + + if (isClosed) { + Standard_Real f, l; + BRep_Tool::Range (E, f, l); + Standard_Boolean isFirst = IsEqual( BRep_Tool::Parameter(V,E), f ); + B.UpdateVertex(SDV, (isFirst ? f : l), E, tol); + SDV.Reverse(); + B.UpdateVertex(SDV, (isFirst ? l : f), E, tol); + } + else + B.UpdateVertex (SDV, BRep_Tool::Parameter(V,E), E, tol); + + } + return SDV; +} + +//======================================================================= +//function : SectionEdgesAD +//purpose : +//======================================================================= + +Handle(BRepAlgo_AsDes) Partition_Inter3d::SectionEdgesAD() const +{ + return mySectionEdgesAD; +} + +//======================================================================= +//function : IsSectionEdge +//purpose : return True if E is an edge of a face and it +// intersects an other face +//======================================================================= + +Standard_Boolean + Partition_Inter3d::IsSectionEdge(const TopoDS_Edge& E) const +{ + return mySectionEdgesAD->HasAscendant(E); +} + +//======================================================================= +//function : HasSectionEdge +//purpose : return True if an edge of F intersects an other +// face or F is intersected by edge of an other face +//======================================================================= + +Standard_Boolean + Partition_Inter3d::HasSectionEdge(const TopoDS_Face& F) const +{ + return mySectionEdgesAD->HasDescendant(F); +} + +//======================================================================= +//function : IsSplitOn +//purpose : return True if NewE is split of OldE on F +//======================================================================= + +Standard_Boolean + Partition_Inter3d::IsSplitOn(const TopoDS_Edge& NewE, + const TopoDS_Edge& OldE, + const TopoDS_Face& F) const +{ + if (! mySectionEdgesAD->HasDescendant(F)) + return Standard_False; + + TopTools_ListIteratorOfListOfShape itE ( mySectionEdgesAD->Descendant(F) ); + for ( ; itE.More(); itE.Next()) { + if ( itE.Value().ShapeType() != TopAbs_EDGE || + ! OldE.IsSame ( itE.Value() )) + continue; + // an edge encountered, its vertices and a split come next + itE.Next(); + if (!itE.More()) break; + const TopoDS_Shape& V3 = itE.Value(); + if (V3.ShapeType() != TopAbs_VERTEX) continue; + itE.Next(); + if (!itE.More()) break; + const TopoDS_Shape& V4 = itE.Value(); + if (V4.ShapeType() != TopAbs_VERTEX) continue; + + TopoDS_Vertex V1, V2; + TopExp::Vertices( OldE, V1, V2); + + if ( V1.IsSame(V2) && + (V1.IsSame(V3) || V1.IsSame(V4)) ) { + // closed old edge; use the split for the test + itE.Next(); + if (!itE.More()) break; + const TopoDS_Edge& split = TopoDS::Edge( itE.Value() ); + // check distance at middle point of NewE + Standard_Real f1,l1, f2,l2; + Handle(Geom2d_Curve) PC1 = BRep_Tool::CurveOnSurface( split, F ,f1,l1); + if (!PC1.IsNull()) { + Handle(Geom2d_Curve) PC2 = BRep_Tool::CurveOnSurface(NewE, F ,f2,l2); + gp_Pnt2d P = PC2->Value( 0.5*(f2+l2) ); + Geom2dAPI_ProjectPointOnCurve proj (P, PC1, f1, l1); + if (proj.NbPoints() && + proj.LowerDistance() <= Precision::Confusion()) + return Standard_True; + } + else { + Handle(Geom_Curve) C1 = BRep_Tool::Curve( split ,f1,l1); + Handle(Geom_Curve) C2 = BRep_Tool::Curve( NewE ,f2,l2); + gp_Pnt P = C2->Value( 0.5*(f2+l2) ); + GeomAPI_ProjectPointOnCurve proj (P, C1, f1, l1); + if (proj.NbPoints() && + proj.LowerDistance() <= Precision::Confusion()) + return Standard_True; + } + } + else { + Standard_Real u3 = BRep_Tool::Parameter( TopoDS::Vertex(V3), OldE); + Standard_Real u4 = BRep_Tool::Parameter( TopoDS::Vertex(V4), OldE); + + Standard_Real f,l, u; + BRep_Tool::Range( NewE, f,l); + u = 0.5*(f+l); + f = Min(u3,u4); + l = Max(u3,u4); + + if (u <= l && u >= f) + return Standard_True; + } + } + return Standard_False; +} + +//======================================================================= +//function : SectionEdgeFaces +//purpose : return faces cut by section edge +//======================================================================= + +const TopTools_ListOfShape& + Partition_Inter3d::SectionEdgeFaces(const TopoDS_Edge& SecE) const +{ + return mySectionEdgesAD->Ascendant( SecE ); +} + +#endif diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Loop.cxx netgen-4.5_new/libsrc/occ/Partition_Loop.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Loop.cxx 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Loop.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -29,12 +29,14 @@ // $Header$ //using namespace std; -#include -#include "Partition_Loop.ixx" #include "utilities.h" +#include + +#include "Partition_Loop.ixx" + #include #include #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Loop.cxx.orig netgen-4.5_new/libsrc/occ/Partition_Loop.cxx.orig --- netgen-4.5_orig/libsrc/occ/Partition_Loop.cxx.orig 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/Partition_Loop.cxx.orig 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,472 @@ +#ifdef OCCGEOMETRY + +// GEOM PARTITION : partition algorithm +// +// Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// 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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org +// +// +// +// File : Partition_Loop.cxx +// Author : Benedicte MARTIN +// Module : GEOM +// $Header$ + +//using namespace std; +#include + +#include "Partition_Loop.ixx" + +#include "utilities.h" + +#include +#include +#include + +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include + +static char* name = new char[100]; +static int nbe = 0; + +//======================================================================= +//function : Partition_Loop +//purpose : +//======================================================================= +Partition_Loop::Partition_Loop() +{ +} + +//======================================================================= +//function : Init +//purpose : +//======================================================================= +void Partition_Loop::Init(const TopoDS_Face& F) +{ + myConstEdges.Clear(); + myNewWires .Clear(); + myNewFaces .Clear(); + myFace = F; +} + +//======================================================================= +//function : AddConstEdge +//purpose : +//======================================================================= +void Partition_Loop::AddConstEdge (const TopoDS_Edge& E) +{ + myConstEdges.Append(E); +} + + +//======================================================================= +//function : FindDelta +//purpose : +//======================================================================= +static Standard_Real FindDelta(TopTools_ListOfShape& LE, + const TopoDS_Face& F) +{ + Standard_Real dist, f, l; + Standard_Real d = Precision::Infinite(); + TopTools_ListIteratorOfListOfShape itl; + + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + Handle(Geom2d_Curve) C = BRep_Tool::CurveOnSurface(E,F,f,l); + gp_Pnt2d p = C->Value(f); + gp_Pnt2d pp = C->Value(l); + Standard_Real d1 = p.Distance(pp); + if (d1 connected by the vertex in the list . +// Is erased of the list. If is too in the list +// with the same orientation, it's erased of the list +//======================================================================= +static Standard_Boolean SelectEdge(const TopoDS_Face& F, + const TopoDS_Edge& CE, + const TopoDS_Vertex& CV, + TopoDS_Edge& NE, + TopTools_ListOfShape& LE) +{ + TopTools_ListIteratorOfListOfShape itl; + NE.Nullify(); + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + if (itl.Value().IsEqual(CE)) { + LE.Remove(itl); + break; + } + } + + if (LE.Extent() > 1) { + //-------------------------------------------------------------- + // Several possible edges. + // - Test the edges differents of CE + //-------------------------------------------------------------- + Standard_Real cf, cl, f, l; + TopoDS_Face FForward = F; + Handle(Geom2d_Curve) Cc, C; + FForward.Orientation(TopAbs_FORWARD); + + Cc = BRep_Tool::CurveOnSurface(CE,FForward,cf,cl); + Standard_Real dist,distmin = 100*BRep_Tool::Tolerance(CV); + Standard_Real uc,u; + if (CE.Orientation () == TopAbs_FORWARD) uc = cl; + else uc = cf; + + gp_Pnt2d P2,PV = Cc->Value(uc); + + Standard_Real delta = FindDelta(LE,FForward); + + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + if (!E.IsSame(CE)) { + C = BRep_Tool::CurveOnSurface(E,FForward,f,l); + if (E.Orientation () == TopAbs_FORWARD) u = f; + else u = l; + P2 = C->Value(u); + dist = PV.Distance(P2); + if (dist <= distmin){ + distmin = dist; + } + + } + } + + Standard_Real anglemax = - PI; + TopoDS_Edge SelectedEdge; + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + if (!E.IsSame(CE)) { + C = BRep_Tool::CurveOnSurface(E,FForward,f,l); + if (E.Orientation () == TopAbs_FORWARD) u = f; + else u = l; + P2 = C->Value(u); + dist = PV.Distance(P2); + if (dist <= distmin + (1./3)*delta){ + gp_Pnt2d PC, P; + gp_Vec2d CTg1, CTg2, Tg1, Tg2; + Cc->D2(uc, PC, CTg1, CTg2); + C->D2(u, P, Tg1, Tg2); + + Standard_Real angle; + + if (CE.Orientation () == TopAbs_REVERSED && E.Orientation () == TopAbs_FORWARD) { + angle = CTg1.Angle(Tg1.Reversed()); + } + else if (CE.Orientation () == TopAbs_FORWARD && E.Orientation () == TopAbs_REVERSED) { + angle = (CTg1.Reversed()).Angle(Tg1); + } + else if (CE.Orientation () == TopAbs_REVERSED && E.Orientation () == TopAbs_REVERSED) { + angle = CTg1.Angle(Tg1); + } + else if (CE.Orientation () == TopAbs_FORWARD && E.Orientation () == TopAbs_FORWARD) { + angle = (CTg1.Reversed()).Angle(Tg1.Reversed()); + } + if (angle >= anglemax) { + anglemax = angle ; + SelectedEdge = E; + } + } + } + } + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + if (E.IsEqual(SelectedEdge)) { + NE = TopoDS::Edge(E); + LE.Remove(itl); + break; + } + } + } + else if (LE.Extent() == 1) { + NE = TopoDS::Edge(LE.First()); + LE.RemoveFirst(); + } + else { + return Standard_False; + } + return Standard_True; +} + +//======================================================================= +//function : SamePnt2d +//purpose : +//======================================================================= +static Standard_Boolean SamePnt2d(TopoDS_Vertex V, + TopoDS_Edge& E1, + TopoDS_Edge& E2, + TopoDS_Face& F) +{ + Standard_Real f1,f2,l1,l2; + gp_Pnt2d P1,P2; + TopoDS_Shape aLocalF = F.Oriented(TopAbs_FORWARD); + TopoDS_Face FF = TopoDS::Face(aLocalF); + Handle(Geom2d_Curve) C1 = BRep_Tool::CurveOnSurface(E1,FF,f1,l1); + Handle(Geom2d_Curve) C2 = BRep_Tool::CurveOnSurface(E2,FF,f2,l2); + if (E1.Orientation () == TopAbs_FORWARD) P1 = C1->Value(f1); + else P1 = C1->Value(l1); + + if (E2.Orientation () == TopAbs_FORWARD) P2 = C2->Value(l2); + else P2 = C2->Value(f2); + Standard_Real Tol = 100*BRep_Tool::Tolerance(V); + Standard_Real Dist = P1.Distance(P2); + return Dist < Tol; +} + +//======================================================================= +//function : PurgeNewEdges +//purpose : +//======================================================================= +static void PurgeNewEdges(TopTools_ListOfShape& ConstEdges, + const TopTools_MapOfOrientedShape& UsedEdges) +{ + TopTools_ListIteratorOfListOfShape it(ConstEdges); + while ( it.More()) { + const TopoDS_Shape& NE = it.Value(); + if (!UsedEdges.Contains(NE)) { + ConstEdges.Remove(it); + } + else { + it.Next(); + } + } +} + +//======================================================================= +//function : StoreInMVE +//purpose : +//======================================================================= +static void StoreInMVE (const TopoDS_Face& F, + TopoDS_Edge& E, + TopTools_DataMapOfShapeListOfShape& MVE ) + +{ + TopoDS_Vertex V1, V2; + TopTools_ListOfShape Empty; + + TopExp::Vertices(E,V1,V2); + if (!MVE.IsBound(V1)) { + MVE.Bind(V1,Empty); + } + MVE(V1).Append(E); + + if (!MVE.IsBound(V2)) { + MVE.Bind(V2,Empty); + } + MVE(V2).Append(E); +} + +//======================================================================= +//function : Perform +//purpose : +//======================================================================= +void Partition_Loop::Perform() +{ + + TopTools_DataMapOfShapeListOfShape MVE; + TopTools_DataMapIteratorOfDataMapOfShapeListOfShape Mapit, Mapit1; + TopTools_ListIteratorOfListOfShape itl; + TopoDS_Vertex V1,V2; + + //----------------------------------- + // Construction map vertex => edges + //----------------------------------- + for (itl.Initialize(myConstEdges); itl.More(); itl.Next()) { + TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + StoreInMVE(myFace,E,MVE); + } + + //---------------------------------------------- + // Construction of all the wires and of all the new faces. + //---------------------------------------------- + TopTools_MapOfOrientedShape UsedEdges; + + while (!MVE.IsEmpty()) { + TopoDS_Vertex VF,CV; + TopoDS_Edge CE,NE,EF; + TopoDS_Wire NW; + BRep_Builder B; + Standard_Boolean End= Standard_False; + + B.MakeWire(NW); + //-------------------------------- + // EF first edge. + //-------------------------------- + Mapit.Initialize(MVE); + EF = CE = TopoDS::Edge(Mapit.Value().First()); + + TopExp::Vertices(CE,V1,V2); + //-------------------------------- + // VF first vertex + //-------------------------------- + if (CE.Orientation() == TopAbs_FORWARD) { + CV = VF = V1; + } + else { + CV = VF = V2; + } + if (!MVE.IsBound(CV)) continue; + for ( itl.Initialize(MVE(CV)); itl.More(); itl.Next()) { + if (itl.Value().IsEqual(CE)) { + MVE(CV).Remove(itl); + break; + } + } + + int i = 0; + while (!End) { + //------------------------------- + // Construction of a wire. + //------------------------------- + TopExp::Vertices(CE,V1,V2); + if (!CV.IsSame(V1)) CV = V1; else CV = V2; + B.Add (NW,CE); + UsedEdges.Add(CE); + + //-------------- + // stop test + //-------------- + if (!MVE.IsBound(CV) || MVE(CV).IsEmpty() || CV.IsSame(VF) ) { + if (CV.IsSame(VF)) { + if (MVE(CV).Extent() == 1 ) MVE.UnBind(CV); + else { + for ( itl.Initialize(MVE(CV)); itl.More(); itl.Next()) { + if (itl.Value().IsEqual(CE)) { + MVE(CV).Remove(itl); + break; + } + } + } + } + End=Standard_True; + } + + //-------------- + // select edge + //-------------- + else { + Standard_Boolean find = SelectEdge(myFace,CE,CV,NE,MVE(CV)); + if (find) { + CE=NE; + if (MVE(CV).IsEmpty()) MVE.UnBind(CV); + if (CE.IsNull() ) { + MESSAGE ( " CE is NULL !!! " ) + End=Standard_True; + } + } + else { + MESSAGE ( " edge doesn't exist " ) + End=Standard_True; + } + } + } + + //----------------------------- + // Test if the wire is closed + //----------------------------- + if (VF.IsSame(CV) && SamePnt2d(VF,EF,CE,myFace)) { + } + else{ + MESSAGE ( "wire not closed" ) + } + myNewWires.Append (NW); + } + + PurgeNewEdges(myConstEdges,UsedEdges); + +} + + +//======================================================================= +//function : NewWires +//purpose : +//======================================================================= +const TopTools_ListOfShape& Partition_Loop::NewWires() const +{ + return myNewWires; +} + +//======================================================================= +//function : NewFaces +//purpose : +//======================================================================= +const TopTools_ListOfShape& Partition_Loop::NewFaces() const +{ + return myNewFaces; +} + +//======================================================================= +//function : WiresToFaces +//purpose : +//======================================================================= +void Partition_Loop::WiresToFaces() +{ + if (!myNewWires.IsEmpty()) { + BRepAlgo_FaceRestrictor FR; + + TopAbs_Orientation OriF = myFace.Orientation(); + TopoDS_Shape aLocalS = myFace.Oriented(TopAbs_FORWARD); + + FR.Init (TopoDS::Face(aLocalS),Standard_False); + TopTools_ListIteratorOfListOfShape it(myNewWires); + for (; it.More(); it.Next()) { + FR.Add(TopoDS::Wire(it.Value())); + } + + FR.Perform(); + + if (FR.IsDone()) { + for (; FR.More(); FR.Next()) { + myNewFaces.Append(FR.Current().Oriented(OriF)); + } + } + } +} + + +#endif diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Loop2d.cxx netgen-4.5_new/libsrc/occ/Partition_Loop2d.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Loop2d.cxx 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Loop2d.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -12,9 +12,11 @@ // $Header$ //using namespace std; -#include "Partition_Loop2d.ixx" + #include "utilities.h" + +#include "Partition_Loop2d.ixx" #include #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Loop2d.cxx.orig netgen-4.5_new/libsrc/occ/Partition_Loop2d.cxx.orig --- netgen-4.5_orig/libsrc/occ/Partition_Loop2d.cxx.orig 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/Partition_Loop2d.cxx.orig 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,1142 @@ +#ifdef OCCGEOMETRY + +// GEOM PARTITION : partition algorithm +// +// Copyright (C) 2003 CEA/DEN, EDF R& D +// +// +// +// File : Partition_Loop2d.cxx +// Author : Benedicte MARTIN +// Module : GEOM +// $Header$ + +//using namespace std; +#include "Partition_Loop2d.ixx" + +#include "utilities.h" +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +//======================================================================= +//function : Partition_Loop2d +//purpose : +//======================================================================= + +Partition_Loop2d::Partition_Loop2d() +{ +} + +//======================================================================= +//function : Init +//purpose : Init with the set of edges must have +// pcurves on . +//======================================================================= + +void Partition_Loop2d::Init(const TopoDS_Face& F) +{ + myConstEdges.Clear(); + myNewWires .Clear(); + myNewFaces .Clear(); + myFace = F; + myFaceOri = myFace.Orientation(); + myFace.Orientation( TopAbs_FORWARD ); +} + +//======================================================================= +//function : AddConstEdge +//purpose : Add as unique edge in the result. +//======================================================================= + +void Partition_Loop2d::AddConstEdge (const TopoDS_Edge& E) +{ +#ifdef DEB + Standard_Real f,l; + Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l); + if (pc.IsNull()) { + INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE"); + } else +#endif + { + myConstEdges.Append(E); + } +} + +void Partition_Loop2d::AddSectionEdge (const TopoDS_Edge& E) +{ +#ifdef DEB + Standard_Real f,l; + Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( E, myFace, f,l); + if (pc.IsNull()) + pc = BRep_Tool::CurveOnSurface( E, myFace, f,l); + gp_Vec2d Tg1; + gp_Pnt2d PC; + pc->D1(0.5*(f+l), PC, Tg1); + if (Tg1.Magnitude() <= gp::Resolution()) { + MESSAGE (""); + } + if (pc.IsNull()) { + INFOS( "AddConstEdge(): EDGE W/O PCURVE on FACE"); + } else +#endif + { + myConstEdges.Append(E); + myConstEdges.Append(E.Reversed()); + mySectionEdges.Add( E ); + } +} + +//======================================================================= +//function : preciseU +//purpose : find u such that the 3D point on theE is just out of tolerance +// of theV +//======================================================================= + +static Standard_Real preciseU (const BRepAdaptor_Surface& theSurf, + const TopoDS_Edge& theE, + const TopoDS_Vertex& theV, + const Handle(Geom2d_Curve)& theC, + const Standard_Boolean theFirstEnd) +{ + Standard_Boolean isForward = ( theE.Orientation () == TopAbs_FORWARD ); + if (theFirstEnd) isForward = !isForward; + + // find the first point in 2d and 3d + Standard_Real f,l; + BRep_Tool::Range( theE, f, l ); + Standard_Real u0 = isForward ? l : f; + gp_Pnt2d aP2d0 = theC->Value( u0 ); + gp_Pnt aPnt0 = theSurf.Value( aP2d0.X(), aP2d0.Y() ); + + // shift in 2d and 3d + Standard_Real du = ( l - f ) / 100, du3d = 0; + if (isForward) + du = -du; + + // target parameter + Standard_Real u; + + while (du3d < ::RealSmall()) + { + // u for test + u = u0 + du; + du *= 10; // for the next iteration: increase du untill du3d is large enough + + // find out how u is far from u0 in 3D + gp_Pnt2d aP2d = theC->Value( u ); + gp_Pnt aPnt = theSurf.Value( aP2d.X(), aP2d.Y() ); + du3d = aPnt0.Distance( aPnt ); + } + + // find u such that the 3D point is just out of tolerance of theV + Standard_Real tolV = BRep_Tool::Tolerance( theV ) + Precision::Confusion(); + u = u0 + du * tolV / du3d; + + // check that u is within the range + if ( isForward ? (u < f) : (u > l) ) + u = u0 + du; + + return u; +} + +//======================================================================= +//function : SelectEdge +//purpose : Find in the list the edge connected with by +// the vertex . +// is removed from the list. If is in +// with the same orientation, it's removed from the list +//======================================================================= + +static Standard_Boolean SelectEdge(const BRepAdaptor_Surface& Surf, + const TopoDS_Edge& CE, + const TopoDS_Vertex& CV, + TopoDS_Edge& NE, + const TopTools_ListOfShape& LE) +{ + NE.Nullify(); + + if (LE.Extent() > 1) { + //-------------------------------------------------------------- + // Several possible edges. + // - Test the edges differents of CE + //-------------------------------------------------------------- + TopoDS_Face FForward = Surf.Face(); + TopoDS_Edge aPrevNE; + + gp_Vec2d CTg1, Tg1, CTg2, Tg2; + gp_Pnt2d PC, P; + + Standard_Real f, l; + Handle(Geom2d_Curve) Cc, C; + Cc = BRep_Tool::CurveOnSurface(CE,FForward,f,l); + + Standard_Boolean isForward = ( CE.Orientation () == TopAbs_FORWARD ); + Standard_Real uc, u, du = Precision::PConfusion(); + uc = isForward ? ( l - du ) : ( f + du ); + Cc->D1(uc, PC, CTg1); + if (!isForward) CTg1.Reverse(); + + Standard_Real anglemin = 3 * PI, tolAng = 1.e-8; + + // select an edge whose first derivative is most left of CTg1 + // ie an angle between Tg1 and CTg1 is least + TopTools_ListIteratorOfListOfShape itl; + for ( itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + if (E.IsSame(CE)) + continue; + if (! CV.IsSame( TopExp::FirstVertex( E, Standard_True ))) + continue; + + isForward = ( E.Orientation () == TopAbs_FORWARD ); + + // get E curve + C = BRep_Tool::CurveOnSurface(E,FForward,f,l); + // get the first derivative Tg1 + u = isForward ? ( f + du ) : ( l - du ); + C->D1(u, P, Tg1); + if (!isForward) Tg1.Reverse(); + + // -PI < angle < PI + Standard_Real angle = Tg1.Angle(CTg1); + + if (PI - Abs(angle) <= tolAng) + { + // an angle is too close to PI; assure that an angle sign really + // reflects an edge position: +PI - an edge is worst, + // -PI - an edge is best. + u = preciseU( Surf, CE, CV, Cc, Standard_False); + gp_Vec2d CTg; + Cc->D1(u, PC, CTg); + if (CE.Orientation() == TopAbs_REVERSED) CTg.Reverse(); + + u = preciseU( Surf, E, CV, C, Standard_True); + C->D1(u, P, Tg1); + if (!isForward) Tg1.Reverse(); + + angle = Tg1.Angle(CTg); + } + + Standard_Boolean isClose = ( Abs( angle - anglemin ) <= tolAng ); + if (angle <= anglemin) { + if (isClose) + aPrevNE = NE; + else + aPrevNE.Nullify(); + anglemin = angle ; + NE = E; + } + else + if (isClose) + aPrevNE = E; + + } + if (!aPrevNE.IsNull()) { + // select one of close edges, the most left one. + Cc = BRep_Tool::CurveOnSurface( NE, FForward, f, l ); + uc = preciseU( Surf, NE, CV, Cc, Standard_True); + Cc->D1(uc, PC, CTg1); + if (NE.Orientation() != TopAbs_FORWARD) CTg1.Reverse(); + + u = preciseU( Surf, aPrevNE, CV, C, Standard_True); + C->D1(u, P, Tg1); + if (aPrevNE.Orientation() != TopAbs_FORWARD) Tg1.Reverse(); + + if ( Tg1.Angle(CTg1) < 0) + NE = aPrevNE; + } + } + else if (LE.Extent() == 1) { + NE = TopoDS::Edge(LE.First()); + } + else { + return Standard_False; + } + return !NE.IsNull(); +} + +//======================================================================= +//function : SamePnt2d +//purpose : +//======================================================================= + +static Standard_Boolean SamePnt2d(const TopoDS_Vertex& V1, + const TopoDS_Edge& E1, + const TopoDS_Vertex& V2, + const TopoDS_Edge& E2, + const TopoDS_Face& F) +{ + Standard_Real f1,f2,l1,l2; + Handle(Geom2d_Curve) C1 = BRep_Tool::CurveOnSurface(E1,F,f1,l1); + Handle(Geom2d_Curve) C2 = BRep_Tool::CurveOnSurface(E2,F,f2,l2); + + gp_Pnt2d P1 = C1->Value( BRep_Tool::Parameter(V1,E1)); + gp_Pnt2d P2 = C2->Value( BRep_Tool::Parameter(V2,E2)); + + Standard_Real Tol = 100 * BRep_Tool::Tolerance(V1); + Standard_Real Dist = P1.Distance(P2); + return Dist < Tol; +} + + +//======================================================================= +//function : StoreInMVE +//purpose : +//======================================================================= + +static void StoreInMVE (const TopoDS_Face& /*F*/, + TopoDS_Edge& E, + TopTools_DataMapOfShapeListOfShape& MVE ) + +{ + TopoDS_Vertex V1, V2; + TopTools_ListOfShape Empty; + + TopExp::Vertices(E,V1,V2); + if (!MVE.IsBound(V1)) { + MVE.Bind(V1,Empty); + } + MVE(V1).Append(E); + + if (!MVE.IsBound(V2)) { + MVE.Bind(V2,Empty); + } + MVE(V2).Append(E); +} + +//======================================================================= +//function : RemoveFromMVE +//purpose : +//======================================================================= + +static void RemoveFromMVE(const TopoDS_Edge& E, + TopTools_DataMapOfShapeListOfShape& MVE) +{ + TopTools_ListIteratorOfListOfShape itl; + TopoDS_Vertex V1,V2; + TopExp::Vertices (E,V1,V2); + if (MVE.IsBound(V1)) + for ( itl.Initialize(MVE(V1)); itl.More(); itl.Next()) { + if (itl.Value().IsEqual(E)) { + MVE(V1).Remove(itl); + break; + } + } + if (MVE.IsBound(V2)) + for ( itl.Initialize(MVE(V2)); itl.More(); itl.Next()) { + if (itl.Value().IsEqual(E)) { + MVE(V2).Remove(itl); + break; + } + } +} +//======================================================================= +//function : addConnected +//purpose : add to all edges reachable from +//======================================================================= + +static void addConnected(const TopoDS_Shape& E, + TopTools_MapOfShape& EM, + TopTools_MapOfShape& VM, + const TopTools_DataMapOfShapeListOfShape& MVE) +{ + // Loop on vertices of E + TopoDS_Iterator itV ( E ); + for ( ; itV.More(); itV.Next()) { + + if ( ! VM.Add ( itV.Value() )) continue; + + // Loop on edges sharing V + TopTools_ListIteratorOfListOfShape itE( MVE( itV.Value() ) ); + for (; itE.More(); itE.Next()) { + if ( EM.Add( itE.Value() )) + addConnected ( itE.Value(), EM, VM, MVE ); + } + } +} +//======================================================================= +//function : canPassToOld +//purpose : +//======================================================================= + +// static Standard_Boolean canPassToOld (const TopoDS_Shape& V, +// TopTools_MapOfShape& UsedShapesMap, +// const TopTools_DataMapOfShapeListOfShape& MVE, +// const TopTools_MapOfShape& SectionEdgesMap) +// { +// TopTools_ListIteratorOfListOfShape itE( MVE(V) ); +// // Loop on edges sharing V +// for (; itE.More(); itE.Next()) { +// if ( !UsedShapesMap.Add( itE.Value() )) +// continue; // already checked + +// if ( !SectionEdgesMap.Contains( itE.Value() )) +// return Standard_True; // WE PASSED + +// TopoDS_Iterator itV( itE.Value() ); +// // Loop on vertices of an edge +// for (; itV.More(); itV.Next()) { +// if ( !UsedShapesMap.Add( itV.Value() )) +// continue; // already checked +// else +// return canPassToOld( itV.Value(), UsedShapesMap, MVE, SectionEdgesMap); +// } +// } +// return Standard_False; +// } + +//======================================================================= +//function : MakeDegenAndSelect +//purpose : Find parameter of intersection of with and +// select an edge with its parameter closest to found one. +// Return new degenerated edge trimming by found parameters +//======================================================================= + +static TopoDS_Edge MakeDegenAndSelect(const TopoDS_Edge& CE, + const TopoDS_Vertex& CV, + TopoDS_Edge& NE, + TopTools_SequenceOfShape& EdgesSeq, + TColStd_SequenceOfReal& USeq, + const TopoDS_Edge& DE) +{ + if (EdgesSeq.Length() < 3) { + if (CE == EdgesSeq.First()) + NE = TopoDS::Edge( EdgesSeq.Last() ); + else + NE = TopoDS::Edge( EdgesSeq.First() ); + return DE; + } + + // find parameter on DE where it intersects CE + + Standard_Real U1; + Standard_Integer i, nb = EdgesSeq.Length(); + for (i=1; i<= nb; ++i) { + if (CE == EdgesSeq(i)) { + U1 = USeq(i); + break; + } + } + + // select NE with param closest to U1 thus finding U2 for a new degen edge + + Standard_Real U2, dU, dUmin = 1.e100; + Standard_Boolean isReversed = ( DE.Orientation() == TopAbs_REVERSED ); + for (i=1; i<= nb; ++i) { + dU = USeq(i) - U1; + if (isReversed ? (dU > 0) : (dU < 0)) + continue; + dU = Abs( dU ); + if ( dU > dUmin || IsEqual( dU, 0.)) + continue; + const TopoDS_Edge& E = TopoDS::Edge ( EdgesSeq(i) ); + if ( ! CV.IsSame( TopExp::FirstVertex( E , Standard_True ))) + continue; + NE = E; + dUmin = dU + Epsilon(dU); + U2 = USeq(i); + } + + // make a new degenerated edge + TopoDS_Edge NewDegen = TopoDS::Edge ( DE.EmptyCopied() ); + + Standard_Real Tol = BRep_Tool::Tolerance( CV ); + TopoDS_Vertex V = CV; + + BRep_Builder B; + V.Orientation( NewDegen.Orientation() ); + B.UpdateVertex( V, U1, NewDegen, Tol); + B.Add ( NewDegen , V ); + + V.Reverse(); + B.UpdateVertex( V, U2, NewDegen, Tol); + B.Add ( NewDegen , V ); + + return NewDegen; +} + +//======================================================================= +//function : prepareDegen +//purpose : Intersect with edges bound to its vertex in +// and store intersection parameter on in +// as well as the edges them-self in . +// Bind to vertex of in +//======================================================================= + +static void prepareDegen (const TopoDS_Edge& DegEdge, + const TopoDS_Face& F, + const TopTools_DataMapOfShapeListOfShape& MVE, + TopTools_SequenceOfShape& EdgesSeq, + TColStd_SequenceOfReal& USeq, + TopTools_DataMapOfShapeInteger& MVDEI, + const Standard_Integer DegEdgeIndex) +{ + const TopoDS_Vertex& V = TopExp::FirstVertex ( DegEdge ); + MVDEI.Bind ( V, DegEdgeIndex ); + + const TopTools_ListOfShape& EdgesList = MVE ( V ); + // if only 2 edges come to degenerated one, no pb in selection and + // no need to intersect them, just simulate asked data + Standard_Boolean doIntersect = ( EdgesList.Extent() > 2 ); + + BRepAdaptor_Curve2d DC, C; + Geom2dInt_GInter InterCC; + Standard_Real Tol = Precision::PConfusion(); + if ( doIntersect ) + DC.Initialize( DegEdge, F ); + + // avoid intersecting twice the same edge + BRepOffset_DataMapOfShapeReal EUMap ( EdgesList.Extent() ); + + Standard_Real U, f, l; + BRep_Tool::Range (DegEdge, f, l); + + TopTools_ListIteratorOfListOfShape itE (EdgesList); + for (; itE.More(); itE.Next()) { + + const TopoDS_Edge& E = TopoDS::Edge ( itE.Value() ); + + if ( !doIntersect) { + U = 0.; // it won't be used + } + else if ( BRep_Tool::IsClosed( E, F )) { + // seam edge: select U among f and l + Standard_Boolean first = Standard_True; + if ( V.IsSame ( TopExp::FirstVertex( E, Standard_True ) )) + first = Standard_False; + if ( DegEdge.Orientation() == TopAbs_REVERSED ) + first = !first; + U = first ? f : l; + } + else if ( EUMap.IsBound( E ) ) { + // same edge already bound + U = EUMap( E ); + } + else { + // intersect 2d curves + C.Initialize( E, F ); + InterCC.Perform ( DC, C , Tol, Tol ); + if (! InterCC.IsDone() || InterCC.NbPoints() == 0) { + MESSAGE ( "NO 2d INTERSECTION ON DEGENERATED EDGE" ); + continue; + } + // hope there is only one point of intersection + U = InterCC.Point( 1 ).ParamOnFirst(); + } + USeq.Append ( U ); + EdgesSeq.Append ( E ); + } +} +//======================================================================= +//function : Perform +//purpose : Make loops. +//======================================================================= + +void Partition_Loop2d::Perform() +{ + + Standard_Integer NbConstEdges = myConstEdges.Extent(); + TopTools_DataMapOfShapeListOfShape MVE(NbConstEdges) , MVE2(NbConstEdges); + TopTools_DataMapIteratorOfDataMapOfShapeListOfShape Mapit; + TopTools_ListIteratorOfListOfShape itl; + TopoDS_Vertex V1,V2; + BRepAdaptor_Surface Surface ( myFace, Standard_False ); + + // degenerated edges and parameters of their 2d intersection with other edges + TopoDS_Edge DE [2]; + TopTools_SequenceOfShape SEID [2]; // seq of edges intersecting degenerated + TColStd_SequenceOfReal SeqU [2]; // n-th U corresponds to n-th edge in SEID + TopTools_DataMapOfShapeInteger MVDEI(2); // map vertex - degenerated edge index + Standard_Integer iDeg = 0; // index of degenerated edge [0,1] + + //--------------------------------------------------------- + // Construction map vertex => edges, find degenerated edges + //--------------------------------------------------------- + for (itl.Initialize(myConstEdges); itl.More(); itl.Next()) { + TopoDS_Edge& E = TopoDS::Edge(itl.Value()); + if ( BRep_Tool::Degenerated( E )) { + if (DE[0].IsNull()) DE[0] = E; + else DE[1] = E; + } + else + StoreInMVE(myFace,E,MVE); + } + + // fill data for degenerated edges + if ( ! DE[0].IsNull() ) + prepareDegen ( DE[0], myFace, MVE, SEID[0], SeqU[0], MVDEI, 0); + if ( ! DE[1].IsNull() ) + prepareDegen ( DE[1], myFace, MVE, SEID[1], SeqU[1], MVDEI, 1); + + + // to detect internal wires + Standard_Boolean isInternCW = 0; + MVE2 = MVE; + + + //------------------------------ + // Construction of all the wires + //------------------------------ + // first, we collect wire edges in WEL list looking for same edges that + // will be then removed possibly exploding a wire into parts; + // second, build wire(s) + + while (!MVE.IsEmpty()) { + + TopoDS_Vertex VF,CV; + TopoDS_Edge CE,NE,EF; + TopoDS_Wire NW; + BRep_Builder B; + Standard_Boolean End = Standard_False; + TopTools_ListOfShape WEL; + + Mapit.Initialize(MVE); + if (Mapit.Value().IsEmpty()) { + MVE.UnBind(Mapit.Key()); + continue; + } + + // EF first edge. + EF = CE = TopoDS::Edge(Mapit.Value().First()); + // VF first vertex + VF = TopExp::FirstVertex( CE, Standard_True); + + isInternCW = Standard_True; + + TopTools_MapOfShape addedEM (NbConstEdges); // map of edges added to WEL + TopTools_MapOfShape doubleEM (NbConstEdges); // edges encountered twice in WEL + + //------------------------------- + // Construction of a wire. + //------------------------------- + while (!End) { + + // only a seam is allowed twice in a wire, the others should be removed + if (addedEM.Add ( CE ) || BRep_Tool::IsClosed( CE, myFace ) ) + WEL.Append( CE ); + else { + doubleEM.Add( CE ); + RemoveFromMVE (CE,MVE2); + TopoDS_Edge CERev = CE; + CERev.Reverse(); + RemoveFromMVE (CERev,MVE2); + } + + RemoveFromMVE (CE,MVE); + + CV = TopExp::LastVertex( CE, Standard_True); + + if (isInternCW && !mySectionEdges.Contains(CE)) + // wire is internal if all edges are section ones + isInternCW = Standard_False; + + if (MVDEI.IsBound( CV )) { // CE comes to the degeneration + iDeg = MVDEI( CV ); + TopoDS_Edge NewDegen; + NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]); + WEL.Append( NewDegen ); + CE = NE; + End = CV.IsSame( VF ); + continue; + } + + //-------------- + // stop test + //-------------- + if (MVE(CV).IsEmpty()) { + End=Standard_True; + MVE.UnBind(CV); + } + else if (CV.IsSame(VF) && SamePnt2d(CV,CE, VF,EF, myFace) ) { + End = Standard_True; + } + else { + //---------------------------- + // select new current edge + //---------------------------- + if (! SelectEdge (Surface,CE,CV,NE,MVE(CV))) { + MESSAGE ( " NOT CLOSED WIRE " ); + End=Standard_True; + } + else + CE = NE; + } + } // while ( !End ) + + + // WEL is built, built wire(s) + + + itl.Initialize( WEL ); + if ( doubleEM.IsEmpty()) { // no double edges + B.MakeWire( NW ); + for (; itl.More(); itl.Next()) + B.Add ( NW, itl.Value()); + if (isInternCW) myInternalWL.Append(NW); + else myNewWires.Append (NW); + } + + else { + // remove double and degenerated edges from WEL + while (itl.More()) { + const TopoDS_Edge& E = TopoDS::Edge ( itl.Value() ); + if ( doubleEM.Contains( E ) || BRep_Tool::Degenerated( E )) + WEL.Remove( itl ); + else + itl.Next(); + } + if ( WEL.IsEmpty()) + continue; + // remove double edges from SEID and SeqU + Standard_Integer i,j; + for (j=0; j<2; ++j) { + for (i=1; i<=SEID[j].Length(); ++i) { + if (doubleEM.Contains( SEID[j].Value(i))) { + SEID[j].Remove( i ); + SeqU[j].Remove( i-- ); + } + } + } + // removal of doulbe edges can explode a wire into parts, + // make new wires of them. + // A Loop like previous one but without 2d check + while ( !WEL.IsEmpty() ) { + CE = TopoDS::Edge( WEL.First() ); + WEL.RemoveFirst(); + B.MakeWire( NW ); + VF = TopExp::FirstVertex ( CE, Standard_True); + + End = Standard_False; + while ( !End) { + B.Add( NW, CE ); + CV = TopExp::LastVertex ( CE, Standard_True); + + if (MVDEI.IsBound( CV )) { // CE comes to the degeneration + iDeg = MVDEI( CV ); + TopoDS_Edge NewDegen; + NewDegen = MakeDegenAndSelect( CE, CV, NE, SEID[iDeg], SeqU[iDeg], DE[iDeg]); + B.Add( NW, NewDegen ); + End = CV.IsSame( VF ); + CE = NE; + if (!NE.IsNull()) { // remove NE from WEL + for (itl.Initialize( WEL ); itl.More(); itl.Next()) + if ( NE == itl.Value()) { + WEL.Remove( itl ); + break; + } + } + } // end degeneration + + else { + if (CV.IsSame( VF )) { + End = Standard_True; + continue; + } + // edges in WEL most often are well ordered + // so try to iterate until the End + Standard_Boolean add = Standard_False; + itl.Initialize(WEL); + while ( itl.More() && !End) { + NE = TopoDS::Edge( itl.Value() ); + if ( CV.IsSame( TopExp::FirstVertex( NE, Standard_True ))) { + WEL.Remove( itl ); + if (add) + B.Add( NW, CE ); + CE = NE; + add = Standard_True; + CV = TopExp::LastVertex( CE, Standard_True); + if (MVDEI.IsBound( CV ) || CV.IsSame( VF )) + break; + } + else + itl.Next(); + } + if (!add) + End = Standard_True; + } + } // !End + + myInternalWL.Append( NW ); + } + } // end building new wire(s) from WEL + + } // end Loop on MVE + + // all wires are built + + + // ============================================================ + // select really internal wires i.e. those from which we can`t + // pass to an old (not section) edge + // ============================================================ + + Standard_Integer nbIW = myInternalWL.Extent(); + if (nbIW == 0) + return; + + if ( myNewWires.Extent() != 1 && nbIW > 1) { + TopTools_MapOfShape outerEM (NbConstEdges); // edges connected to non-section ones + TopTools_MapOfShape visitedVM (NbConstEdges); + for ( itl.Initialize( myConstEdges ); itl.More(); itl.Next()) { + if ( ! mySectionEdges.Contains( itl.Value() )) + addConnected (itl.Value(), outerEM, visitedVM, MVE2); + } + // if an edge of a wire is in , the wire is not internal + TopExp_Explorer expIWE; + TopTools_ListIteratorOfListOfShape itIW ( myInternalWL ); + while (itIW.More()) { + expIWE.Init ( itIW.Value() , TopAbs_EDGE ); + if ( outerEM.Contains( expIWE.Current() )) { + myNewWires.Append ( itIW.Value() ); + myInternalWL.Remove( itIW ); // == itIW.Next() + } + else + itIW.Next(); + } + } +} +//======================================================================= +//function : isHole +//purpose : +//======================================================================= + +static Standard_Boolean isHole (const TopoDS_Wire& W, + const TopoDS_Face& F) +{ + BRep_Builder B; + TopoDS_Shape newFace = F.EmptyCopied(); + B.Add(newFace,W.Oriented(TopAbs_FORWARD)); + BRepTopAdaptor_FClass2d classif (TopoDS::Face(newFace), + Precision::PConfusion()); + return (classif.PerformInfinitePoint() == TopAbs_IN); +} + +//======================================================================= +//function : IsInside +//purpose : check if W1 is inside W2. Suppose W2 is not a hole !!!! +//======================================================================= + +static Standard_Boolean isInside(const TopoDS_Face& F, + const TopoDS_Wire& W1, + const TopoDS_Wire& W2) +{ + // make a face with wire W2 + BRep_Builder B; + TopoDS_Shape aLocalShape = F.EmptyCopied(); + TopoDS_Face newFace = TopoDS::Face(aLocalShape); + B.Add(newFace,W2); + + // get any 2d point of W1 + TopExp_Explorer exp(W1,TopAbs_EDGE); + if (BRep_Tool::Degenerated( TopoDS::Edge( exp.Current() ))) + exp.Next(); + const TopoDS_Edge& e = TopoDS::Edge(exp.Current()); + Standard_Real f,l; + Handle(Geom2d_Curve) C2d = BRep_Tool::CurveOnSurface(e,F,f,l); + gp_Pnt2d pt2d(C2d->Value( 0.5 * ( f + l ))); + + BRepTopAdaptor_FClass2d classif(newFace,Precision::PConfusion()); + return (classif.Perform(pt2d) == TopAbs_IN); +} + +//======================================================================= +//function : NewWires +//purpose : Returns the list of wires performed. +// can be an empty list. +//======================================================================= + +const TopTools_ListOfShape& Partition_Loop2d::NewWires() const +{ + return myNewWires; +} + +//======================================================================= +//function : NewFaces +//purpose : Returns the list of faces. +//Warning : The method as to be called before. +// can be an empty list. +//======================================================================= + +const TopTools_ListOfShape& Partition_Loop2d::NewFaces() const +{ + return myNewFaces; +} + +//======================================================================= +//function : findEqual +//purpose : move wires form to pairs of wires build of the same edges +//======================================================================= + +static void findEqual (TopTools_ListOfShape& WL, + TopTools_DataMapOfShapeShape& EqWM, + const TopoDS_Face& F) +{ + TopTools_ListIteratorOfListOfShape it1, it2; + Standard_Integer i,j; + TColStd_MapOfInteger IndMap; + for (it1.Initialize(WL), i=1; it1.More(); it1.Next(), i++) { + + if (IndMap.Contains(i)) continue; + const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value()); + + for (it2.Initialize(WL), j=1; it2.More(); it2.Next(), j++) { + + if (j <= i || IndMap.Contains(j)) continue; + + TopTools_IndexedMapOfShape EdgesMap; + TopExp::MapShapes (Wire1, TopAbs_EDGE, EdgesMap); + + const TopoDS_Shape& Wire2 = it2.Value(); + TopoDS_Iterator itE ( Wire2); + for (; itE.More(); itE.Next()) { + if ( !EdgesMap.Contains( itE.Value()) ) + break; + } + if (!itE.More()) { // all edges are same + if (isHole( Wire1, F)) { + EqWM.Bind ( Wire1, Wire2 ); + } + else { + EqWM.Bind ( Wire2, Wire1 ); + } + IndMap.Add(i); + IndMap.Add(j); + break; + } + } + } + // clear WL + it1.Initialize(WL); + i=1; + while (it1.More()) { + if (IndMap.Contains(i)) + WL.Remove(it1); // next node becomes current and with Next() we would miss it + else + it1.Next(); + i++; + } +} + +//======================================================================= +//function : classify +//purpose : bind to a wire a list of internal wires +//======================================================================= + +static void classify(const TopTools_DataMapOfShapeShape& EqWM, + BRepAlgo_AsDes& OuterInner, + const TopoDS_Face& F) +{ + TopTools_DataMapIteratorOfDataMapOfShapeShape it1, it2; + + for (it1.Initialize(EqWM); it1.More(); it1.Next()) { + // find next after it1.Value() + for (it2.Initialize(EqWM); it2.More(); it2.Next()) + if (it1.Value().IsSame( it2.Value() )) + { + it2.Next(); + break; + } + for ( ; it2.More(); it2.Next()) { + const TopoDS_Wire& Wire1 = TopoDS::Wire( it1.Value() ); + const TopoDS_Wire& Wire2 = TopoDS::Wire( it2.Value() ); + if (isInside(F, Wire1, Wire2)) + OuterInner.Add (Wire2, Wire1); + else if (isInside(F, Wire2, Wire1)) + OuterInner.Add (Wire1, Wire2); + } + } +} +//======================================================================= +//function : WiresToFaces +//purpose : Build faces from the wires result. +// serves to find original edge by new +// one.
contains edges resulting from face +// intersections +//======================================================================= + +void Partition_Loop2d::WiresToFaces(const BRepAlgo_Image& ) +{ + Standard_Integer nbW = myNewWires.Extent() + myInternalWL.Extent(); + if (nbW==0) + return; + + BRepAlgo_FaceRestrictor FR; + FR.Init (myFace,Standard_False); + + // FaceRestrictor is instable in rather simple cases + // (ex. a single face of bellecoque.brep splited by 10 planes: + // sometimes 1-2 faces are missing ). + // So we use it as less as possible: no holes -> make faces by hands + + + // are there holes in myFace ? + Standard_Boolean hasOldHoles = Standard_False; + TopoDS_Iterator itOldW (myFace); + if ( itOldW.More()) { + const TopoDS_Wire& FirstOldWire = TopoDS::Wire( itOldW.Value() ); + itOldW.Next(); + hasOldHoles = itOldW.More() || isHole( FirstOldWire, myFace); + } + if (myInternalWL.IsEmpty() && !hasOldHoles) { + // each wire bounds one face + BRep_Builder B; + TopTools_ListIteratorOfListOfShape itNW (myNewWires); + for (; itNW.More(); itNW.Next()) { + TopoDS_Face NF = TopoDS::Face ( myFace.EmptyCopied() ); + B.Add ( NF, itNW.Value() ); + NF.Orientation( myFaceOri); + myNewFaces.Append ( NF ); + } + return; + } + + // FaceRestrictor can't classify wires build on all the same edges + // and gives incorrect result in such cases (ex. a plane cut into 2 parts by cylinder) + // We must make faces of equal wires separately. One of equal wires makes a + // hole in a face and should come together with outer wires of face. + // The other of a wires pair bounds a face that may have holes in turn. + + // Find equal wires among internal wires + TopTools_DataMapOfShapeShape EqWM; // key is a hole part of a pair of equal wires + findEqual (myInternalWL, EqWM, myFace); + + if (!EqWM.IsEmpty()) { // there are equal wires + + if (hasOldHoles) + myInternalWL.Append( myNewWires ); // an old wire can be inside an equal wire + + // classify equal wire pairs + BRepAlgo_AsDes OuterInner; + classify (EqWM,OuterInner,myFace); + + // make face of most internal of equal wires and its inner wires + while ( !EqWM.IsEmpty()) { + + TopTools_ListOfShape prevHolesL; // list of hole-part of previous most internal equal wires + + // find most internal wires among pairs (key - hole, value - outer part) + TopTools_DataMapIteratorOfDataMapOfShapeShape it(EqWM); + Standard_Integer nbEqW = EqWM.Extent(); // protection against infinite loop + for ( ; it.More(); it.Next()) { + + TopoDS_Wire outerW = TopoDS::Wire ( it.Value() ); + if ( OuterInner.HasDescendant( outerW ) && // has internal + ! OuterInner.Descendant( outerW ).IsEmpty() ) + continue; + + FR.Add( outerW ); + + // add internal wires that are inside of outerW + TopTools_ListIteratorOfListOfShape itIW (myInternalWL); + while ( itIW.More()) { + TopoDS_Wire IW = TopoDS::Wire ( itIW.Value() ); + if ( isInside (myFace, IW, outerW)) { + FR.Add (IW); + myInternalWL.Remove( itIW ); // == itIW.Next() !!! + } + else + itIW.Next(); + } + + // the hole-part of current pair of equal wires will be in the next new face + prevHolesL.Append ( it.Key() ); + + } // Loop on map of equal pairs searching for innermost wires + + // make faces + FR.Perform(); + if (FR.IsDone()) { + for (; FR.More(); FR.Next()) + myNewFaces.Append(FR.Current()); + } + + FR.Clear(); + + // add hole-parts to FaceRestrictor, + // remove them from the EqWM, + // remove found wires as internal of resting classified wires + Standard_Boolean clearOuterInner = ( prevHolesL.Extent() < EqWM.Extent() ); + TopTools_ListIteratorOfListOfShape itPrev (prevHolesL); + for (; itPrev.More(); itPrev.Next()) { + TopoDS_Wire& Hole = TopoDS::Wire ( itPrev.Value() ); + FR.Add ( Hole ); + if (clearOuterInner) { + const TopoDS_Wire& outerW = TopoDS::Wire ( EqWM.Find( Hole ) ); + // Loop on wires including outerW + TopTools_ListIteratorOfListOfShape itO( OuterInner.Ascendant( outerW )); + for (; itO.More(); itO.Next()) { + TopTools_ListOfShape& innerL = OuterInner.ChangeDescendant( itO.Value() ); + TopTools_ListIteratorOfListOfShape itI (innerL); + // Loop on internal wires of current including wire + for (; itI.More(); itI.Next()) + if ( outerW.IsSame( itI.Value() )) { + innerL.Remove( itI ); break; + } + } + } + EqWM.UnBind ( Hole ); + } + + if (nbEqW == EqWM.Extent()) + { + // error: pb with wires classification +#ifdef DEB + MESSAGE("Partition_Loop2d::WiresToFaces(), pb with wires classification"); +#endif + break; + } + + } // while (!EqWM.IsEmpty) + + } // if !EqWM.IsEmpty() + + myNewWires.Append ( myInternalWL ); + + TopTools_ListIteratorOfListOfShape itW (myNewWires); + for (; itW.More(); itW.Next()) { + TopoDS_Wire& W = TopoDS::Wire ( itW.Value() ); + FR.Add(W); + } + FR.Perform(); + for (; FR.IsDone() && FR.More(); FR.Next()) + myNewFaces.Append(FR.Current()); + + + TopTools_ListIteratorOfListOfShape itNF (myNewFaces); + for (; itNF.More(); itNF.Next()) + itNF.Value().Orientation( myFaceOri ); +} + +#endif diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Loop3d.cxx netgen-4.5_new/libsrc/occ/Partition_Loop3d.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Loop3d.cxx 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Loop3d.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -10,6 +10,11 @@ // Module : GEOM //using namespace std; + + + +#include "utilities.h" + #include "Partition_Loop3d.ixx" #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Spliter.cxx netgen-4.5_new/libsrc/occ/Partition_Spliter.cxx --- netgen-4.5_orig/libsrc/occ/Partition_Spliter.cxx 2005-07-11 10:33:27.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/Partition_Spliter.cxx 2010-06-23 13:19:48.000000000 +0400 @@ -29,14 +29,15 @@ // $Header$ //using namespace std; + +#include "utilities.h" + #include "Partition_Inter2d.hxx" #include "Partition_Inter3d.hxx" #include "Partition_Loop2d.hxx" #include "Partition_Loop3d.hxx" #include "Partition_Spliter.ixx" -#include "utilities.h" - #include #include #include diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/Partition_Spliter.cxx.orig netgen-4.5_new/libsrc/occ/Partition_Spliter.cxx.orig --- netgen-4.5_orig/libsrc/occ/Partition_Spliter.cxx.orig 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/Partition_Spliter.cxx.orig 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,2166 @@ +#ifdef OCCGEOMETRY + +// GEOM PARTITION : partition algorithm +// +// Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN, +// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS +// +// 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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org +// +// +// +// File : Partition_Spliter.cxx +// Author : Benedicte MARTIN +// Module : GEOM +// $Header$ + +//using namespace std; +#include "Partition_Inter2d.hxx" +#include "Partition_Inter3d.hxx" +#include "Partition_Loop2d.hxx" +#include "Partition_Loop3d.hxx" +#include "Partition_Spliter.ixx" + +#include "utilities.h" + +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#ifdef DEB +//# define PART_PERF +#endif + +#ifdef PART_PERF +# include +#endif + +//======================================================================= +//function : isClosed +//purpose : check id a shape is closed, ie is a solid or a closed shell +//======================================================================= + +static Standard_Boolean isClosed(const TopoDS_Shape& theShape) +{ + Standard_Boolean isClosed = (theShape.ShapeType() == TopAbs_SOLID); + + if (!isClosed && theShape.ShapeType() == TopAbs_SHELL) { + TopTools_IndexedDataMapOfShapeListOfShape MEF; + TopExp::MapShapesAndAncestors(theShape, TopAbs_EDGE, TopAbs_FACE, MEF); + for (Standard_Integer i=1; isClosed && i<=MEF.Extent(); ++i) + isClosed = ( MEF(i).Extent() != 1 ); + } + + return isClosed; +} + +//======================================================================= +//function : Partition_Spliter +//purpose : constructor +//======================================================================= + +Partition_Spliter::Partition_Spliter() +{ + myAsDes = new BRepAlgo_AsDes; + Clear(); +} + +//======================================================================= +//function : AddTool +//purpose : add cutting tool that will _NOT_ be in result +//======================================================================= + +void Partition_Spliter::AddTool(const TopoDS_Shape& S) +{ + if (S.ShapeType() < TopAbs_SOLID) { // compound or compsolid + TopoDS_Iterator it (S); + for (; it.More(); it.Next()) + { + AddTool( it.Value()); + myFaceShapeMap.Bind( it.Value(), S ); // to know compound by shape + } + return; + } + + for (TopExp_Explorer exp(S,TopAbs_FACE); exp.More(); exp.Next()) + { + myMapTools.Add(exp.Current()); + myFaceShapeMap.Bind( exp.Current(), S ); + } + if (isClosed( S )) + myClosedShapes.Add( S ); +} + +//======================================================================= +//function : AddShape +//purpose : add object Shape to be splited +//======================================================================= + +void Partition_Spliter::AddShape(const TopoDS_Shape& S) +{ + if (S.ShapeType() < TopAbs_SOLID) { // compound or compsolid + TopoDS_Iterator it (S); + for (; it.More(); it.Next()) + { + AddShape( it.Value()); + myFaceShapeMap.Bind( it.Value(), S ); // to know compound by shape + } + return; + } + + TopExp_Explorer exp(S,TopAbs_FACE); + if (!exp.More()) { // do not split edges and vertices + //myBuilder.Add( myShape, S ); + return; + } + + Standard_Integer nbFacesBefore = myMapFaces.Extent(); // not to add twice the same S + for (; exp.More(); exp.Next()) { + const TopoDS_Shape & aFace = exp.Current(); + if ( ! myFaceShapeMap.IsBound( aFace )) // keep shape of tool face added as object + myFaceShapeMap.Bind( aFace, S ); + if (myMapFaces.Add( aFace )) + myImagesFaces.SetRoot( aFace ); + } + + if (nbFacesBefore == myMapFaces.Extent()) + return; + + // solids must be processed before all + if (S.ShapeType() == TopAbs_SOLID) + myListShapes.Prepend(S); + else + myListShapes.Append(S); + + if (isClosed( S )) + myClosedShapes.Add( S ); + +} + +//======================================================================= +//function : Shape +//purpose : return resulting compound +//======================================================================= + +TopoDS_Shape Partition_Spliter::Shape() const +{ + return myShape; +} + +//======================================================================= +//function : Clear +//purpose : clear fields +//======================================================================= + +void Partition_Spliter::Clear() +{ + myDoneStep = TopAbs_SHAPE; + + myListShapes.Clear(); + myMapFaces.Clear(); + myMapTools.Clear(); + myEqualEdges.Clear(); + myNewSection.Clear(); + myClosedShapes.Clear(); + mySharedFaces.Clear(); + myWrappingSolid.Clear(); + myFaceShapeMap.Clear(); + + myInternalFaces.Clear(); + myIntNotClFaces.Clear(); + + myAsDes->Clear(); + myImagesFaces.Clear(); + myImagesEdges.Clear(); + myImageShape.Clear(); + + // myInter3d = Partition_Inter3d(myAsDes); + Partition_Inter3d hinter3d (myAsDes); + myInter3d = hinter3d; + + myAddedFacesMap.Clear(); + +} + +//======================================================================= +//function : Compute +//purpose : produce a result +//======================================================================= + +void Partition_Spliter::Compute(const TopAbs_ShapeEnum Limit) +{ + if ((Limit != TopAbs_SHAPE && myDoneStep == Limit) || + (Limit == TopAbs_SHAPE && myDoneStep == TopAbs_SOLID)) + return; + + myBuilder.MakeCompound( myShape ); + + TopTools_MapIteratorOfMapOfShape it; + TopTools_ListIteratorOfListOfShape itl; + TopExp_Explorer exp; + +#ifdef PART_PERF + OSD_Chronometer aCron; +#endif + + if (myDoneStep > TopAbs_VERTEX) { + + TopTools_ListOfShape aListFaces; + aListFaces = myImagesFaces.Roots(); + for (it.Initialize(myMapTools); it.More(); it.Next()) + aListFaces.Append(it.Key()); + +#ifdef PART_PERF + aCron.Start(); +#endif + + //----------------------------------------------- + // Intersection between faces + //----------------------------------------------- + // result is in myAsDes as a map Face - list of new edges; + // special care is done for section edges, same domain faces and vertices: + // data about them is inside myInter3d + myInter3d.CompletPart3d(aListFaces, myFaceShapeMap); + +#ifdef PART_PERF + MESSAGE("+++ CompletPart3d()"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + //----------------------------------------------- + // Intersection of edges + //----------------------------------------------- + + // add tool faces which must be reconstructed to myMapFaces too + FindToolsToReconstruct(); + +#ifdef PART_PERF + MESSAGE("+++ FindToolsToReconstruct()"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + + // add existing vertices to edges of object faces in myAsDes + TopTools_MapOfShape DoneEM; + for ( it.Initialize(myMapFaces); it.More(); it.Next()) { + const TopoDS_Shape& F = it.Key(); + TopoDS_Face FForward = TopoDS::Face(F.Oriented(TopAbs_FORWARD)); + for (exp.Init(FForward,TopAbs_EDGE); exp.More(); exp.Next()) { + const TopoDS_Edge& E = TopoDS::Edge( exp.Current() ); + myAsDes->Add(FForward,E); + if (DoneEM.Add(E)) { + TopoDS_Iterator itV(E); + for (; itV.More(); itV.Next()) { + const TopoDS_Vertex& V = TopoDS::Vertex( itV.Value()); + myAsDes->Add(E, myInter3d.ReplaceSameDomainV( V, E )); + } + } + } + } + + // intersect edges that are descendants of a face in myAsDes + TopTools_MapOfShape& Modif = myInter3d.TouchedFaces(); + for ( it.Initialize(Modif); it.More(); it.Next()) { + const TopoDS_Face& F = TopoDS::Face(it.Key()); + Partition_Inter2d::CompletPart2d (myAsDes, F, myInter3d.NewEdges()); + } + // now myAsDes contains also new vertices made at edge intersection as + // descendant of edges both new and old + + myDoneStep = TopAbs_VERTEX; + +#ifdef PART_PERF + MESSAGE("+++ CompletPart2d()"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + } // if (myDoneStep > TopAbs_VERTEX) + + if (Limit == TopAbs_VERTEX) { + // add new vertices to myShape + for ( it.Initialize( myInter3d.NewEdges() ); it.More(); it.Next()) { + if (! myAsDes->HasDescendant( it.Key() )) + continue; + itl.Initialize( myAsDes->Descendant( it.Key() )); + for (; itl.More(); itl.Next()) + myBuilder.Add ( myShape, itl.Value() ); + } + return; + } + + + if (myDoneStep > TopAbs_EDGE) { + + //----------------------------------------------- + //----------------------------------------------- + // ------- Reconstruction of all the edges.------ + //----------------------------------------------- + //----------------------------------------------- + + // ============== + // cut new edges + // ============== + TopTools_ListOfShape LSE; // all edge splits + for ( it.Initialize(myInter3d.NewEdges()); it.More(); it.Next()) { + + TopoDS_Vertex V1,V2; + TopoDS_Edge EE = TopoDS::Edge(it.Key()); + + TopTools_ListOfShape aListV, aListF; + aListV = myAsDes->Descendant(EE); // intersection vertices + aListF = myAsDes->Ascendant(EE); // intersected faces + + if (aListV.IsEmpty()) + continue; // new edge does not intersect any other edge + + // Add end vertices to new edges only if + // one face is Tool and the other is Shape + Standard_Boolean isTool1 = ! myMapFaces.Contains( aListF.First() ); + Standard_Boolean isTool2 = ! myMapFaces.Contains( aListF.Last() ); + if (isTool1 || isTool2) + { + TopExp::Vertices(EE,V1,V2); + Standard_Real Tol = Max (BRep_Tool::Tolerance( V1 ), + BRep_Tool::Tolerance( V2 )); + + gp_Pnt P1 = BRep_Tool::Pnt(V1); + gp_Pnt P2 = BRep_Tool::Pnt(V2); + Standard_Boolean AddV1 = Standard_True; + Standard_Boolean AddV2 = Standard_True; + + // add only if there is no intersection at end vertex + for (itl.Initialize(aListV); itl.More(); itl.Next()) { + const TopoDS_Vertex& Ve = TopoDS::Vertex(itl.Value()) ; + Standard_Real Tol2 = Max ( Tol, BRep_Tool::Tolerance( Ve )); + Tol2 *= Tol2; + gp_Pnt P = BRep_Tool::Pnt(Ve); + if (AddV1 && P.SquareDistance(P1) <= Tol2) + AddV1 = Standard_False; + + if (AddV2 && P.SquareDistance(P2) <= Tol2) + AddV2 = Standard_False; + } + + if (AddV1) { + aListV.Append(V1); + myAsDes->Add(EE,V1); + } + + if (AddV2) { + aListV.Append(V2); + myAsDes->Add(EE,V2); + } + } + + // cut new edges + Standard_Integer NbV=aListV.Extent() ; + if (NbV>1 || (NbV==1 && V1.IsSame(V2)) ) { + TopTools_ListOfShape LNE; + MakeEdges (EE,aListV, LNE); + myImagesEdges.Bind(EE,LNE); + LSE.Append( LNE ); + } + } + + // ============== + // cut old edges + // ============== + for ( it.Initialize(myMapFaces); it.More(); it.Next()) { + for (exp.Init( it.Key(), TopAbs_EDGE); exp.More(); exp.Next()) { + const TopoDS_Edge& EE = TopoDS::Edge( exp.Current() ); + if ( myImagesEdges.HasImage( EE )) + continue; + TopTools_ListOfShape LNE; + const TopTools_ListOfShape& aListVV = myAsDes->Descendant(EE); + MakeEdges (EE, aListVV, LNE); + myImagesEdges.Bind(EE,LNE); + LSE.Append( LNE ); + } + } +#ifdef PART_PERF + MESSAGE("+++ Cut Edges"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + + // process same domain section edges + MergeEqualEdges( LSE ); + + myDoneStep = TopAbs_EDGE; + +#ifdef PART_PERF + MESSAGE("+++ MergeEqualEdges()"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + } // if (myDoneStep > TopAbs_EDGE) + + if (Limit == TopAbs_EDGE) { + // add splits of old edges + TopTools_ListIteratorOfListOfShape itNE; + for (itl.Initialize( myListShapes );itl.More();itl.Next()) { + if (myMapTools.Contains( itl.Value() )) + continue; // skip tool faces + for ( exp.Init( itl.Value(), TopAbs_EDGE ); exp.More(); exp.Next()) { + itNE.Initialize( myImagesEdges.Image( exp.Current() )); + for ( ; itNE.More(); itNE.Next()) + myBuilder.Add ( myShape, itNE.Value() ); + } + } + // add splits of new edges + for ( it.Initialize( myInter3d.NewEdges() ); it.More(); it.Next()) { + itNE.Initialize( myImagesEdges.Image( it.Key() )); + for (; itNE.More(); itNE.Next()) + myBuilder.Add ( myShape, itNE.Value() ); + } + return; + } + + + //----------------------------------------------- + // split faces + //----------------------------------------------- + + if (myDoneStep > TopAbs_FACE) { + + for (itl.Initialize(myListShapes);itl.More();itl.Next()) { + TopoDS_Shape FacesComp = MakeFaces ( itl.Value()); + // there is a cunning here: myImagesFaces keeps faces made by Loop2d + // but some of them may be replaced with splits of same domain face + // and myImageShape keeps ultimate result + myImageShape.Bind( itl.Value(), FacesComp ); + } + + myDoneStep = TopAbs_FACE; +#ifdef PART_PERF + MESSAGE("+++ MakeFaces()"); + aCron.Show( cout ); + aCron.Reset(); + aCron.Start(); +#endif + } + + if (Limit == TopAbs_WIRE || + Limit == TopAbs_FACE) { + for (itl.Initialize(myListShapes);itl.More();itl.Next()) { + if ( myMapTools.Contains( itl.Value() )) + continue; // no result needed for a tool face + const TopoDS_Shape& FacesComp = myImageShape.Image( itl.Value() ).First(); + for ( exp.Init( FacesComp, Limit); exp.More(); exp.Next()) + myBuilder.Add ( myShape, exp.Current()); + } + return; + } + + + //----------------------------------------------- + // split and add solids and shells + //----------------------------------------------- + + Standard_Boolean makeSolids = (Limit == TopAbs_SHAPE || + Limit < TopAbs_SHELL); + for (itl.Initialize(myListShapes);itl.More();itl.Next()) + { + const TopoDS_Shape & S = itl.Value(); + if (S.ShapeType() > TopAbs_SHELL) + continue; + + TopTools_ListOfShape NSL; // new shape list + MakeShells (S , NSL); + if (makeSolids && S.ShapeType() == TopAbs_SOLID ) + MakeSolids( S, NSL ); + + // store new shells or solids + TopTools_ListIteratorOfListOfShape itNSL (NSL); + for ( ; itNSL.More(); itNSL.Next()) + myBuilder.Add (myShape, itNSL.Value()); + } +#ifdef PART_PERF + MESSAGE("+++ MakeShells()"); + aCron.Show( cout ); +#endif + + //----------------------------------------------- + // add split faces + //----------------------------------------------- + + for (itl.Initialize(myListShapes);itl.More();itl.Next()) + { + const TopoDS_Shape & S = itl.Value(); + if (S.ShapeType() != TopAbs_FACE || + myMapTools.Contains( S )) + continue; + TopoDS_Iterator itS( myImageShape.Image(S).First() ); + for (; itS.More(); itS.Next()) + if (! myAddedFacesMap.Contains( itS.Value() )) + myBuilder.Add (myShape, itS.Value()); + } + + myDoneStep = makeSolids ? TopAbs_SOLID : TopAbs_SHELL; + +} + +//======================================================================= +//function : MakeSolids +//purpose : make solids out of Shells +//======================================================================= + +void Partition_Spliter::MakeSolids(const TopoDS_Shape & theSolid, + TopTools_ListOfShape & theShellList) +{ + // for a solid wrapping other shells or solids without intersection, + // it is necessary to find shells making holes in it + + TopTools_ListOfShape aNewSolids; // result + TopTools_ListOfShape aHoleShells; + TopoDS_Shape anInfinitePointShape; + + Standard_Boolean isWrapping = myWrappingSolid.Contains( theSolid ); + if (!isWrapping && !theShellList.IsEmpty()) + { + // check if theSolid initially has internal shells + TopoDS_Iterator aShellExp (theSolid); + aShellExp.Next(); + isWrapping = aShellExp.More(); + } + + TopTools_ListIteratorOfListOfShape aShellIt(theShellList); + for ( ; aShellIt.More(); aShellIt.Next()) + { + const TopoDS_Shape & aShell = aShellIt.Value(); + + // check if a shell is a hole + if (isWrapping && IsInside (anInfinitePointShape, aShell)) + aHoleShells.Append( aShell ); + else + { + // make a solid from a shell + TopoDS_Solid Solid; + myBuilder.MakeSolid( Solid ); + myBuilder.Add (Solid, aShell); + + aNewSolids.Append (Solid); + } + } + + // find an outer shell most close to each hole shell + TopTools_DataMapOfShapeShape aInOutMap; + for (aShellIt.Initialize( aHoleShells ); aShellIt.More(); aShellIt.Next()) + { + const TopoDS_Shape & aHole = aShellIt.Value(); + TopTools_ListIteratorOfListOfShape aSolisIt (aNewSolids); + for ( ; aSolisIt.More(); aSolisIt.Next()) + { + const TopoDS_Shape & aSolid = aSolisIt.Value(); + if (! IsInside( aHole, aSolid )) + continue; + + if ( aInOutMap.IsBound (aHole)) + { + const TopoDS_Shape & aSolid2 = aInOutMap( aHole ); + if ( IsInside( aSolid, aSolid2 )) + { + aInOutMap.UnBind( aHole ); + aInOutMap.Bind ( aHole, aSolid ); + } + } + else + aInOutMap.Bind ( aHole, aSolid ); + } + + // add aHole to a solid + if (aInOutMap.IsBound( aHole )) + myBuilder.Add ( aInOutMap( aHole ), aHole ); + } + + theShellList.Clear(); + theShellList.Append( aNewSolids ); +} + +//======================================================================= +//function : FindFacesInside +//purpose : return compound of faces of other shapes that are +// inside . +// is an object shape. +// makes avoid faces that do not form a +// closed shell +// makes return already added faces +//======================================================================= + +TopoDS_Shape Partition_Spliter::FindFacesInside(const TopoDS_Shape& theShape, + const Standard_Boolean CheckClosed, + const Standard_Boolean All) +{ + // ================================================ + // check if internal faces have been already found + // ================================================ + TopExp_Explorer expl; + if (myInternalFaces.IsBound( theShape )) + { + TopoDS_Shape aIntFComp = myInternalFaces.Find ( theShape ); + TopoDS_Shape aIntRemFComp = myIntNotClFaces.Find ( theShape ); + + expl.Init( aIntRemFComp, TopAbs_FACE); + if (CheckClosed || !expl.More()) + return aIntFComp; + + TopoDS_Compound C; + myBuilder.MakeCompound( C ); + // add removed faces + for (; expl.More(); expl.Next()) + myBuilder.Add( C, expl.Current() ); + // add good internal faces + for (expl.Init( aIntFComp, TopAbs_FACE); expl.More(); expl.Next()) + myBuilder.Add( C, expl.Current() ); + return C; + } + + // =================================== + // get data for internal faces search + // =================================== + + // compound of split faces of theShape + const TopoDS_Shape& CSF = myImageShape.Image(theShape).First(); + + TopTools_MapOfShape MSE, MFP; + TopTools_DataMapOfShapeListOfShape DMSEFP; + TopTools_MapIteratorOfMapOfShape itm; + TopTools_ListOfShape EmptyL; + + // MSE filling: map of new section edges of CSF + for (expl.Init(CSF,TopAbs_EDGE); expl.More(); expl.Next()) { + const TopoDS_Shape & resE = expl.Current() ; + if (myNewSection.Contains( resE )) // only new edges + MSE.Add(resE); + } + + // DMEF: map edge of CSF - faces of CSF + TopTools_IndexedDataMapOfShapeListOfShape DMEF; + TopExp::MapShapesAndAncestors(CSF, TopAbs_EDGE, TopAbs_FACE, DMEF); + + // Fill + // 1. MFP - a map of faces to process: map of resulting faces except + // those of theShape; we`ll add to C those of them which are inside CSF + // 2. DMSEFP - edge of MSE => faces of MFP + TopTools_ListIteratorOfListOfShape itl; + for (itl.Initialize(myListShapes);itl.More();itl.Next()) { + const TopoDS_Shape& aShape = itl.Value(); + if ( theShape.IsSame( aShape )) continue; + // fill maps + // iterate on split faces of aShape + TopoDS_Iterator itF ( myImageShape.Image(aShape).First() ); + for ( ; itF.More(); itF.Next()) { + const TopoDS_Shape& sf = itF.Value(); + MFP.Add(sf); + // iterate on edges of split faces of aShape, + // add to DMSEFP edges that are new + for (expl.Init( sf, TopAbs_EDGE ); expl.More(); expl.Next()) { + TopoDS_Shape se = expl.Current(); + if ( MSE.Contains(se)) {// section edge + if (!DMSEFP.IsBound(se)) + DMSEFP.Bind(se,EmptyL); + DMSEFP(se).Append(sf); + } + } + } + } + + // add tool faces having section edges on faces of theShape to MFP and DMSEFP; + // such tool faces need not to be reconstructed and so they are not in myListShapes + for (itm.Initialize(myMapTools); itm.More(); itm.Next()) + { + const TopoDS_Shape & aToolFace = itm.Key(); + if (myMapFaces.Contains( aToolFace )) + continue; + MFP.Add(aToolFace); + for (expl.Init( aToolFace, TopAbs_EDGE ); expl.More(); expl.Next()) { + TopoDS_Shape se = expl.Current(); + if ( MSE.Contains( se )) {// section edge + if (!DMSEFP.IsBound( se )) + DMSEFP.Bind( se, EmptyL ); + DMSEFP( se ).Append( aToolFace ); + } + } + } + + + // =========================== + // find faces inside theShape + // =========================== + + Standard_Boolean skipAlreadyAdded = Standard_False; + Standard_Boolean GoodOri, inside; + Standard_Real dot; + TopTools_ListOfShape KeepFaces; + TopTools_DataMapIteratorOfDataMapOfShapeListOfShape Mapit; + + // iterate on section edges, check faces of other shapes + // sharing section edges and put internal faces to KeepFaces + for (Mapit.Initialize(DMSEFP); Mapit.More() ; Mapit.Next() ) { + // a new edge of theShape + const TopoDS_Edge& E = TopoDS::Edge (Mapit.Key()); + // an original edge of which E is a split + const TopoDS_Edge& OrigE = TopoDS::Edge ( myImagesEdges.Root( E )); + // does OrigE itself splits a face + Standard_Boolean isSectionE = myInter3d.IsSectionEdge ( OrigE ); + + // split faces of other shapes sharing E + TopTools_ListOfShape& LSF = DMSEFP.ChangeFind(E); + itl.Initialize( LSF ); + while (itl.More()) { + // a split faces of other shape + TopoDS_Face aFace1 = TopoDS::Face(itl.Value()); + // remove aFace1 form DMSEFP and MFP + LSF.Remove( itl ); // == itl.Next(); + if (!MFP.Remove( aFace1 )) + continue; // was not is MFP ( i.e already checked) + // check if aFace1 was already added to 2 shells + if (!All && + myAddedFacesMap.Contains( aFace1 ) && + myAddedFacesMap.Contains( aFace1.Reversed() )) { + skipAlreadyAdded = Standard_True; + continue; + } + + // find another face which originates from the same face as aFace1: + // usually aFace2 is internal if aFace1 is not and vice versa + + TopoDS_Shape anOrigFace = aFace1; + if (myImagesFaces.IsImage(aFace1)) + anOrigFace = myImagesFaces.Root(aFace1); + TopoDS_Shape aFace2; + if ( !isSectionE ) { + while (itl.More()) { + aFace2 = itl.Value(); + if (!MFP.Contains( aFace2 )) { + LSF.Remove( itl ); + continue; + } + if (anOrigFace.IsSame( myImagesFaces.Root( aFace2 ))) + break; + itl.Next(); + } + if (itl.More()) { // aFace2 found, remove it from maps + LSF.Remove( itl ); + MFP.Remove(aFace2); + } + else + aFace2.Nullify(); + itl.Initialize( LSF ); + } + + // check that anOrigFace is not same domain with CSF faces it intersects + + const TopTools_ListOfShape& FL = DMEF.FindFromKey(E); //faces of CSF sharing E + const TopoDS_Shape& origF1 = myImagesFaces.Root(FL.First()); + const TopoDS_Shape& origF2 = myImagesFaces.Root(FL.Last()); + Standard_Boolean sameDom1 = anOrigFace.IsSame( origF1 ); + Standard_Boolean sameDom2 = anOrigFace.IsSame( origF2 ); + if (!(sameDom1 || sameDom2) && myInter3d.HasSameDomainF( anOrigFace )) { + sameDom1 = myInter3d.IsSameDomainF( anOrigFace, origF1); + if (origF1 == origF2) + sameDom2 = sameDom1; + else + myInter3d.IsSameDomainF( anOrigFace, origF2); + } + if (sameDom1 && sameDom2) + continue; + if ((sameDom1 || sameDom2)) { + inside = Partition_Loop3d::IsInside (E, + TopoDS::Face(FL.First()), + TopoDS::Face(FL.Last()), + 1, dot, GoodOri); + if (inside || (dot + Precision::Angular() >= 1.0)) + continue; // E is convex between origF1 and origF2 or they are tangent + } + + + // keep one of found faces + + //face of CSF sharing E + const TopoDS_Shape& aShapeFace = sameDom1 ? FL.Last() : FL.First(); + // analyse aFace1 state + inside = Partition_Loop3d::IsInside (E, TopoDS::Face(aShapeFace), aFace1, + 1, dot, GoodOri); + if (inside && isSectionE) + { + // aFace1 must be tested with both adjacent faces of CSF + const TopoDS_Shape& aShapeFace2 = sameDom1 ? FL.First() : FL.Last(); + if (aShapeFace2 != aShapeFace) + inside = Partition_Loop3d::IsInside (E, TopoDS::Face(aShapeFace2), aFace1, + 1, dot, GoodOri); + } + + // store internal face + if (inside) + KeepFaces.Append(aFace1); + + else if (!aFace2.IsNull()) + { + if (dot + Precision::Angular() >= 1.0) + { + // aFace2 state is not clear, it will be analysed alone, + // put it back to the maps + MFP.Add( aFace2 ); + LSF.Append( aFace2 ); + } + else + KeepFaces.Append(aFace2); + } + } + } + + // =================================================== + // add not distributed faces connected with KeepFaces + // =================================================== + + // ultimate list of internal faces + TopTools_ListOfShape KeptFaces; + + // add to MFP not split tool faces as well, they may be connected with + // tool faces interfering with theShape + for ( itm.Initialize(myMapTools); itm.More(); itm.Next() ) { + const TopoDS_Shape& aToolFace = itm.Key(); + if (!myImageShape.HasImage(aToolFace)) + MFP.Add (aToolFace); + } + + if (MFP.IsEmpty()) + KeptFaces.Append (KeepFaces); + + while (!KeepFaces.IsEmpty()) + { + // KeepEdges : map of edges of faces kept last time + TopTools_IndexedMapOfShape KeepEdges; + for ( itl.Initialize(KeepFaces); itl.More(); itl.Next() ) { + TopExp::MapShapes( itl.Value(), TopAbs_EDGE, KeepEdges); + KeptFaces.Append( itl.Value() ); + } + + KeepFaces.Clear(); + + // keep faces connected with already kept faces by KeepEdges + for ( itm.Initialize(MFP); itm.More(); itm.Next() ) { + const TopoDS_Shape& FP = itm.Key(); + for (expl.Init(FP,TopAbs_EDGE); expl.More(); expl.Next()) { + const TopoDS_Shape& se = expl.Current(); + if (!MSE.Contains(se) && KeepEdges.Contains(se) ) { + KeepFaces.Append(FP); + MFP.Remove(FP); + break; + } + } + } + } + + // =============================================================== + // here MFP contains faces outer of theShape and those of shapes + // which do not interfere with theShape at all and between which + // there may be those wrapped by theShape and whose faces may be + // needed to be returned as well + // =============================================================== + + Standard_Boolean isSolid = (theShape.ShapeType() == TopAbs_SOLID); + if (All || isSolid) // All is for sub-result removal + { + // loop on not used faces; checked faces will be removed from MFP + // during the loop + for ( itm.Initialize( MFP ); itm.More(); itm.Next() ) { + const TopoDS_Shape & aFace = itm.Key(); + + // a shape which aFace originates from + TopoDS_Shape anOrigShape = GetOriginalShape( aFace ); + + // find out if all split faces of anOrigShape are not in MFP + // and by the way remove them from MFP + Standard_Boolean isAllOut = Standard_True; + TopoDS_Shape aSplitFaces = anOrigShape; + if (myImageShape.HasImage(anOrigShape)) + aSplitFaces = myImageShape.Image(anOrigShape).First(); + + TopTools_ListOfShape aSplitFaceL; // faces candidate to be kept + for (expl.Init( aSplitFaces, TopAbs_FACE ); expl.More(); expl.Next()) + { + const TopoDS_Shape & aSpFace = expl.Current(); + // a tool face which became object has image but the whole tool shape has not + if (myImageShape.HasImage( aSpFace )) + { + TopExp_Explorer exF (myImageShape.Image( aSpFace ).First(), TopAbs_FACE ); + for ( ; exF.More(); exF.Next() ) + { + aSplitFaceL.Append( exF.Current() ); + if ( ! MFP.Remove( exF.Current() ) && isAllOut ) + // a shared face might be removed from MFP during a prev loop + isAllOut = mySharedFaces.Contains( exF.Current() ); + } + } + else + { + aSplitFaceL.Append( aSpFace ); + if ( ! MFP.Remove( aSpFace ) && isAllOut) + // a shared face might be removed from MFP during a prev loop + isAllOut = mySharedFaces.Contains( aSpFace ); + } + } + itm.Initialize( MFP ); // iterate remaining faces + if ( !isAllOut ) + continue; + + // classify anOrigShape against theShape + if (IsInside (anOrigShape, theShape)) + { + if (isSolid && myClosedShapes.Contains( anOrigShape )) + // to make a special care at solid reconstruction + myWrappingSolid.Add ( theShape ); + + // keep faces of an internal shape anOrigShape + KeptFaces.Append( aSplitFaceL ); + } + } + } + + // ==================================================== + // check if kept faces form a shell without free edges + // ==================================================== + + DMEF.Clear(); // edge - kept faces + MFP.Clear(); // reuse it for wrong faces + if (CheckClosed) { + for (itl.Initialize(KeptFaces); itl.More(); itl.Next() ) + TopExp::MapShapesAndAncestors(itl.Value(), TopAbs_EDGE, TopAbs_FACE, DMEF); + + Standard_Integer i, nb = DMEF.Extent(); + Standard_Boolean isClosed = Standard_False; + while (!isClosed) { + isClosed = Standard_True; + for (i=1; isClosed && i<=nb; ++i) { + const TopoDS_Shape& E = DMEF.FindKey( i ); + if (! BRep_Tool::Degenerated( TopoDS::Edge( E )) && + ! MSE.Contains( E )) + isClosed = ( DMEF(i).Extent() != 1 ); + } + if (!isClosed) { + const TopoDS_Shape& F = DMEF.FindFromIndex( i-1 ).First(); // bad face + MFP.Add( F ); + // remove bad face from DMEF + for (expl.Init( F, TopAbs_EDGE); expl.More(); expl.Next()) { + const TopoDS_Shape& E = expl.Current(); + TopTools_ListOfShape& FL = DMEF.ChangeFromKey( E ); + for (itl.Initialize( FL ); itl.More(); itl.Next() ) { + if ( F.IsSame( itl.Value() )) { + FL.Remove( itl ); + break; + } + } + } + } + } + } + + // ============== + // make a result + // ============== + + TopoDS_Compound C; + // compound of removed internal faces + TopoDS_Compound CNotCl; + + myBuilder.MakeCompound(C); + myBuilder.MakeCompound(CNotCl); + + // add to compounds + for (itl.Initialize(KeptFaces); itl.More(); itl.Next() ) + { + TopoDS_Shape & aIntFace = itl.Value(); + if (! MFP.Contains( aIntFace )) + myBuilder.Add( C, aIntFace); + else + myBuilder.Add( CNotCl, aIntFace); + } + + if (!skipAlreadyAdded && CheckClosed) + { + myInternalFaces.Bind( theShape, C ); + myIntNotClFaces.Bind( theShape, CNotCl ); + } + + return C; +} + +//======================================================================= +//function : MakeShell +//purpose : split S into compound of shells +//======================================================================= + +void Partition_Spliter::MakeShells(const TopoDS_Shape& S, + TopTools_ListOfShape& NS) +{ + Partition_Loop3d ShellMaker; + // get compound of split faces of S + const TopoDS_Shape& FacesComp = myImageShape.Image(S).First(); + ShellMaker.AddConstFaces( FacesComp ); + // add split faces inside S + if (myClosedShapes.Contains( S )) { + TopoDS_Shape InternalFacesComp = FindFacesInside(S, Standard_True); + ShellMaker.AddSectionFaces( InternalFacesComp ); + } + + NS = ShellMaker.MakeShells( myAddedFacesMap ); + + // Add faces added to new shell to myAddedFacesMap: + // avoid rebuilding twice commont part of 2 solids. + TopTools_ListIteratorOfListOfShape itS(NS); + while ( itS.More()) { + TopExp_Explorer expF (itS.Value(), TopAbs_FACE); + for (; expF.More(); expF.Next()) + myAddedFacesMap.Add (expF.Current()); + + itS.Next(); + } +} + +//======================================================================= +//function : findEqual +//purpose : compare edges of EL1 against edges of EL2, +// Result is in EMM binding EL1 edges to list of equal edges. +// Edges are considered equall only if they have same vertices. +// ==True makes consider same edges as equal +// Put in all equal edges +//======================================================================= + +static void findEqual (const TopTools_ListOfShape& EL1, + const TopTools_ListOfShape& EL2, + const Standard_Boolean addSame, + TopTools_DataMapOfShapeListOfShape& EEM, + TopTools_MapOfShape& AllEqMap) +{ + // map vertices to edges for EL2 + TopTools_DataMapOfShapeListOfShape VEM; + TopTools_ListIteratorOfListOfShape itE1, itE2(EL2); + TopoDS_Iterator itV; + TopTools_ListOfShape emptyL; + for (; itE2.More(); itE2.Next()) { + for (itV.Initialize( itE2.Value() ); itV.More(); itV.Next()) { + const TopoDS_Shape& V = itV.Value(); + if (! VEM.IsBound( V ) ) + VEM.Bind( V, emptyL); + VEM( V ).Append( itE2.Value()); + } + } + + gp_Vec D1, D2; + gp_Pnt P; + Standard_Real f,l,u,tol; + Handle(Geom_Curve) C1, C2; + Extrema_ExtPC Extrema; + TopoDS_Vertex V1, V2, V3, V4; + + AllEqMap.Clear(); + + for (itE1.Initialize(EL1); itE1.More(); itE1.Next()) { + const TopoDS_Edge& E1 = TopoDS::Edge( itE1.Value()); + if (BRep_Tool::Degenerated( E1 ) || AllEqMap.Contains (E1)) + continue; + TopExp::Vertices( E1, V1, V2 ); + + if (VEM.IsBound(V1)) + itE2.Initialize( VEM(V1) ); + for (; itE2.More(); itE2.Next()) { + const TopoDS_Edge& E2 = TopoDS::Edge( itE2.Value()); + if (BRep_Tool::Degenerated( E2 ) || AllEqMap.Contains (E2)) + continue; + + if (E1.IsSame(E2)) { + if (!addSame) + continue; + } + else { + TopExp::Vertices( E2, V3, V4); + if (!V2.IsSame(V4) && !V2.IsSame(V3)) + continue; + // E1 and E2 have same vertices + // check D1 at end points. + C2 = BRep_Tool::Curve( E2, f,l); + C1 = BRep_Tool::Curve( E1, f,l); + u = BRep_Tool::Parameter(V1,E1); + C1->D1(u, P, D1); + u = BRep_Tool::Parameter(V1.IsSame(V3) ? V3 : V4, E2); + C2->D1(u, P, D2); + D1.Normalize(); D2.Normalize(); + if (Abs(D1*D2) + Precision::Angular() < 1.0) + continue; + if (! V1.IsSame(V2)) { + u = BRep_Tool::Parameter(V2,E1); + C1->D1(u, P, D1); + u = BRep_Tool::Parameter(V2.IsSame(V3) ? V3 : V4, E2); + C2->D1(u, P, D2); + D1.Normalize(); D2.Normalize(); + if (Abs(D1*D2) + Precision::Angular() < 1.0) + continue; + } + // check distance at a couple of internal points + tol = Max(BRep_Tool::Tolerance(E1), + BRep_Tool::Tolerance(E2)); + GeomAdaptor_Curve AC1(C1); + Extrema.Initialize(AC1,f,l); + Standard_Boolean ok = Standard_True, hasMin = Standard_False; + BRep_Tool::Range( E2, f, l); + Standard_Integer i=1, nbi=3; + for (; iValue( f+(l-f)*i/nbi )); + Standard_Integer j=1, nbj=Extrema.NbExt(); + for (; j<=nbj && ok; ++j) { + if (Extrema.IsMin(j)) { + hasMin = Standard_True; + ok = Extrema.Value(j) <= tol; + } + } + } + if ( !hasMin || !ok) + continue; + } + // bind E2 to E1 in EEM + if (!EEM.IsBound(E1)) { + EEM.Bind (E1, emptyL); + AllEqMap.Add (E1); + } + EEM(E1).Append(E2); + AllEqMap.Add (E2); + } + } +} + +//======================================================================= +//function : MakeFaces +//purpose : split faces of S, return compound of new faces +//======================================================================= + +TopoDS_Shape Partition_Spliter::MakeFaces (const TopoDS_Shape& S) +{ + TopoDS_Compound C; + myBuilder.MakeCompound(C); + + TopTools_ListIteratorOfListOfShape itl, itNE; + + TopExp_Explorer exp(S,TopAbs_FACE); + for (; exp.More(); exp.Next()) { + + const TopoDS_Face& F = TopoDS::Face(exp.Current()); + + TopTools_ListOfShape LNF; + + if (myImagesFaces.HasImage( F )) { + myImagesFaces.LastImage( F, LNF ); + TopAbs_Orientation oriF = F.Orientation(); + for ( itl.Initialize( LNF ); itl.More(); itl.Next()) + itl.Value().Orientation( oriF ); + } + else { + + Partition_Loop2d loops; + loops.Init(F); + + TopTools_IndexedMapOfShape EM; + TopExp::MapShapes( F, TopAbs_EDGE, EM); + + TopTools_MapOfShape AddedEqualM, EqualSeamM; + Standard_Boolean needRebuild = Standard_False; + + // add splits to loops + + // LE: old edges + new not splitted edges + const TopTools_ListOfShape& LE = myAsDes->Descendant(F); + for (itl.Initialize(LE); itl.More(); itl.Next()) { + const TopoDS_Edge& E = TopoDS::Edge( itl.Value() ); + + Standard_Boolean isSectionE = myInter3d.IsSectionEdge(E); + Standard_Boolean isNewE = !EM.Contains( E ); + + // LSE: list of split edges + TopTools_ListOfShape LSE; + myImagesEdges.LastImage(E,LSE); // splits of E or E itself + + for (itNE.Initialize(LSE); itNE.More(); itNE.Next()) { + + TopoDS_Edge NE = TopoDS::Edge( itNE.Value() ); + Standard_Boolean isSameE = NE.IsSame ( E ); + + if ( isNewE || isSectionE || !isSameE) { + if (AddedEqualM.Contains( NE )) { + // a seam must be twice in a loop + if (!BRep_Tool::IsClosed( E, F ) || !EqualSeamM.Add( NE )) + continue; + } + + if (isNewE) { + if (isSectionE) { + if ( ! myInter3d.IsSplitOn( NE, E, F) ) + continue; + } + else { + TopoDS_Vertex V1,V2; + TopExp::Vertices(NE,V1,V2); + const TopTools_ListOfShape& EL1 = myAsDes->Ascendant(V1); + const TopTools_ListOfShape& EL2 = myAsDes->Ascendant(V2); + if ( EL1.Extent() < 2 && EL2.Extent() < 2 ) + continue; + } + } + else { + NE.Orientation( E.Orientation()); + if (!isSameE) { + // orient NE because it may be a split of other edge + Standard_Real f,l,u; + Handle(Geom_Curve) C3d = BRep_Tool::Curve( E,f,l ); + Handle(Geom_Curve) NC3d = BRep_Tool::Curve( NE,f,l); + if ( C3d != NC3d) { + gp_Vec D1, ND1; gp_Pnt P; + TopoDS_Vertex V = TopExp::FirstVertex(NE); + u = BRep_Tool::Parameter(V,NE); + NC3d->D1 (u, P, ND1); + u = BRep_Tool::Parameter(V,E); + C3d ->D1 (u, P, D1); + if (ND1.Dot(D1) < 0) + NE.Reverse(); + } + } + } + if (myEqualEdges.Contains( NE )) + AddedEqualM.Add( NE ); + + needRebuild = Standard_True; + } + + if (isNewE || isSectionE) + myNewSection.Add( NE ); + + if (isNewE) + loops.AddSectionEdge(NE); + else + loops.AddConstEdge(NE); + } + } + + //------------------- + // Build the faces. + //------------------- + + if (needRebuild) { + + loops.Perform(); + loops.WiresToFaces(myImagesEdges); + + LNF = loops.NewFaces(); + + myImagesFaces.Bind(F,LNF); + + // replace the result faces that have already been built + // during same domain faces reconstruction done earlier + if (myInter3d.HasSameDomainF( F )) + { + // build map edge to same domain faces: EFM + TopTools_IndexedDataMapOfShapeListOfShape EFM; + TopTools_MapOfShape SDFM; // avoid doubling + itl.Initialize( myInter3d.SameDomain( F )); + for (; itl.More(); itl.Next()) { + if ( !myImagesFaces.HasImage( itl.Value() )) + continue; + // loop on splits of a SD face + TopTools_ListIteratorOfListOfShape itNF; + itNF.Initialize (myImagesFaces.Image( itl.Value() )); + for ( ; itNF.More(); itNF.Next()) { + TopoDS_Shape SDF = itNF.Value(); + if (myImagesFaces.HasImage( SDF )) // already replaced + SDF = myImagesFaces.Image( SDF ).First(); + if (SDFM.Add (SDF)) + TopExp::MapShapesAndAncestors(SDF, TopAbs_EDGE, TopAbs_FACE, EFM); + } + } + // do replace faces in the LNF + TopTools_ListOfShape LOF; + if ( !EFM.IsEmpty() ) + itl.Initialize( LNF ); + while (itl.More()) { + const TopoDS_Shape& NF = itl.Value(); + TopExp_Explorer expE ( NF, TopAbs_EDGE ); + const TopoDS_Edge& E = TopoDS::Edge (expE.Current()); + if (EFM.Contains( E )) { + const TopTools_ListOfShape& SDFL = EFM.FindFromKey( E ); + TopoDS_Shape SDF = SDFL.First(); + Standard_Boolean GoodOri; + Standard_Real dot; + Partition_Loop3d::IsInside (E, TopoDS::Face(NF), TopoDS::Face(SDF), + 1, dot, GoodOri); + if (dot < 0) + { + // NF and SDF are on different side of E + if (SDFL.Extent() == 1) { + itl.Next(); + continue; + } + else + SDF = SDFL.Last(); // next face must be on the same side + } + gp_Vec V1 = Partition_Loop3d::Normal( E, TopoDS::Face( NF )); + gp_Vec V2 = Partition_Loop3d::Normal( E, TopoDS::Face( SDF )); + if (V1*V2 < 0) + SDF.Reverse(); + + if (!myImagesFaces.HasImage( NF )) + myImagesFaces.Bind( NF, SDF ); + + // mySharedFaces is used in FindFacesInside() + mySharedFaces.Add( SDF ); + + LOF.Prepend ( SDF ); + LNF.Remove (itl); + } + else + itl.Next(); + } + + LNF.Append (LOF); + } + } // if (needRebuild) + + else { + LNF.Append( F ); + myImagesFaces.Bind(F,LNF); + } + } // if (myImagesFaces.HasImage( F )) + + // fill the resulting compound + for (itl.Initialize(LNF); itl.More(); itl.Next()) + myBuilder.Add ( C, itl.Value()); + + } // loop on faces of S + + return C; +} + + +//======================================================================= +//function : Tri +//purpose : +//======================================================================= + +static void Tri(const TopoDS_Edge& E, + TopTools_SequenceOfShape& Seq, + const Partition_Inter3d & theInter3d) +{ + Standard_Boolean Invert = Standard_True; + Standard_Real U1,U2; + TopoDS_Vertex V1,V2; + + while (Invert) { + Invert = Standard_False; + for ( Standard_Integer i = 1; i < Seq.Length(); i++) { + + V1 = TopoDS::Vertex(Seq.Value(i)); + V2 = TopoDS::Vertex(Seq.Value(i+1)); + + V1.Orientation(TopAbs_INTERNAL); + V2.Orientation(TopAbs_INTERNAL); + + U1 = BRep_Tool::Parameter(V1,E); + U2 = BRep_Tool::Parameter(V2,E); + + if (IsEqual(U1,U2)) { + if (theInter3d.ReplaceSameDomainV( V1, E ).IsSame( V1 )) + Seq.Remove(i+1); // remove V2 + else + Seq.Remove(i); + i--; + continue; + } + if (U2 < U1) { + Seq.Exchange(i,i+1); + Invert = Standard_True; + } + } + } +} + +//======================================================================= +//function : MakeEdges +//purpose : cut E by vertices VOnE, return list of new edges NE +//======================================================================= + +void Partition_Spliter::MakeEdges (const TopoDS_Edge& E, + const TopTools_ListOfShape& VOnE, + TopTools_ListOfShape& NE ) const +{ + TopoDS_Edge WE = E; + WE.Orientation(TopAbs_FORWARD); + + Standard_Real U1,U2, f, l; + TopoDS_Vertex V1,V2,VF,VL; + + BRep_Tool::Range(WE,f,l); + TopExp::Vertices(WE,VF,VL); + + if (VOnE.Extent() < 3) { // do not rebuild not cut edge + if (( VF.IsSame( VOnE.First() ) && VL.IsSame( VOnE.Last() )) || + VL.IsSame( VOnE.First() ) && VF.IsSame( VOnE.Last() ) ) { + NE.Append( E ); + return; + } + } + + TopTools_SequenceOfShape SV; + TopTools_ListIteratorOfListOfShape itv(VOnE); + TopTools_MapOfOrientedShape VM( VOnE.Extent() ); + for (; itv.More(); itv.Next()) + if ( VM.Add( itv.Value() )) + SV.Append(itv.Value()); + + Tri( WE, SV, myInter3d ); + + if (SV.Length() < 3) { // do not rebuild not cut edge + if (( VF.IsSame( SV.First() ) && VL.IsSame( SV.Last() )) || + VL.IsSame( SV.First() ) && VF.IsSame( SV.Last() ) ) { + NE.Append( E ); + return; + } + } + + Standard_Integer iVer, NbVer = SV.Length(); + + + //---------------------------------------------------------------- + // Construction of the new edges . + //---------------------------------------------------------------- + + if (VF.IsSame(VL)) { // closed edge + if (NbVer==1) + SV.Append( SV.First() ); + else if (!SV.First().IsSame(SV.Last())) { + Standard_Boolean isFirst=0; + Standard_Real minDU = 1.e10; + TopoDS_Vertex endV = Partition_Inter2d::FindEndVertex(VOnE, f,l, E, isFirst,minDU); + if (endV.IsSame(SV.First())) + SV.Append(endV); + else if (endV.IsSame(SV.Last())) + SV.Prepend(endV); + else + MESSAGE ("END VERTEX IS IN SEQUNCE MIDDLE"); + } + NbVer = SV.Length(); + } + + for (iVer=1; iVer < NbVer; iVer++) { + V1 = TopoDS::Vertex(SV(iVer)); + V2 = TopoDS::Vertex(SV(iVer+1)); + + TopoDS_Shape NewEdge = WE.EmptyCopied(); + V1.Orientation(TopAbs_FORWARD); + myBuilder.Add (NewEdge,V1); + V2.Orientation(TopAbs_REVERSED); + myBuilder.Add (NewEdge,V2); + + if (iVer==1) + U1 = f; + else { + V1.Orientation(TopAbs_INTERNAL); + U1=BRep_Tool::Parameter(V1,WE); + } + if (iVer+1 == NbVer) + U2 = l; + else { + V2.Orientation(TopAbs_INTERNAL); + U2=BRep_Tool::Parameter(V2,WE); + } + if (Abs(U1-U2) <= Precision::PConfusion()) { + MESSAGE( "MakeEdges(), EQUAL PARAMETERS OF DIFFERENT VERTICES"); + continue; + } + TopoDS_Edge EE=TopoDS::Edge(NewEdge); + myBuilder.Range (EE,U1,U2); + + TopoDS_Edge NEdge = TopoDS::Edge(NewEdge); + myBuilder.SameParameter(NEdge,Standard_False); + + Standard_Real tol = 1.0e-2; + Standard_Boolean flag = BRep_Tool::SameParameter(NEdge); + if (!flag) { + BRepLib::SameParameter(NEdge,tol); + } + NE.Append(NEdge.Oriented(E.Orientation())); + } +} + +//======================================================================= +//function : MergeEqualEdges +//purpose : find equal edges, choose ones to keep and make +// them have pcurves on all faces they are shared by +//======================================================================= + +void Partition_Spliter::MergeEqualEdges (const TopTools_ListOfShape& LSE) +{ + // find equal edges + // map: edge - equal edges + TopTools_DataMapOfShapeListOfShape EEM( LSE.Extent() ); + findEqual (LSE, LSE, 0, EEM, myEqualEdges); + + TopTools_ListOfShape EEL; // list of equal edges + TopTools_DataMapIteratorOfDataMapOfShapeListOfShape itM (EEM); + for ( ; itM.More(); itM.Next()) { + EEL = itM.Value(); + EEL.Append( itM.Key() ); + + // choose an edge to keep, section edges have priority + TopoDS_Edge EKeep; + TopTools_ListIteratorOfListOfShape itEE (EEL); + for (; itEE.More(); itEE.Next()) { + EKeep = TopoDS::Edge( itEE.Value() ); + const TopoDS_Edge& EKeepOrig = TopoDS::Edge( myImagesEdges.Root( EKeep )); + if (myInter3d.IsSectionEdge( EKeepOrig )) + break; + } + + // update edge images and build pcurves + Standard_Real f,l, tol; + for (itEE.Initialize (EEL); itEE.More(); itEE.Next()) { + const TopoDS_Edge& E = TopoDS::Edge( itEE.Value() ); + if ( E.IsSame( EKeep )) + continue; + + // 1. build pcurves of the kept edge on faces where replaced edges exist + const TopoDS_Edge& EReplOrig = TopoDS::Edge( myImagesEdges.Root( E )); + TopTools_ListOfShape FL; + FL = myAsDes->Ascendant( EReplOrig ); + Standard_Integer iFace, iFirstSectionFace = FL.Extent() + 1; + // add faces where the replaced edge is a section edge + if (myInter3d.IsSectionEdge( EReplOrig )) { + TopTools_ListIteratorOfListOfShape seIt; + seIt.Initialize( myInter3d.SectionEdgeFaces ( EReplOrig )); + for ( ; seIt.More(); seIt.Next()) + FL.Append( seIt.Value() ); + } + // loop on faces + TopTools_ListIteratorOfListOfShape itF (FL); + for ( iFace = 1 ; itF.More(); itF.Next(), ++iFace ) { + const TopoDS_Face& F = TopoDS::Face( itF.Value()); + + Handle(Geom2d_Curve) pc = BRep_Tool::CurveOnSurface( EKeep, F, f,l); + if (pc.IsNull()) { + Handle(Geom_Curve) C3d = BRep_Tool::Curve( EKeep, f, l); + C3d = new Geom_TrimmedCurve( C3d, f,l); + pc = TopOpeBRepTool_CurveTool::MakePCurveOnFace (F,C3d,tol); + if (pc.IsNull()) { + MESSAGE (" CANT BUILD PCURVE "); + } + myBuilder.UpdateEdge( EKeep, pc, F, tol); + } + + if (iFace >= iFirstSectionFace || + !BRep_Tool::IsClosed( EReplOrig, F )) + continue; + + // build the second pcurve for a seam + TopoDS_Vertex V = TopExp::FirstVertex( EKeep ); + Standard_Real Ukeep = BRep_Tool::Parameter( V, EKeep ); + Standard_Real Urepl = BRep_Tool::Parameter( V, E ); + + TopoDS_Edge EReplRev = E; + EReplRev.Reverse(); + Handle(Geom2d_Curve) pcRepl1 = BRep_Tool::CurveOnSurface( E, F, f,l); + Handle(Geom2d_Curve) pcRepl2 = BRep_Tool::CurveOnSurface( EReplRev, F, f,l); + + gp_Pnt2d p1r, p2r, pk; + p1r = pcRepl1->Value( Urepl ); + p2r = pcRepl2->Value( Urepl ); + pk = pc->Value( Ukeep ); + + // suppose that pk is equal to either p1r or p2r + Standard_Boolean isUPeriod = + ( Abs( p1r.X() - p2r.X() ) > Abs( p1r.Y() - p2r.Y() )); + Standard_Boolean is1Equal; + if (isUPeriod) + is1Equal = ( Abs( p1r.X() - pk.X() ) < Abs( p2r.X() - pk.X() )); + else + is1Equal = ( Abs( p1r.Y() - pk.Y() ) < Abs( p2r.Y() - pk.Y() )); + + Handle(Geom2d_Curve) pc2 = Handle(Geom2d_Curve)::DownCast + ( pc->Translated( pk, is1Equal ? p2r : p1r ) ); + + if (E.Orientation() == TopAbs_REVERSED) + is1Equal = !is1Equal; + + if (is1Equal) + myBuilder.UpdateEdge( EKeep, pc, pc2, F, tol); + else + myBuilder.UpdateEdge( EKeep, pc2, pc, F, tol); + + } // loop on a Faces where a replaced edge exists + + + // 2. update edge images according to replacement + if (myImagesEdges.HasImage( E )) + myImagesEdges.Remove( E ); + myImagesEdges.Bind( E, EKeep ); + + } // loop on a list of equal edges EEL + } // loop on a map of equal edges EEM +} + +//======================================================================= +//function : KeepShapesInside +//purpose : remove shapes that are outside of S from resul +//======================================================================= + +void Partition_Spliter::KeepShapesInside (const TopoDS_Shape& S) +{ + TopoDS_Iterator it; + if (S.ShapeType() < TopAbs_SOLID) { // compound or compsolid + for (it.Initialize( S ); it.More(); it.Next()) + KeepShapesInside( it.Value()); + return; + } + + Standard_Boolean isTool = Standard_False; + if (!myImageShape.HasImage( S )) { + isTool = CheckTool( S ); + if (!isTool) return; + } + + // build map of internal faces + TopTools_IndexedMapOfShape MIF; + TopoDS_Shape IntFacesComp = FindFacesInside( S, Standard_False, Standard_True); + TopExp::MapShapes( IntFacesComp, TopAbs_FACE, MIF ); + + TopoDS_Compound C; + myBuilder.MakeCompound(C); + + TopAbs_ShapeEnum anInternalShapeType = TopAbs_SHAPE; + if (!MIF.IsEmpty()) + { + // leave in the result only those shapes having a face in MIF + for (it.Initialize( myShape ); it.More(); it.Next()) { + const TopoDS_Shape & aResShape = it.Value(); + TopExp_Explorer expResF( aResShape, TopAbs_FACE ); + for (; expResF.More(); expResF.Next()) { + if ( MIF.Contains( expResF.Current())) { + myBuilder.Add( C, aResShape ); + if (aResShape.ShapeType() < anInternalShapeType) + anInternalShapeType = aResShape.ShapeType(); + break; + } + } + } + } + + // may be S was not split by internal faces then it is missing + // in myShape, add it + if (!isTool && + (anInternalShapeType > TopAbs_SOLID || S.ShapeType() > TopAbs_SOLID)) + { + TopTools_IndexedMapOfShape MSF; // map of split faces of S + TopExp::MapShapes( myImageShape.Image(S).First(), TopAbs_FACE, MSF); + + // find a shape having all faces in MSF + for (it.Initialize( myShape ); it.More(); it.Next()) { + TopExp_Explorer expResF( it.Value(), TopAbs_FACE ); + for (; expResF.More(); expResF.Next()) { + if (! MSF.Contains( expResF.Current())) + break; + } + if (! expResF.More()) { + myBuilder.Add( C, it.Value() ); + break; + } + } + } + + myShape = C; +} + +//======================================================================= +//function : RemoveShapesInside +//purpose : remove shapes that are inside S from resul +//======================================================================= + +void Partition_Spliter::RemoveShapesInside (const TopoDS_Shape& S) +{ + TopoDS_Iterator it; + if (S.ShapeType() < TopAbs_SOLID) { // compound or compsolid + for (it.Initialize( S ); it.More(); it.Next()) + RemoveShapesInside( it.Value()); + return; + } + Standard_Boolean isTool = Standard_False; + if (!myImageShape.HasImage( S )) { + isTool = CheckTool( S ); + if (!isTool) return; + } + + TopoDS_Shape IntFacesComp = FindFacesInside( S, Standard_False, Standard_True); + TopTools_IndexedMapOfShape MIF; // map of internal faces + TopExp::MapShapes( IntFacesComp, TopAbs_FACE, MIF); + + if (MIF.IsEmpty()) return; + + // add to MIF split faces of S + if (myImageShape.HasImage(S)) + TopExp::MapShapes( myImageShape.Image(S).First(), TopAbs_FACE, MIF); + + // leave in the result only those shapes not having all face in MIF + + TopoDS_Compound C; + myBuilder.MakeCompound(C); + + // RMF : faces of removed shapes that encounter once + TopTools_MapOfShape RFM; + + for (it.Initialize( myShape ); it.More(); it.Next()) { + + TopExp_Explorer expResF( it.Value(), TopAbs_FACE ); + for (; expResF.More(); expResF.Next()) + if (!MIF.Contains( expResF.Current())) + break; + + if (expResF.More()) + // add shape to result + myBuilder.Add( C, it.Value() ); + else + // add faces of a removed shape to RFM + for (expResF.ReInit(); expResF.More(); expResF.Next()) { + const TopoDS_Shape& F = expResF.Current(); + if ( ! RFM.Remove ( F )) + RFM.Add( F ); + } + } + + if (!isTool) { + + // rebuild S, it must remain in the result + + Standard_Boolean isClosed = Standard_False; + switch (S.ShapeType()) { + case TopAbs_SOLID : + isClosed = Standard_True; break; + case TopAbs_SHELL: { + TopTools_IndexedDataMapOfShapeListOfShape MEF; + TopExp::MapShapesAndAncestors(S, TopAbs_EDGE, TopAbs_FACE, MEF); + Standard_Integer i; + for (i=1; isClosed && i<=MEF.Extent(); ++i) + isClosed = ( MEF(i).Extent() != 1 ); + break; + } + default: + isClosed = Standard_False; + } + if (isClosed) { + + // add to a new shape external faces of removed shapes, ie those in RFM + + TopoDS_Shell Shell; + myBuilder.MakeShell( Shell ); + + // exclude redundant internal face with edges encounterd only once + TopTools_IndexedDataMapOfShapeListOfShape MEF; + TopTools_MapIteratorOfMapOfShape itF (RFM); + for ( ; itF.More(); itF.Next()) + TopExp::MapShapesAndAncestors(itF.Key(), TopAbs_EDGE, TopAbs_FACE, MEF); + + // add only faces forming a closed shell + for (itF.Reset() ; itF.More(); itF.Next()) + { + TopExp_Explorer expE (itF.Key(), TopAbs_EDGE); + for (; expE.More(); expE.Next()) + if (MEF.FindFromKey(expE.Current()).Extent() == 1) + break; + if (!expE.More()) + myBuilder.Add( Shell, itF.Key()); + } + + if (S.ShapeType() == TopAbs_SOLID) { + TopoDS_Solid Solid; + myBuilder.MakeSolid( Solid ); + myBuilder.Add (Solid, Shell); + myBuilder.Add (C, Solid); + } + else + myBuilder.Add (C, Shell); + } + else { + if (myImageShape.HasImage( S )) { + for (it.Initialize( myImageShape.Image(S).First()); it.More(); it.Next()) + myBuilder.Add (C, it.Value()); + } + } + } + + myShape = C; +} + +//======================================================================= +//function : CheckTool +//purpose : Return True if is a tool shape. Prepare tool +// faces of for the search of internal faces. +//======================================================================= + +Standard_Boolean Partition_Spliter::CheckTool(const TopoDS_Shape& S) +{ + // suppose S has not an image + + Standard_Boolean isTool = Standard_False; + TopoDS_Compound C; + myBuilder.MakeCompound( C ); + + TopExp_Explorer expF( S, TopAbs_FACE); + for (; expF.More(); expF.Next()) { + + const TopoDS_Face& F = TopoDS::Face( expF.Current() ); + if (myMapTools.Contains( F )) + isTool = Standard_True; + else + continue; + + if (myImagesFaces.HasImage( F )) { + // F has been reconstructed + TopAbs_Orientation Fori = F.Orientation(); + TopTools_ListOfShape LNF; + myImagesFaces.LastImage( F, LNF); + TopTools_ListIteratorOfListOfShape itF (LNF); + for ( ; itF.More(); itF.Next()) + myBuilder.Add( C, itF.Value().Oriented(Fori) ); + continue; + } + + Standard_Boolean hasSectionE = myInter3d.HasSectionEdge( F ); + Standard_Boolean hasNewE = myAsDes->HasDescendant( F ); + if (!hasSectionE && !hasNewE) + { + // F intersects nothing + myBuilder.Add( C, F ); + continue; + } + + // make an image for F + + TopoDS_Face NF = F; + NF.Orientation(TopAbs_FORWARD); + NF = TopoDS::Face( NF.EmptyCopied() ); // make a copy + TopoDS_Wire NW; + myBuilder.MakeWire( NW ); + + // add edges, as less as possible + TopTools_ListOfShape NEL; + TopTools_ListIteratorOfListOfShape itNE; + if (hasSectionE) { + // add section edges + TopExp_Explorer expE; + for ( ; expE.More(); expE.Next()) { + if (! myImagesEdges.HasImage( expE.Current() )) + continue; + myImagesEdges.LastImage( expE.Current(), NEL ); + for ( itNE.Initialize( NEL ); itNE.More(); itNE.Next()) + myBuilder.Add ( NW, itNE.Value()); + } + } + if (hasNewE) { + // add new adges + NEL = myAsDes->Descendant( F ); + for ( itNE.Initialize( NEL ); itNE.More(); itNE.Next()) { + TopTools_ListOfShape SEL; // splits + myImagesEdges.LastImage( itNE.Value(), SEL ); + TopTools_ListIteratorOfListOfShape itSE (SEL); + for ( ; itSE.More(); itSE.Next()) + myBuilder.Add ( NW, itSE.Value()); + } + } + myBuilder.Add( NF, NW ); + myBuilder.Add (C, NF); + + NF.Orientation( F.Orientation() ); // NF is most probably invalid + myImagesFaces.Bind (F, NF); + } + if (isTool) + myImageShape.Bind (S, C); + + return isTool; +} + +//======================================================================= +//function : IsInside +//purpose : Return True if the first vertex of S1 inside S2. +// If S1.IsNull(), check infinite point against S2. +//======================================================================= + +Standard_Boolean Partition_Spliter::IsInside (const TopoDS_Shape& theS1, + const TopoDS_Shape& theS2) +{ + BRepClass3d_SolidClassifier aClassifier( theS2 ); + + TopExp_Explorer expl( theS1, TopAbs_VERTEX ); + if (!expl.More()) + aClassifier.PerformInfinitePoint( ::RealSmall()); + else + { + const TopoDS_Vertex & aVertex = TopoDS::Vertex( expl.Current() ); + aClassifier.Perform (BRep_Tool::Pnt( aVertex ), + BRep_Tool::Tolerance( aVertex )); + } + + return ( aClassifier.State() == TopAbs_IN ); +} + +//======================================================================= +//function : GetOriginalShape +//purpose : Return the shape aShape originates from. aShape +// should be a face or more complex result shape +//======================================================================= + +TopoDS_Shape Partition_Spliter::GetOriginalShape(const TopoDS_Shape& theShape) const +{ + TopoDS_Shape anOrigShape; + + TopExp_Explorer expl( theShape, TopAbs_FACE); + if (expl.More()) + { + + TopoDS_Shape aFace = expl.Current(); + if (myImagesFaces.IsImage( aFace )) + aFace = myImagesFaces.Root( aFace ); + anOrigShape = myFaceShapeMap.Find( aFace ); + } + return anOrigShape; +} + +//======================================================================= +//function : FindToolsToReconstruct +//purpose : find and store as objects tools which interfere +// with solids or are inside solids without +// an interference +//======================================================================= + +void Partition_Spliter::FindToolsToReconstruct() +{ + if (myMapTools.IsEmpty()) + return; + + Standard_Integer nbFoundTools = 0; + + // build edge - face map in order to detect interference with section edges + TopTools_IndexedDataMapOfShapeListOfShape EFM; + TopTools_MapIteratorOfMapOfShape aMapIt; + for (aMapIt.Initialize(myMapTools); aMapIt.More(); aMapIt.Next()) + TopExp::MapShapesAndAncestors( aMapIt.Key(), TopAbs_EDGE, TopAbs_FACE, EFM); + for (aMapIt.Initialize(myMapFaces); aMapIt.More(); aMapIt.Next()) + TopExp::MapShapesAndAncestors( aMapIt.Key(), TopAbs_EDGE, TopAbs_FACE, EFM); + + TopTools_MapOfShape aCurrentSolids, aCheckedShapes; + + // faces cut by new edges + TopTools_MapOfShape & aSectionFaces = myInter3d.TouchedFaces(); + + // keep solids interfering with each other in aCurrentSolids map + // and add tool faces intersecting solids as object shapes + + TopTools_ListIteratorOfListOfShape itS, itF, itCF, itE; + for (itS.Initialize( myListShapes ); itS.More(); itS.Next()) { + TopExp_Explorer expSo (itS.Value(), TopAbs_SOLID); + for (; expSo.More(); expSo.Next()) { + + // check if a solid has been already processed + const TopoDS_Shape & aSo = expSo.Current(); + if (!aCheckedShapes.Add( aSo )) + continue; + aCurrentSolids.Add( aSo ); + + // faces to check + TopTools_ListOfShape aFacesToCheck; + TopExp_Explorer exp( aSo, TopAbs_FACE ); + for ( ; exp.More(); exp.Next()) + aFacesToCheck.Append ( exp.Current()); + + // add other shapes interefering with a solid. + // iterate faces to check while appending new ones + for (itCF.Initialize (aFacesToCheck) ; itCF.More(); itCF.Next()) + { + const TopoDS_Shape& aCheckFace = itCF.Value(); +// if (!aCheckedShapes.Add( aCheckFace )) +// continue; + + // find faces interfering with aCheckFace + TopTools_ListOfShape anIntFaces; + + // ** 1. faces intersecting aCheckFace with creation of new edges on it + if ( myAsDes->HasDescendant( aCheckFace )) + { + // new edges on aCheckFace + const TopTools_ListOfShape& NEL = myAsDes->Descendant( aCheckFace ); + for (itE.Initialize( NEL); itE.More(); itE.Next()) + { + const TopoDS_Shape & aNewEdge = itE.Value(); + if (!aCheckedShapes.Add( aNewEdge )) + continue; + + // faces interfering by aNewEdge + itF.Initialize (myAsDes->Ascendant( aNewEdge )); + for (; itF.More(); itF.Next()) + if (aCheckFace != itF.Value()) + anIntFaces.Append( itF.Value() ); + + // ** 2. faces having section edge aNewEdge on aFacesToCheck + if (EFM.Contains( aNewEdge)) + { + itF.Initialize ( EFM.FindFromKey (itE.Value())); + for (; itF.More(); itF.Next()) + if (aCheckFace != itF.Value()) + anIntFaces.Append( itF.Value() ); + } + } + } + + // ** 3. faces cut by edges of aCheckFace + TopExp_Explorer expE (aCheckFace, TopAbs_EDGE); + for ( ; expE.More(); expE.Next()) + { + const TopoDS_Shape & aCheckEdge = expE.Current(); + if (aCheckedShapes.Add( aCheckEdge ) && + myInter3d.IsSectionEdge( TopoDS::Edge( aCheckEdge ))) + { + itF.Initialize( myInter3d.SectionEdgeFaces( TopoDS::Edge( aCheckEdge ))); + for (; itF.More(); itF.Next()) + if (aCheckFace != itF.Value()) + anIntFaces.Append( itF.Value() ); + } + } + + // process faces interfering with aCheckFace and shapes they + // belong to + for (itF.Initialize (anIntFaces); itF.More(); itF.Next()) + { + const TopoDS_Shape & F = itF.Value(); + if (! aCheckedShapes.Add( F )) + continue; + + Standard_Boolean isTool = myMapTools.Contains( F ); + if (isTool && + myFaceShapeMap( aCheckFace ).ShapeType() == TopAbs_SOLID ) + { + // a tool interfering with a solid + if (aSectionFaces.Contains( F )) + AddShape( F ); + ++ nbFoundTools; + if (nbFoundTools == myMapTools.Extent()) + return; + } + + const TopoDS_Shape & S = myFaceShapeMap( F ); + if (aCheckedShapes.Add( S )) + { + // a new shape interefering with aCurrentSolids is found + if (!isTool && S.ShapeType() == TopAbs_SOLID) + aCurrentSolids.Add ( S ); + // add faces to aFacesToCheck list + for ( exp.Init( S, TopAbs_FACE ); exp.More(); exp.Next()) + aFacesToCheck.Append ( exp.Current() ); + } + } + } // loop on aFacesToCheck + + // Here aCurrentSolids contains all solids interfering with each other. + // aCheckedShapes contains all faces belonging to shapes included + // in or interfering with aCurrentSolids or previously checked solids. + // Test if tool faces that do not interefere with other shapes are + // wrapped by any of aCurrentSolids + + TopTools_MapIteratorOfMapOfShape aSolidIt (aCurrentSolids); + for ( ; aSolidIt.More(); aSolidIt.Next()) + { + const TopoDS_Shape & aSolid = aSolidIt.Key(); + TopTools_MapOfShape aCheckedTools( myMapTools.Extent() ); + + TopTools_MapIteratorOfMapOfShape aToolIt (myMapTools); + for ( ; aToolIt.More(); aToolIt.Next()) + { + const TopoDS_Shape & aToolFace = aToolIt.Key(); + if (aCheckedShapes.Contains( aToolFace ) || // already found + aCheckedTools.Contains( aToolFace )) // checked against aSolid + continue; + + const TopoDS_Shape & aToolShape = myFaceShapeMap( aToolFace ); + TopExp_Explorer aToolFaceIt( aToolShape, TopAbs_FACE ); + + Standard_Boolean isInside = IsInside( aToolShape, aSolid ); + for ( ; aToolFaceIt.More(); aToolFaceIt.Next() ) + { + const TopoDS_Shape & aTool = aToolFaceIt.Current(); + aCheckedTools.Add( aTool ); + if (isInside) + { + if (aSectionFaces.Contains( aTool )) + AddShape( aTool ); + ++ nbFoundTools; + if (nbFoundTools == myMapTools.Extent()) + return; + aCheckedShapes.Add( aTool ); + } + } + } + } + + } // loop on solid shapes + } +} + +#endif diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occconstruction.cpp netgen-4.5_new/libsrc/occ/occconstruction.cpp --- netgen-4.5_orig/libsrc/occ/occconstruction.cpp 2005-12-06 18:15:53.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/occconstruction.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -28,8 +28,8 @@ #include #include #include -#include -#include +//#include +//#include #include #include namespace netgen diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occgenmesh.cpp netgen-4.5_new/libsrc/occ/occgenmesh.cpp --- netgen-4.5_orig/libsrc/occ/occgenmesh.cpp 2006-02-07 13:12:48.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/occgenmesh.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -28,7 +28,7 @@ return Point<3> (p.X(), p.Y(), p.Z()); } - void DivideEdge (TopoDS_Edge & edge, + static void DivideEdge (TopoDS_Edge & edge, ARRAY & ps, ARRAY & params, Mesh & mesh) @@ -49,23 +49,19 @@ hvalue[0] = 0; pnt = c->Value(s0); - double olddist = 0; - double dist = 0; - - for (int i = 1; i <= DIVIDEEDGESECTIONS; i++) + int i; + for (i = 1; i <= DIVIDEEDGESECTIONS; i++) { oldpnt = pnt; pnt = c->Value(s0+(i/double(DIVIDEEDGESECTIONS))*(s1-s0)); + double dist = pnt.Distance(oldpnt); hvalue[i] = hvalue[i-1] + 1.0/mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z()))* - pnt.Distance(oldpnt); + dist; //(*testout) << "mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z())) " << mesh.GetH(Point3d(pnt.X(), pnt.Y(), pnt.Z())) // << " pnt.Distance(oldpnt) " << pnt.Distance(oldpnt) << endl; - - olddist = dist; - dist = pnt.Distance(oldpnt); } // nsubedges = int(ceil(hvalue[DIVIDEEDGESECTIONS])); @@ -74,7 +70,7 @@ ps.SetSize(nsubedges-1); params.SetSize(nsubedges+1); - int i = 1; + i = 1; int i1 = 0; do { @@ -112,7 +108,7 @@ static void FindEdges (OCCGeometry & geom, Mesh & mesh) { - char * savetask = multithread.task; + const char * savetask = multithread.task; multithread.task = "Edge meshing"; (*testout) << "edge meshing" << endl; @@ -124,6 +120,7 @@ (*testout) << "nedges = " << nedges << endl; double eps = 1e-6 * geom.GetBoundingBox().Diam(); + double eps2 = eps * eps; for (int i = 1; i <= nvertices; i++) { @@ -133,7 +130,7 @@ bool exists = 0; if (merge_solids) for (PointIndex pi = 1; pi <= mesh.GetNP(); pi++) - if ( Dist2 (mesh[pi], Point<3>(mp)) < eps*eps) + if ( Dist2 (mesh[pi], Point<3>(mp)) < eps2) { exists = 1; break; @@ -163,6 +160,7 @@ { TopoDS_Face face = TopoDS::Face(exp1.Current()); int facenr = geom.fmap.FindIndex(face); + if ( facenr < 1 ) continue; if (face2solid[0][facenr-1] == 0) face2solid[0][facenr-1] = solidnr; @@ -184,6 +182,9 @@ int facenr = 0; int edgenr = 0; + // EAP, IMP [SALOME platform 0013410]. + // take into account nb of already meshed edges + edgenr = mesh.GetNSeg(); (*testout) << "faces = " << geom.fmap.Extent() << endl; int curr = 0; @@ -232,6 +233,11 @@ continue; } + // EAP, IMP [SALOME platform 0013410]. + // Do not divide already meshed edges + if ( geom.emap.FindIndex(edge) < 1 ) + continue; + if (geom.vmap.FindIndex(TopExp::FirstVertex (edge)) == geom.vmap.FindIndex(TopExp::LastVertex (edge))) { @@ -276,8 +282,8 @@ pnums.Last() = -1; for (PointIndex pi = 1; pi < first_ep; pi++) { - if (Dist2 (mesh[pi], fp) < eps*eps) pnums[0] = pi; - if (Dist2 (mesh[pi], lp) < eps*eps) pnums.Last() = pi; + if (Dist2 (mesh[pi], fp) < eps2) pnums[0] = pi; + if (Dist2 (mesh[pi], lp) < eps2) pnums.Last() = pi; } } @@ -287,7 +293,7 @@ bool exists = 0; int j; for (j = first_ep; j <= mesh.GetNP(); j++) - if ((mesh.Point(j)-Point<3>(mp[i-1])).Length() < eps) + if (Dist2(mesh.Point(j), Point<3>(mp[i-1])) < eps2) { exists = 1; break; @@ -394,7 +400,7 @@ int i, j, k; int changed; - char * savetask = multithread.task; + const char * savetask = multithread.task; multithread.task = "Surface meshing"; geom.facemeshstatus = 0; @@ -751,7 +757,7 @@ multithread.task = savetask; } - double ComputeH (double kappa) + static double ComputeH (double kappa) { double hret; kappa *= mparam.curvaturesafety; @@ -779,7 +785,7 @@ double nq = n*q; Point<3> p = p0 + 0.5*n; - double lambda = (p-l.p0)*n / nq; + double lambda = (fabs(nq) > 1e-10 ? (p-l.p0)*n / nq : -1); if (lambda >= 0 && lambda <= 1) { @@ -799,55 +805,55 @@ - void RestrictHTriangle (gp_Pnt2d & par0, gp_Pnt2d & par1, gp_Pnt2d & par2, - BRepLProp_SLProps * prop, Mesh & mesh, const double maxside, int depth, double h = 0) + static void RestrictHTriangle (gp_Pnt2d & par0, gp_Pnt2d & par1, gp_Pnt2d & par2, + BRepAdaptor_Surface& surf, Mesh & mesh, const double maxside, int depth, double h = 0) { - + BRepLProp_SLProps prop(surf, 2, 1e-5); gp_Pnt2d parmid; parmid.SetX(0.3*(par0.X()+par1.X()+par2.X())); parmid.SetY(0.3*(par0.Y()+par1.Y()+par2.Y())); - if (depth == 0) + //if (depth == 0) { double curvature = 0; - prop->SetParameters (parmid.X(), parmid.Y()); - if (!prop->IsCurvatureDefined()) + prop.SetParameters (parmid.X(), parmid.Y()); + if (!prop.IsCurvatureDefined()) { (*testout) << "curvature not defined!" << endl; return; } - curvature = max(fabs(prop->MinCurvature()), - fabs(prop->MaxCurvature())); + curvature = max(fabs(prop.MinCurvature()), + fabs(prop.MaxCurvature())); - prop->SetParameters (par0.X(), par0.Y()); - if (!prop->IsCurvatureDefined()) + prop.SetParameters (par0.X(), par0.Y()); + if (!prop.IsCurvatureDefined()) { (*testout) << "curvature not defined!" << endl; return; } - curvature = max(curvature,max(fabs(prop->MinCurvature()), - fabs(prop->MaxCurvature()))); + curvature = max(curvature,max(fabs(prop.MinCurvature()), + fabs(prop.MaxCurvature()))); - prop->SetParameters (par1.X(), par1.Y()); - if (!prop->IsCurvatureDefined()) + prop.SetParameters (par1.X(), par1.Y()); + if (!prop.IsCurvatureDefined()) { (*testout) << "curvature not defined!" << endl; return; } - curvature = max(curvature,max(fabs(prop->MinCurvature()), - fabs(prop->MaxCurvature()))); + curvature = max(curvature,max(fabs(prop.MinCurvature()), + fabs(prop.MaxCurvature()))); - prop->SetParameters (par2.X(), par2.Y()); - if (!prop->IsCurvatureDefined()) + prop.SetParameters (par2.X(), par2.Y()); + if (!prop.IsCurvatureDefined()) { (*testout) << "curvature not defined!" << endl; return; } - curvature = max(curvature,max(fabs(prop->MinCurvature()), - fabs(prop->MaxCurvature()))); + curvature = max(curvature,max(fabs(prop.MinCurvature()), + fabs(prop.MaxCurvature()))); //(*testout) << "curvature " << curvature << endl; @@ -886,51 +892,47 @@ pm1.SetX(0.5*(par0.X()+par2.X())); pm1.SetY(0.5*(par0.Y()+par2.Y())); pm2.SetX(0.5*(par1.X()+par0.X())); pm2.SetY(0.5*(par1.Y()+par0.Y())); - RestrictHTriangle (pm0, pm1, pm2, prop, mesh, 0.5*maxside, depth+1, h); - RestrictHTriangle (par0, pm1, pm2, prop, mesh, 0.5*maxside, depth+1, h); - RestrictHTriangle (par1, pm0, pm2, prop, mesh, 0.5*maxside, depth+1, h); - RestrictHTriangle (par2, pm1, pm0, prop, mesh, 0.5*maxside, depth+1, h); + RestrictHTriangle (pm0, pm1, pm2, surf, mesh, 0.5*maxside, depth+1, h); + RestrictHTriangle (par0, pm1, pm2, surf, mesh, 0.5*maxside, depth+1, h); + RestrictHTriangle (par1, pm0, pm2, surf, mesh, 0.5*maxside, depth+1, h); + RestrictHTriangle (par2, pm1, pm0, surf, mesh, 0.5*maxside, depth+1, h); } else { gp_Pnt pnt; Point3d p3d; - prop->SetParameters (parmid.X(), parmid.Y()); - pnt = prop->Value(); + surf.D0(parmid.X(), parmid.Y(), pnt); p3d = Point3d(pnt.X(), pnt.Y(), pnt.Z()); mesh.RestrictLocalH (p3d, h); - prop->SetParameters (par0.X(), par0.Y()); - pnt = prop->Value(); + surf.D0(par0.X(), par0.Y(), pnt); p3d = Point3d(pnt.X(), pnt.Y(), pnt.Z()); mesh.RestrictLocalH (p3d, h); - prop->SetParameters (par1.X(), par1.Y()); - pnt = prop->Value(); + surf.D0(par1.X(), par1.Y(), pnt); p3d = Point3d(pnt.X(), pnt.Y(), pnt.Z()); mesh.RestrictLocalH (p3d, h); - prop->SetParameters (par2.X(), par2.Y()); - pnt = prop->Value(); + surf.D0(par2.X(), par2.Y(), pnt); p3d = Point3d(pnt.X(), pnt.Y(), pnt.Z()); mesh.RestrictLocalH (p3d, h); - (*testout) << "p = " << p3d << ", h = " << h << ", maxside = " << maxside << endl; + //(*testout) << "p = " << p3d << ", h = " << h << ", maxside = " << maxside << endl; /* (*testout) << pnt.X() << " " << pnt.Y() << " " << pnt.Z() << endl; - prop->SetParameters (par0.X(), par0.Y()); - pnt = prop->Value(); + prop.SetParameters (par0.X(), par0.Y()); + pnt = prop.Value(); (*testout) << pnt.X() << " " << pnt.Y() << " " << pnt.Z() << endl; - prop->SetParameters (par1.X(), par1.Y()); - pnt = prop->Value(); + prop.SetParameters (par1.X(), par1.Y()); + pnt = prop.Value(); (*testout) << pnt.X() << " " << pnt.Y() << " " << pnt.Z() << endl; - prop->SetParameters (par2.X(), par2.Y()); - pnt = prop->Value(); + prop.SetParameters (par2.X(), par2.Y()); + pnt = prop.Value(); (*testout) << pnt.X() << " " << pnt.Y() << " " << pnt.Z() << endl; */ } @@ -970,7 +972,7 @@ if (mparam.uselocalh) { - char * savetask = multithread.task; + const char * savetask = multithread.task; multithread.percent = 0; mesh->SetLocalH (bb.PMin(), bb.PMax(), mparam.grading); @@ -1075,7 +1077,6 @@ if (triangulation.IsNull()) continue; BRepAdaptor_Surface sf(face, Standard_True); - BRepLProp_SLProps prop(sf, 2, 1e-5); int ntriangles = triangulation -> NbTriangles(); for (int j = 1; j <= ntriangles; j++) @@ -1096,7 +1097,7 @@ maxside = max (maxside, p[1].Distance(p[2])); //cout << "\rFace " << i << " pos11 ntriangles " << ntriangles << " maxside " << maxside << flush; - RestrictHTriangle (par[0], par[1], par[2], &prop, *mesh, maxside, 0); + RestrictHTriangle (par[0], par[1], par[2], sf, *mesh, maxside, 0); //cout << "\rFace " << i << " pos12 ntriangles " << ntriangles << flush; } } diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occgeom.cpp netgen-4.5_new/libsrc/occ/occgeom.cpp --- netgen-4.5_orig/libsrc/occ/occgeom.cpp 2006-01-25 16:35:50.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/occgeom.cpp 2010-06-23 12:38:22.000000000 +0400 @@ -7,6 +7,8 @@ #include "ShapeAnalysis_ShapeContents.hxx" #include "ShapeAnalysis_CheckSmallFace.hxx" #include "ShapeAnalysis_DataMapOfShapeListOfReal.hxx" +#include +#include #include "BRepAlgoAPI_Fuse.hxx" #include "BRepCheck_Analyzer.hxx" #include "BRepLib.hxx" @@ -16,11 +18,19 @@ #include "Partition_Spliter.hxx" //#include "VrmlAPI.hxx" //#include "StlAPI.hxx" +#include namespace netgen { + OCCGeometry::~OCCGeometry() + { + NCollection_DataMap::Iterator it(fclsmap); + for (; it.More(); it.Next()) + delete it.Value(); + } + void OCCGeometry :: PrintNrShapes () { TopExp_Explorer e; @@ -947,13 +957,13 @@ void OCCGeometry :: BuildVisualizationMesh () { - - cout << "Preparing visualization (deflection = " << vispar.occdeflection << ") ... " << flush; + double vispar_occdeflection = 0.01; + cout << "Preparing visualization (deflection = " << vispar_occdeflection << ") ... " << flush; BRepTools::Clean (shape); //WriteOCC_STL("test.stl"); - BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh (shape, vispar.occdeflection, true); + BRepMesh_IncrementalMesh::BRepMesh_IncrementalMesh (shape, vispar_occdeflection, true); cout << "done" << endl; @@ -973,8 +983,27 @@ } + void OCCGeometry::GetFaceTools(int surfi, Handle(ShapeAnalysis_Surface)& proj, + BRepTopAdaptor_FClass2d*& cls) const + { + //MSV: organize caching projector in the map + if (fprjmap.IsBound(surfi)) + { + proj = fprjmap.Find(surfi); + cls = fclsmap.Find(surfi); + } + else + { + const TopoDS_Face& aFace = TopoDS::Face(fmap(surfi)); + Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace); + proj = new ShapeAnalysis_Surface(aSurf); + fprjmap.Bind(surfi, proj); + cls = new BRepTopAdaptor_FClass2d(aFace,Precision::Confusion()); + fclsmap.Bind(surfi, cls); + } + } - void OCCGeometry :: Project (int surfi, Point<3> & p) const + bool OCCGeometry :: Project (int surfi, Point<3> & p, double& u, double& v) const { static int cnt = 0; if (++cnt % 1000 == 0) cout << "Project cnt = " << cnt << endl; @@ -983,18 +1012,22 @@ //(*testout) << "before " << pnt.X() << " "<< pnt.Y() << " "<< pnt.Z() << " " << endl; - GeomAPI_ProjectPointOnSurf proj(pnt, BRep_Tool::Surface(TopoDS::Face(fmap(surfi)))); - if (proj.NbPoints() == 0) - { - cout << "Projection fails" << endl; - } - else - { - pnt = proj.NearestPoint(); - //(*testout) << "after " << pnt.X() << " "<< pnt.Y() << " "<< pnt.Z() << " " << endl; + Handle(ShapeAnalysis_Surface) proj; + BRepTopAdaptor_FClass2d *cls; + GetFaceTools(surfi, proj, cls); - p = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); - } + gp_Pnt2d p2d = proj->ValueOfUV(pnt, Precision::Confusion()); + if (cls->Perform(p2d) == TopAbs_OUT) + { + //cout << "Projection fails" << endl; + return false; + } + pnt = proj->Value(p2d); + p2d.Coord(u, v); + //(*testout) << "after " << pnt.X() << " "<< pnt.Y() << " "<< pnt.Z() << " " << endl; + + p = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); + return true; } @@ -1002,54 +1035,20 @@ { gp_Pnt p(ap(0), ap(1), ap(2)); - Handle(Geom_Surface) surface = BRep_Tool::Surface(TopoDS::Face(fmap(surfi))); + Handle(ShapeAnalysis_Surface) proj; + BRepTopAdaptor_FClass2d *cls; + GetFaceTools(surfi, proj, cls); - gp_Pnt x = surface->Value (u,v); - - if (p.SquareDistance(x) <= sqr(PROJECTION_TOLERANCE)) return true; - - gp_Vec du, dv; - - surface->D1(u,v,x,du,dv); - - int count = 0; - - gp_Pnt xold; - gp_Vec n; - double det, lambda, mu; - - do { - count++; - - n = du^dv; - - det = Det3 (n.X(), du.X(), dv.X(), - n.Y(), du.Y(), dv.Y(), - n.Z(), du.Z(), dv.Z()); - - if (det < 1e-15) return false; - - lambda = Det3 (n.X(), p.X()-x.X(), dv.X(), - n.Y(), p.Y()-x.Y(), dv.Y(), - n.Z(), p.Z()-x.Z(), dv.Z())/det; - - mu = Det3 (n.X(), du.X(), p.X()-x.X(), - n.Y(), du.Y(), p.Y()-x.Y(), - n.Z(), du.Z(), p.Z()-x.Z())/det; - - u += lambda; - v += mu; - - xold = x; - surface->D1(u,v,x,du,dv); - - } while (xold.SquareDistance(x) > sqr(PROJECTION_TOLERANCE) && count < 50); - - // (*testout) << "FastProject count: " << count << endl; - - if (count == 50) return false; + gp_Pnt2d p2d = proj->NextValueOfUV(gp_Pnt2d(u,v), p, Precision::Confusion()); + if (cls->Perform(p2d) == TopAbs_OUT) + { + //cout << "Projection fails" << endl; + return false; + } - ap = Point<3> (x.X(), x.Y(), x.Z()); + p = proj->Value(p2d); + p2d.Coord(u, v); + ap = Point<3> (p.X(), p.Y(), p.Z()); return true; } @@ -1190,16 +1189,16 @@ return occgeo; } - char * shapesname[] = + const char * shapesname[] = {" ", "CompSolids", "Solids", "Shells", "Faces", "Wires", "Edges", "Vertices"}; - char * shapename[] = + const char * shapename[] = {" ", "CompSolid", "Solid", "Shell", "Face", "Wire", "Edge", "Vertex"}; - char * orientationstring[] = + const char * orientationstring[] = {"+", "-"}; void OCCGeometry :: RecursiveTopologyTree (const TopoDS_Shape & sh, diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occgeom.hpp netgen-4.5_new/libsrc/occ/occgeom.hpp --- netgen-4.5_orig/libsrc/occ/occgeom.hpp 2006-01-25 16:35:50.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/occgeom.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -15,8 +15,6 @@ #include "Geom_Curve.hxx" #include "Geom2d_Curve.hxx" #include "Geom_Surface.hxx" -#include "GeomAPI_ProjectPointOnSurf.hxx" -#include "GeomAPI_ProjectPointOnCurve.hxx" #include "BRepTools.hxx" #include "TopExp.hxx" #include "BRepBuilderAPI_MakeVertex.hxx" @@ -41,8 +39,6 @@ #include "Geom_Curve.hxx" #include "Geom2d_Curve.hxx" #include "Geom_Surface.hxx" -#include "GeomAPI_ProjectPointOnSurf.hxx" -#include "GeomAPI_ProjectPointOnCurve.hxx" #include "TopoDS_Wire.hxx" #include "BRepTools_WireExplorer.hxx" #include "BRepTools.hxx" @@ -69,7 +65,7 @@ #include "IGESToBRep_Reader.hxx" #include "Interface_Static.hxx" #include "GeomAPI_ExtremaCurveCurve.hxx" -#include "Standard_ErrorHandler.hxx" +//#include "Standard_ErrorHandler.hxx" #include "Standard_Failure.hxx" #include "ShapeUpgrade_ShellSewing.hxx" #include "ShapeFix_Shape.hxx" @@ -84,11 +80,15 @@ #include "STEPControl_Writer.hxx" #include "StlAPI_Writer.hxx" #include "STEPControl_StepModelType.hxx" +#include + +class Handle_ShapeAnalysis_Surface; +class BRepTopAdaptor_FClass2d; namespace netgen { -#include "../visualization/vispar.hpp" + //#include "../visualization/vispar.hpp" // class VisualizationParameters; // extern VisualizationParameters vispar; @@ -159,6 +159,8 @@ class OCCGeometry { Point<3> center; + mutable NCollection_DataMap fprjmap; + mutable NCollection_DataMap fclsmap; public: TopoDS_Shape shape; @@ -189,6 +191,7 @@ vmap.Clear(); } + ~OCCGeometry(); void BuildFMap(); @@ -204,10 +207,12 @@ Point<3> Center() { return center; } - void Project (int surfi, Point<3> & p) const; + bool Project (int surfi, Point<3> & p, double& u, double& v) const; bool FastProject (int surfi, Point<3> & ap, double& u, double& v) const; - + void GetFaceTools(int surfi, Handle(ShapeAnalysis_Surface)& proj, + BRepTopAdaptor_FClass2d*& cls) const; + OCCSurface GetSurface (int surfi) { cout << "OCCGeometry::GetSurface using PLANESPACE" << endl; diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occmeshsurf.cpp netgen-4.5_new/libsrc/occ/occmeshsurf.cpp --- netgen-4.5_orig/libsrc/occ/occmeshsurf.cpp 2006-01-25 16:36:26.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/occmeshsurf.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -5,6 +5,8 @@ #include #include #include +#include +#include namespace netgen @@ -411,11 +413,16 @@ } - void MeshOptimize2dOCCSurfaces :: ProjectPoint (INDEX surfind, Point3d & p) const + bool MeshOptimize2dOCCSurfaces :: ProjectPoint (INDEX surfind, Point3d & p, PointGeomInfo& gi) const { Point<3> hp = p; - geometry.Project (surfind, hp); + bool ok; + if (gi.trignum > 0) + ok = geometry.FastProject (surfind, hp, gi.u, gi.v); + else + ok = geometry.Project (surfind, hp, gi.u, gi.v); p = hp; + return ok; } void MeshOptimize2dOCCSurfaces :: ProjectPoint2 (INDEX surfind, INDEX surfind2, @@ -506,38 +513,6 @@ } - int MeshOptimize2dOCCSurfaces :: - CalcPointGeomInfo(int surfind, PointGeomInfo& gi, const Point3d& p) const - { - Standard_Real u,v; - - gp_Pnt pnt(p.X(), p.Y(), p.Z()); - - Handle(Geom_Surface) occface; - occface = BRep_Tool::Surface(TopoDS::Face(geometry.fmap(surfind))); - - GeomAPI_ProjectPointOnSurf proj(pnt, occface); - - if (proj.NbPoints() < 1) - { - cout << "ERROR: OCCSurface :: GetNormalVector: GeomAPI_ProjectPointOnSurf failed!" - << endl; - cout << p << endl; - return 0; - } - - proj.LowerDistanceParameters (u, v); - - gi.u = u; - gi.v = v; - return 1; - } - - - - - - OCCRefinementSurfaces :: OCCRefinementSurfaces (const OCCGeometry & ageometry) : Refinement(), geometry(ageometry) { @@ -627,10 +602,11 @@ if (!geometry.FastProject (surfi, hnewp, u, v)) { cout << "Fast projection to surface fails! Using OCC projection" << endl; - geometry.Project (surfi, hnewp); + double u, v; + geometry.Project (surfi, hnewp, u, v); } - newgi.trignum = 1; + newgi.trignum = surfi; } newp = hnewp; @@ -653,14 +629,17 @@ hnewp = Point<3> (pnt.X(), pnt.Y(), pnt.Z()); newp = hnewp; newgi = ap1; - }; + } void OCCRefinementSurfaces :: ProjectToSurface (Point<3> & p, int surfi) { if (surfi > 0) - geometry.Project (surfi, p); - }; + { + double u, v; + geometry.Project (surfi, p, u, v); + } + } void OCCRefinementSurfaces :: ProjectToSurface (Point<3> & p, int surfi, PointGeomInfo & gi) { @@ -668,9 +647,10 @@ if (!geometry.FastProject (surfi, p, gi.u, gi.v)) { cout << "Fast projection to surface fails! Using OCC projection" << endl; - geometry.Project (surfi, p); + double u, v; + geometry.Project (surfi, p, u, v); } - }; + } diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/occmeshsurf.hpp netgen-4.5_new/libsrc/occ/occmeshsurf.hpp --- netgen-4.5_orig/libsrc/occ/occmeshsurf.hpp 2005-06-09 18:51:10.000000000 +0400 +++ netgen-4.5_new/libsrc/occ/occmeshsurf.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -151,7 +151,7 @@ MeshOptimize2dOCCSurfaces (const OCCGeometry & ageometry); /// - virtual void ProjectPoint (INDEX surfind, Point3d & p) const; + virtual bool ProjectPoint (INDEX surfind, Point3d & p, PointGeomInfo& gi) const; /// virtual void ProjectPoint2 (INDEX surfind, INDEX surfind2, Point3d & p) const; /// @@ -159,9 +159,6 @@ /// virtual void GetNormalVector(INDEX surfind, const Point3d & p, PointGeomInfo & gi, Vec3d & n) const; - - virtual int CalcPointGeomInfo(int surfind, PointGeomInfo& gi, const Point3d& p3) const; - }; diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/occ/utilities.h netgen-4.5_new/libsrc/occ/utilities.h --- netgen-4.5_orig/libsrc/occ/utilities.h 2005-02-11 14:35:43.000000000 +0300 +++ netgen-4.5_new/libsrc/occ/utilities.h 2010-06-23 13:19:48.000000000 +0400 @@ -33,6 +33,7 @@ #include #include +#include #include // #include "SALOME_Log.hxx" diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/stlgeom/meshstlsurface.cpp netgen-4.5_new/libsrc/stlgeom/meshstlsurface.cpp --- netgen-4.5_orig/libsrc/stlgeom/meshstlsurface.cpp 2006-01-11 19:08:20.000000000 +0300 +++ netgen-4.5_new/libsrc/stlgeom/meshstlsurface.cpp 2010-06-23 13:19:48.000000000 +0400 @@ -946,20 +946,23 @@ } -void MeshOptimizeSTLSurface :: ProjectPoint (INDEX surfind, Point3d & p) const +bool MeshOptimizeSTLSurface :: ProjectPoint (INDEX surfind, Point3d & p, PointGeomInfo& gi) const { Point<3> hp = p; - if (!geom.Project (hp)) + if (gi.trignum > 0) + ((STLGeometry&)geom).SelectChartOfTriangle (gi.trignum); + if (!(gi.trignum = geom.Project (hp))) { PrintMessage(7,"project failed"); - if (!geom.ProjectOnWholeSurface(hp)) + if (!(gi.trignum = geom.ProjectOnWholeSurface(hp))) { PrintMessage(7, "project on whole surface failed"); } } p = hp; // geometry.GetSurface(surfind)->Project (p); + return gi.trignum > 0; } void MeshOptimizeSTLSurface :: ProjectPoint2 (INDEX surfind, INDEX surfind2, Point3d & p) const @@ -970,20 +973,6 @@ */ } -int MeshOptimizeSTLSurface :: CalcPointGeomInfo(PointGeomInfo& gi, const Point3d& p3) const -{ - Point<3> hp = p3; - gi.trignum = geom.Project (hp); - - if (gi.trignum) - { - return 1; - } - - return 0; - -} - void MeshOptimizeSTLSurface :: GetNormalVector(INDEX surfind, const Point3d & p, Vec3d & n) const { n = geom.GetChartNormalVector(); diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/stlgeom/meshstlsurface.hpp netgen-4.5_new/libsrc/stlgeom/meshstlsurface.hpp --- netgen-4.5_orig/libsrc/stlgeom/meshstlsurface.hpp 2004-09-30 17:13:56.000000000 +0400 +++ netgen-4.5_new/libsrc/stlgeom/meshstlsurface.hpp 2010-06-23 13:19:48.000000000 +0400 @@ -79,12 +79,10 @@ virtual void SelectSurfaceOfPoint (const Point3d & p, const PointGeomInfo & gi); /// - virtual void ProjectPoint (INDEX surfind, Point3d & p) const; + virtual bool ProjectPoint (INDEX surfind, Point3d & p, PointGeomInfo& gi) const; /// virtual void ProjectPoint2 (INDEX surfind, INDEX surfind2, Point3d & p) const; /// - virtual int CalcPointGeomInfo(PointGeomInfo& gi, const Point3d& p3) const; - /// virtual void GetNormalVector(INDEX surfind, const Point3d & p, Vec3d & n) const; }; diff -Naur --exclude=CVS netgen-4.5_orig/libsrc/stlgeom/stlgeommesh.cpp netgen-4.5_new/libsrc/stlgeom/stlgeommesh.cpp --- netgen-4.5_orig/libsrc/stlgeom/stlgeommesh.cpp 2004-08-10 03:39:45.000000000 +0400 +++ netgen-4.5_new/libsrc/stlgeom/stlgeommesh.cpp 2010-06-23 12:39:38.000000000 +0400 @@ -1437,7 +1437,7 @@ if (!optstring || strlen(optstring) == 0) { - mparam.optimize2d = "smcm"; + mparam.optimize2d = (char*)"smcm"; } else { @@ -1453,7 +1453,7 @@ mparam.grading); mesh -> LoadLocalMeshSize (mparam.meshsizefilename); mesh -> CalcLocalHFromSurfaceCurvature (stlparam.resthsurfmeshcurvfac); - mparam.optimize2d = "cmsmSm"; + mparam.optimize2d = (char*)"cmsmSm"; STLSurfaceOptimization (*stlgeometry, *mesh, mparam); #ifdef STAT_STREAM (*statout) << GetTime() << " & "; @@ -1559,7 +1559,7 @@ if (!optstring || strlen(optstring) == 0) { - mparam.optimize3d = "cmdmstm"; + mparam.optimize3d = (char*)"cmdmstm"; } else { diff -Naur --exclude=CVS netgen-4.5_orig/makeForSalome.sh netgen-4.5_new/makeForSalome.sh --- netgen-4.5_orig/makeForSalome.sh 1970-01-01 03:00:00.000000000 +0300 +++ netgen-4.5_new/makeForSalome.sh 2010-06-23 13:19:48.000000000 +0400 @@ -0,0 +1,31 @@ +#! /bin/sh +cp ngtcltk/ngnewdelete.* libsrc/interface/ + +if test `uname -m` == "x86_64" ; then + MACHINE=LINUX64 +else + MACHINE=LINUX +fi +export MACHINE +make -C libsrc/csg +make -C libsrc/general +make -C libsrc/geom2d +make -C libsrc/gprim +make -C libsrc/interface +make -C libsrc/linalg +make -C libsrc/meshing +make -C libsrc/opti +make -C libsrc/stlgeom +make -C libsrc/occ + +if [ ! -d install ] ; then + mkdir install +fi + +cp -r lib install/ + +if [ ! -d install/include ] ; then + mkdir install/include +fi + +cp libsrc/interface/nglib.h libsrc/general/*.hpp libsrc/csg/*.hpp libsrc/geom2d/*.hpp \