[occt.git] / src / IntStart / IntStart_SearchInside.gxx

// Copyright (c) 1995-1999 Matra Datavision
// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.

#ifndef OCCT_DEBUG
#define No_Standard_RangeError
#define No_Standard_OutOfRange
#endif


#include <math_FunctionSetRoot.hxx>
#include <Precision.hxx>
#include <gp_Pnt2d.hxx>

#include <TopAbs_State.hxx>

IntStart_SearchInside::IntStart_SearchInside (): done(Standard_False)
{}

IntStart_SearchInside::IntStart_SearchInside (TheFunction& Func,
                                              const ThePSurface& PS,
					      const Handle(TheTopolTool)& T,
					      const Standard_Real Epsilon) {

  Perform(Func,PS,T,Epsilon);
}


//=======================================================================
//function : Perform
//purpose  : Search all inside points
//=======================================================================

void IntStart_SearchInside::Perform (TheFunction& Func,
				     const ThePSurface& PS,
				     const Handle(TheTopolTool)& T,
				     const Standard_Real Epsilon) {

  done = Standard_False;
  list.Clear();
  Standard_Real aBinf[2], aBsup[2], aUVap[2], atoler[2];
  math_Vector Binf(aBinf,1,2), Bsup(aBsup,1,2), UVap(aUVap,1,2), toler(atoler,1,2);
  gp_Pnt psol;
  Standard_Boolean testpnt;
  Standard_Integer i,j,nbpt;
  TopAbs_State situ;
  Standard_Real umin,umax,vmin,vmax;
  Binf(1) = umin = ThePSurfaceTool::FirstUParameter(PS);
  Binf(2) = vmin = ThePSurfaceTool::FirstVParameter(PS);
  Bsup(1) = umax = ThePSurfaceTool::LastUParameter(PS);
  Bsup(2) = vmax = ThePSurfaceTool::LastVParameter(PS);
  
  Standard_Integer NbsampleU= T->NbSamplesU();
  Standard_Integer NbsampleV= T->NbSamplesV();
  Standard_Integer Nbsample = T->NbSamples();
  
  Standard_Real du = Bsup(1)-Binf(1);
  Standard_Real dv = Bsup(2)-Binf(2);  
  du/=(Standard_Real)NbsampleU*0.5;
  dv/=(Standard_Real)NbsampleV*0.5;
  
  Standard_Real toler1 = toler(1) = ThePSurfaceTool::UResolution(PS,Precision::Confusion());
  Standard_Real toler2 = toler(2) = ThePSurfaceTool::VResolution(PS,Precision::Confusion());
  Standard_Real Maxtoler1toler2 = toler1;
  if(toler2>Maxtoler1toler2) Maxtoler1toler2 = toler2;

  //-- lbr le 15 mai 97
  //-- on interdit aux points d'etre trop prets des restrictions 
  Maxtoler1toler2*=1000;
  if(Maxtoler1toler2>du*0.001) Maxtoler1toler2=du*0.001;
  if(Maxtoler1toler2>dv*0.001) Maxtoler1toler2=dv*0.001;
  
  
  Func.Set(PS);
  Standard_Real Tol = Func.Tolerance();

  math_FunctionSetRoot Rsnld(Func,toler);

  //-- lbr le 15 mai 97 
  umin+=du*0.01;
  vmin+=dv*0.01;
  umax-=du*0.01;
  vmax-=dv*0.01;

  //-- lbr le 30 octobre 97 : 
  //-- Si une surface vient tangenter 2 edges proche d un coin
  //-- il faut faire attention qu un point de depart soit trouve au 
  //-- voisinage du coin. Car ds le cas contraire, le cheminement ne
  //-- pourra pas passer au travers des frontieres : 
  //--
  //-- typiquement I est un cylindre (conge)
  //--
  //--                 PPPPPPPPPPPPPPPPPPPP*PPPPPPPPPPPPPPPP
  //--                 P               I     I
  //--                 P           I           I
  //--                 P        I
  //--                 P      #  il faut trouver un point ici
  //--                 P    I
  //--                 P  I
  //--                 PI
  //--                 *                    I
  //--                 PI                 I 
  //--                 P  I I I I I I I I
  //--

  
  for (i=1; i <= Nbsample+12; i++) {
    gp_Pnt2d s2d;
    gp_Pnt s3d;
    Standard_Boolean nepastester=Standard_False;
    if(i<=Nbsample) { 
      T->SamplePoint(i,s2d,s3d);      
      UVap(1)=s2d.X(); UVap(2)=s2d.Y();
      
      Standard_Real u1,v1,u2,v2;
      u1 = Binf(1) = Max(umin,UVap(1)-du);
      v1 = Binf(2) = Max(vmin,UVap(2)-dv);
      u2 = Bsup(1) = Min(umax,UVap(1)+du);
      v2 = Bsup(2) = Min(vmax,UVap(2)+dv);
      
      
      //-- gp_Pnt Pmilieu = ThePSurfaceTool::Value(PS,0.5*(u1+u2),0.5*(v1+v2));
      gp_Pnt Pextrm1 = ThePSurfaceTool::Value(PS,u1,v1);
      gp_Pnt Pextrm2 = ThePSurfaceTool::Value(PS,u2,v2);
      Standard_Real aValf[1];
      math_Vector Valf(aValf,1,1);
      Func.Value(UVap,Valf);
      Standard_Real rvalf = Valf(1);
      Standard_Real DistPP = Pextrm1.SquareDistance(Pextrm2);
      if(rvalf*rvalf > 3.0*DistPP) { 
	nepastester=Standard_True;
      }      
    }
    else { 
                 if(i==Nbsample+1) { s2d.SetCoord(umin+du*0.02,vmin+dv*0.02); }
      else       if(i==Nbsample+2) { s2d.SetCoord(umax-du*0.02,vmin+dv*0.02);  }
      else       if(i==Nbsample+3) { s2d.SetCoord(umin+du*0.02,vmax-dv*0.02); }
      else       if(i==Nbsample+4) { s2d.SetCoord(umax-du*0.02,vmax-dv*0.02);  }

      else       if(i==Nbsample+5) { s2d.SetCoord(umin+du*0.02,vmin+dv*0.02); }
      else       if(i==Nbsample+6) { s2d.SetCoord(umax-du*0.02,vmin+dv*0.02);  }
      else       if(i==Nbsample+7) { s2d.SetCoord(umin+du*0.02,vmax-dv*0.02); }
      else       if(i==Nbsample+8) { s2d.SetCoord(umax-du*0.02,vmax-dv*0.02);  }

      else       if(i==Nbsample+9) { s2d.SetCoord(umin+du*0.005,vmin+dv*0.005); }
      else       if(i==Nbsample+10){ s2d.SetCoord(umax-du*0.005,vmin+dv*0.005);  }
      else       if(i==Nbsample+11){ s2d.SetCoord(umin+du*0.005,vmax-dv*0.005); }
      else                         { s2d.SetCoord(umax-du*0.005,vmax-dv*0.005);  }

      UVap(1)=s2d.X(); UVap(2)=s2d.Y();

      Binf(1) = Max(umin,UVap(1)-du);
      Binf(2) = Max(vmin,UVap(2)-dv);
      Bsup(1) = Min(umax,UVap(1)+du);
      Bsup(2) = Min(vmax,UVap(2)+dv);      
    }
    

    if(nepastester==Standard_False) { 
      Rsnld.Perform(Func,UVap,Binf,Bsup);
      if (Rsnld.IsDone()) {
	if (Abs(Func.Root()) <= Tol) {
	  if (!Func.IsTangent()) {
	    psol = Func.Point();
	    Rsnld.Root(UVap);
	    // On regarde si le point trouve est bien un nouveau point.	  
	    j = 1;
	    nbpt = list.Length();
	    testpnt = (j <= nbpt);
	    
	    while (testpnt) {
	      const IntSurf_InteriorPoint&  IPj = list(j);
	      const gp_Pnt& Pj = IPj.Value();
	      if (   (Abs(Pj.X()-psol.X()) <= Epsilon)
		  && (Abs(Pj.Y()-psol.Y()) <= Epsilon)
		  && (Abs(Pj.Z()-psol.Z()) <= Epsilon)
		  && (Abs(UVap(1)-IPj.UParameter()) <= toler1)
		  && (Abs(UVap(2)-IPj.VParameter()) <= toler2) ) {
		testpnt = Standard_False;
	      }
	      else {
		j = j+1;
		testpnt = (j <= nbpt);
	      }
	    }
	    if (j > nbpt) {
	      //	    situ = TheSITool::Classify(PS,UVap(1),UVap(2));
	      situ = T->Classify(gp_Pnt2d(UVap(1),UVap(2)),
				 Maxtoler1toler2,Standard_False);   //-- ,Standard_False pour ne pas recadrer on Periodic
	      if (situ == TopAbs_IN) {
		list.Append(IntSurf_InteriorPoint(psol,UVap(1),UVap(2),
						  Func.Direction3d(),
						  Func.Direction2d()));
	      }
	    }
	  }
	}
      }
    }
  }
  //-- printf("\n Total : %d    Rejet : %d   RatioPointCalc : %g   nbpt =%d\n",REJET_OK+REJET_KO,REJET_OK,(double)(REJET_KO)/(double)(REJET_OK+REJET_KO),list.Length());
  done = Standard_True;
}


//=======================================================================
//function : Perform
//purpose  : Test the given inside point
//=======================================================================

void IntStart_SearchInside::Perform (TheFunction& Func,
				     const ThePSurface& PS,
				     const Standard_Real UStart,
				     const Standard_Real VStart)
{
  done = Standard_False;
  list.Clear();
  math_Vector Binf(1,2), Bsup(1,2), toler(1,2);

  Binf(1) = ThePSurfaceTool::FirstUParameter(PS);
  Binf(2) = ThePSurfaceTool::FirstVParameter(PS);
  Bsup(1) = ThePSurfaceTool::LastUParameter(PS);
  Bsup(2) = ThePSurfaceTool::LastVParameter(PS);

  toler(1) = ThePSurfaceTool::UResolution(PS,Precision::Confusion());
  toler(2) = ThePSurfaceTool::VResolution(PS,Precision::Confusion());

  if (UStart-Binf(1) > -toler(1) && UStart-Bsup(1) < toler(1) &&
      VStart-Binf(2) > -toler(2) && VStart-Bsup(2) < toler(2)) {

    Func.Set(PS);
    math_Vector UVap(1,2);
    UVap(1)=UStart; UVap(2)=VStart;

    math_FunctionSetRoot Rsnld(Func,toler);
    Rsnld.Perform(Func,UVap,Binf,Bsup);
    if (Rsnld.IsDone()) {
      Standard_Real tol = Func.Tolerance();
      Standard_Real valf = Func.Root();
      if (Abs(valf) <= tol && !Func.IsTangent()) {
	const gp_Pnt& psol = Func.Point();
	Rsnld.Root(UVap);
	IntSurf_InteriorPoint intp (psol,UVap(1),UVap(2),
				    Func.Direction3d(),Func.Direction2d());
	list.Append(intp);
      }
    }
  }

  done = Standard_True;
}
Commit	Line	Data
b311480e	1	// Copyright (c) 1995-1999 Matra Datavision
973c2be1	2	// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e	3	//
973c2be1	4	// This file is part of Open CASCADE Technology software library.
b311480e	5	//
d5f74e42	6	// This library is free software; you can redistribute it and/or modify it under
d5f74e42	7	// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1	8	// by the Free Software Foundation, with special exception defined in the file
	9	// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
	10	// distribution for complete text of the license and disclaimer of any warranty.
b311480e	11	//
973c2be1	12	// Alternatively, this file may be used under the terms of Open CASCADE
973c2be1	13	// commercial license or contractual agreement.
7fd59977	14
0797d9d3	15	#ifndef OCCT_DEBUG
7fd59977	16	#define No_Standard_RangeError
	17	#define No_Standard_OutOfRange
	18	#endif
	19
	20
	21
	22	#include <math_FunctionSetRoot.hxx>
	23	#include <Precision.hxx>
	24	#include <gp_Pnt2d.hxx>
	25
	26	#include <TopAbs_State.hxx>
	27
	28	IntStart_SearchInside::IntStart_SearchInside (): done(Standard_False)
	29	{}
	30
	31	IntStart_SearchInside::IntStart_SearchInside (TheFunction& Func,
	32	const ThePSurface& PS,
	33	const Handle(TheTopolTool)& T,
	34	const Standard_Real Epsilon) {
	35
	36	Perform(Func,PS,T,Epsilon);
	37	}
	38
	39
	40	//=======================================================================
	41	//function : Perform
	42	//purpose : Search all inside points
	43	//=======================================================================
	44
	45	void IntStart_SearchInside::Perform (TheFunction& Func,
	46	const ThePSurface& PS,
	47	const Handle(TheTopolTool)& T,
	48	const Standard_Real Epsilon) {
	49
	50	done = Standard_False;
	51	list.Clear();
1ef32e96 RL	52	Standard_Real aBinf[2], aBsup[2], aUVap[2], atoler[2];
1ef32e96 RL	53	math_Vector Binf(aBinf,1,2), Bsup(aBsup,1,2), UVap(aUVap,1,2), toler(atoler,1,2);
7fd59977	54	gp_Pnt psol;
	55	Standard_Boolean testpnt;
	56	Standard_Integer i,j,nbpt;
	57	TopAbs_State situ;
	58	Standard_Real umin,umax,vmin,vmax;
	59	Binf(1) = umin = ThePSurfaceTool::FirstUParameter(PS);
	60	Binf(2) = vmin = ThePSurfaceTool::FirstVParameter(PS);
	61	Bsup(1) = umax = ThePSurfaceTool::LastUParameter(PS);
	62	Bsup(2) = vmax = ThePSurfaceTool::LastVParameter(PS);
	63
	64	Standard_Integer NbsampleU= T->NbSamplesU();
	65	Standard_Integer NbsampleV= T->NbSamplesV();
	66	Standard_Integer Nbsample = T->NbSamples();
	67
	68	Standard_Real du = Bsup(1)-Binf(1);
	69	Standard_Real dv = Bsup(2)-Binf(2);
	70	du/=(Standard_Real)NbsampleU*0.5;
	71	dv/=(Standard_Real)NbsampleV*0.5;
	72
	73	Standard_Real toler1 = toler(1) = ThePSurfaceTool::UResolution(PS,Precision::Confusion());
	74	Standard_Real toler2 = toler(2) = ThePSurfaceTool::VResolution(PS,Precision::Confusion());
	75	Standard_Real Maxtoler1toler2 = toler1;
	76	if(toler2>Maxtoler1toler2) Maxtoler1toler2 = toler2;
	77
	78	//-- lbr le 15 mai 97
	79	//-- on interdit aux points d'etre trop prets des restrictions
	80	Maxtoler1toler2*=1000;
	81	if(Maxtoler1toler2>du0.001) Maxtoler1toler2=du0.001;
	82	if(Maxtoler1toler2>dv0.001) Maxtoler1toler2=dv0.001;
	83
	84
	85	Func.Set(PS);
	86	Standard_Real Tol = Func.Tolerance();
	87
	88	math_FunctionSetRoot Rsnld(Func,toler);
	89
7fd59977	90	//-- lbr le 15 mai 97
	91	umin+=du*0.01;
	92	vmin+=dv*0.01;
	93	umax-=du*0.01;
	94	vmax-=dv*0.01;
	95
	96	//-- lbr le 30 octobre 97 :
	97	//-- Si une surface vient tangenter 2 edges proche d un coin
	98	//-- il faut faire attention qu un point de depart soit trouve au
	99	//-- voisinage du coin. Car ds le cas contraire, le cheminement ne
	100	//-- pourra pas passer au travers des frontieres :
	101	//--
	102	//-- typiquement I est un cylindre (conge)
	103	//--
	104	//-- PPPPPPPPPPPPPPPPPPPP*PPPPPPPPPPPPPPPP
	105	//-- P I I
	106	//-- P I I
	107	//-- P I
	108	//-- P # il faut trouver un point ici
	109	//-- P I
	110	//-- P I
	111	//-- PI
	112	//-- * I
	113	//-- PI I
	114	//-- P I I I I I I I I
	115	//--
	116
	117
	118	for (i=1; i <= Nbsample+12; i++) {
	119	gp_Pnt2d s2d;
	120	gp_Pnt s3d;
	121	Standard_Boolean nepastester=Standard_False;
	122	if(i<=Nbsample) {
	123	T->SamplePoint(i,s2d,s3d);
	124	UVap(1)=s2d.X(); UVap(2)=s2d.Y();
	125
	126	Standard_Real u1,v1,u2,v2;
	127	u1 = Binf(1) = Max(umin,UVap(1)-du);
	128	v1 = Binf(2) = Max(vmin,UVap(2)-dv);
	129	u2 = Bsup(1) = Min(umax,UVap(1)+du);
	130	v2 = Bsup(2) = Min(vmax,UVap(2)+dv);
	131
	132
	133	//-- gp_Pnt Pmilieu = ThePSurfaceTool::Value(PS,0.5(u1+u2),0.5(v1+v2));
	134	gp_Pnt Pextrm1 = ThePSurfaceTool::Value(PS,u1,v1);
	135	gp_Pnt Pextrm2 = ThePSurfaceTool::Value(PS,u2,v2);
1ef32e96 RL	136	Standard_Real aValf[1];
1ef32e96 RL	137	math_Vector Valf(aValf,1,1);
7fd59977	138	Func.Value(UVap,Valf);
	139	Standard_Real rvalf = Valf(1);
	140	Standard_Real DistPP = Pextrm1.SquareDistance(Pextrm2);
	141	if(rvalfrvalf > 3.0DistPP) {
7fd59977	142	nepastester=Standard_True;
	143	}
	144	}
	145	else {
	146	if(i==Nbsample+1) { s2d.SetCoord(umin+du0.02,vmin+dv0.02); }
	147	else if(i==Nbsample+2) { s2d.SetCoord(umax-du0.02,vmin+dv0.02); }
	148	else if(i==Nbsample+3) { s2d.SetCoord(umin+du0.02,vmax-dv0.02); }
	149	else if(i==Nbsample+4) { s2d.SetCoord(umax-du0.02,vmax-dv0.02); }
	150
	151	else if(i==Nbsample+5) { s2d.SetCoord(umin+du0.02,vmin+dv0.02); }
	152	else if(i==Nbsample+6) { s2d.SetCoord(umax-du0.02,vmin+dv0.02); }
	153	else if(i==Nbsample+7) { s2d.SetCoord(umin+du0.02,vmax-dv0.02); }
	154	else if(i==Nbsample+8) { s2d.SetCoord(umax-du0.02,vmax-dv0.02); }
	155
	156	else if(i==Nbsample+9) { s2d.SetCoord(umin+du0.005,vmin+dv0.005); }
	157	else if(i==Nbsample+10){ s2d.SetCoord(umax-du0.005,vmin+dv0.005); }
	158	else if(i==Nbsample+11){ s2d.SetCoord(umin+du0.005,vmax-dv0.005); }
	159	else { s2d.SetCoord(umax-du0.005,vmax-dv0.005); }
	160
	161	UVap(1)=s2d.X(); UVap(2)=s2d.Y();
	162
	163	Binf(1) = Max(umin,UVap(1)-du);
	164	Binf(2) = Max(vmin,UVap(2)-dv);
	165	Bsup(1) = Min(umax,UVap(1)+du);
	166	Bsup(2) = Min(vmax,UVap(2)+dv);
	167	}
	168
	169
	170	if(nepastester==Standard_False) {
7fd59977	171	Rsnld.Perform(Func,UVap,Binf,Bsup);
	172	if (Rsnld.IsDone()) {
	173	if (Abs(Func.Root()) <= Tol) {
	174	if (!Func.IsTangent()) {
	175	psol = Func.Point();
	176	Rsnld.Root(UVap);
	177	// On regarde si le point trouve est bien un nouveau point.
	178	j = 1;
	179	nbpt = list.Length();
	180	testpnt = (j <= nbpt);
	181
	182	while (testpnt) {
	183	const IntSurf_InteriorPoint& IPj = list(j);
	184	const gp_Pnt& Pj = IPj.Value();
	185	if ( (Abs(Pj.X()-psol.X()) <= Epsilon)
	186	&& (Abs(Pj.Y()-psol.Y()) <= Epsilon)
	187	&& (Abs(Pj.Z()-psol.Z()) <= Epsilon)
	188	&& (Abs(UVap(1)-IPj.UParameter()) <= toler1)
	189	&& (Abs(UVap(2)-IPj.VParameter()) <= toler2) ) {
	190	testpnt = Standard_False;
	191	}
	192	else {
	193	j = j+1;
	194	testpnt = (j <= nbpt);
	195	}
	196	}
	197	if (j > nbpt) {
	198	// situ = TheSITool::Classify(PS,UVap(1),UVap(2));
	199	situ = T->Classify(gp_Pnt2d(UVap(1),UVap(2)),
	200	Maxtoler1toler2,Standard_False); //-- ,Standard_False pour ne pas recadrer on Periodic
	201	if (situ == TopAbs_IN) {
	202	list.Append(IntSurf_InteriorPoint(psol,UVap(1),UVap(2),
	203	Func.Direction3d(),
	204	Func.Direction2d()));
	205	}
	206	}
	207	}
	208	}
	209	}
	210	}
	211	}
	212	//-- printf("\n Total : %d Rejet : %d RatioPointCalc : %g nbpt =%d\n",REJET_OK+REJET_KO,REJET_OK,(double)(REJET_KO)/(double)(REJET_OK+REJET_KO),list.Length());
	213	done = Standard_True;
	214	}
	215
	216
	217	//=======================================================================
	218	//function : Perform
	219	//purpose : Test the given inside point
	220	//=======================================================================
	221
	222	void IntStart_SearchInside::Perform (TheFunction& Func,
	223	const ThePSurface& PS,
	224	const Standard_Real UStart,
	225	const Standard_Real VStart)
	226	{
	227	done = Standard_False;
	228	list.Clear();
	229	math_Vector Binf(1,2), Bsup(1,2), toler(1,2);
	230
	231	Binf(1) = ThePSurfaceTool::FirstUParameter(PS);
	232	Binf(2) = ThePSurfaceTool::FirstVParameter(PS);
	233	Bsup(1) = ThePSurfaceTool::LastUParameter(PS);
	234	Bsup(2) = ThePSurfaceTool::LastVParameter(PS);
235
236	toler(1) = ThePSurfaceTool::UResolution(PS,Precision::Confusion());
237	toler(2) = ThePSurfaceTool::VResolution(PS,Precision::Confusion());
238
239	if (UStart-Binf(1) > -toler(1) && UStart-Bsup(1) < toler(1) &&
240	VStart-Binf(2) > -toler(2) && VStart-Bsup(2) < toler(2)) {
241
242	Func.Set(PS);
243	math_Vector UVap(1,2);
244	UVap(1)=UStart; UVap(2)=VStart;
245
246	math_FunctionSetRoot Rsnld(Func,toler);
247	Rsnld.Perform(Func,UVap,Binf,Bsup);
248	if (Rsnld.IsDone()) {
249	Standard_Real tol = Func.Tolerance();
250	Standard_Real valf = Func.Root();
251	if (Abs(valf) <= tol && !Func.IsTangent()) {
252	const gp_Pnt& psol = Func.Point();
253	Rsnld.Root(UVap);
254	IntSurf_InteriorPoint intp (psol,UVap(1),UVap(2),
255	Func.Direction3d(),Func.Direction2d());
256	list.Append(intp);
257	}
258	}
259	}
260
261	done = Standard_True;
262	}