#include <GeometryTest.hxx>
#include <GeomAPI_ProjectPointOnCurve.hxx>
#include <GeomAPI_ProjectPointOnSurf.hxx>
+#include <Extrema_GenLocateExtPS.hxx>
#include <GeomAPI_ExtremaCurveCurve.hxx>
#include <GeomAPI_ExtremaCurveSurface.hxx>
#include <GeomAPI_ExtremaSurfaceSurface.hxx>
#include <TColStd_Array2OfReal.hxx>
#include <Precision.hxx>
#include <stdio.h>
+
#ifdef _WIN32
Standard_IMPORT Draw_Viewer dout;
#endif
//purpose :
//=======================================================================
+static void showProjSolution(Draw_Interpretor& di,
+ const Standard_Integer i,
+ const gp_Pnt& P, const gp_Pnt& P1,
+ const Standard_Real U, const Standard_Real V,
+ const Standard_Boolean isSurface)
+{
+ char name[100];
+ Sprintf(name, "%s%d", "ext_", i);
+ di << name << " ";
+ char* temp = name; // portage WNT
+ if (P.Distance(P1) > Precision::Confusion())
+ {
+ Handle(Geom_Line) L = new Geom_Line(P, gp_Vec(P, P1));
+ Handle(Geom_TrimmedCurve) CT =
+ new Geom_TrimmedCurve(L, 0., P.Distance(P1));
+ DrawTrSurf::Set(temp, CT);
+ }
+ else
+ {
+ DrawTrSurf::Set(temp, P1);
+ if (isSurface)
+ di << " Point on surface ";
+ else
+ di << " Point on curve ";
+ }
+ if (isSurface)
+ di << " Parameters: " << U << " " << V << "\n";
+ else
+ di << " parameter " << i << " = " << U << "\n";
+}
+
+//=======================================================================
+//function : proj
+//purpose :
+//=======================================================================
+
static Standard_Integer proj (Draw_Interpretor& di, Standard_Integer n, const char** a)
{
if ( n < 5)
{
- cout << " Use proj curve/surf x y z [extrema algo: g(grad)/t(tree)]" << endl;
+ std::cout << " Use proj curve/surf x y z [{extrema algo: g(grad)/t(tree)}|{u v}]" << std::endl;
return 1;
}
gp_Pnt P(Draw::Atof(a[2]),Draw::Atof(a[3]),Draw::Atof(a[4]));
- char name[100];
-
Handle(Geom_Curve) GC = DrawTrSurf::GetCurve(a[1]);
Handle(Geom_Surface) GS;
Extrema_ExtAlgo aProjAlgo = Extrema_ExtAlgo_Grad;
if (GS.IsNull())
return 1;
- Standard_Real U1, U2, V1, V2;
- GS->Bounds(U1,U2,V1,V2);
+ if (n <= 6)
+ {
+ Standard_Real U1, U2, V1, V2;
+ GS->Bounds(U1,U2,V1,V2);
- GeomAPI_ProjectPointOnSurf proj(P,GS,U1,U2,V1,V2,aProjAlgo);
+ GeomAPI_ProjectPointOnSurf proj(P,GS,U1,U2,V1,V2,aProjAlgo);
+ if (!proj.IsDone())
+ {
+ di << "projection failed.";
+ return 0;
+ }
- Standard_Real UU,VV;
- for ( Standard_Integer i = 1; i <= proj.NbPoints(); i++)
- {
- gp_Pnt P1 = proj.Point(i);
- if ( P.Distance(P1) > Precision::Confusion())
+ Standard_Real UU,VV;
+ for ( Standard_Integer i = 1; i <= proj.NbPoints(); i++)
{
- Handle(Geom_Line) L = new Geom_Line(P,gp_Vec(P,P1));
- Handle(Geom_TrimmedCurve) CT =
- new Geom_TrimmedCurve(L, 0., P.Distance(P1));
- Sprintf(name,"%s%d","ext_",i);
- char* temp = name; // portage WNT
- DrawTrSurf::Set(temp, CT);
- di << name << " ";
+ gp_Pnt P1 = proj.Point(i);
+ proj.Parameters(i, UU, VV);
+ showProjSolution(di, i, P, P1, UU, VV, Standard_True);
}
- else
+ }
+ else if (n == 7)
+ {
+ const gp_XY aP2d(Draw::Atof(a[5]), Draw::Atof(a[6]));
+ GeomAdaptor_Surface aGAS(GS);
+ Extrema_GenLocateExtPS aProjector(aGAS, Precision::PConfusion(), Precision::PConfusion());
+ aProjector.Perform(P, aP2d.X(), aP2d.Y());
+ if (!aProjector.IsDone())
{
- Sprintf(name,"%s%d","ext_",i);
- di << name << " ";
- char* temp = name; // portage WNT
- DrawTrSurf::Set(temp, P1);
- proj.Parameters(i,UU,VV);
- di << " Le point est sur la surface." << "\n";
- di << " Ses parametres sont: UU = " << UU << "\n";
- di << " VV = " << VV << "\n";
+ di << "projection failed.";
+ return 0;
}
+
+ const Extrema_POnSurf& aP = aProjector.Point();
+ Standard_Real UU, VV;
+ aP.Parameter(UU, VV);
+ showProjSolution(di, 1, P, aP.Value(), UU, VV, Standard_True);
}
}
else
if(proj.NbPoints() == 0)
{
- cout << "No project point was found." << endl;
+ std::cout << "No project point was found." << std::endl;
return 0;
}
{
gp_Pnt P1 = proj.Point(i);
Standard_Real UU = proj.Parameter(i);
- di << " parameter " << i << " = " << UU << "\n";
- if ( P.Distance(P1) > Precision::Confusion())
- {
- Handle(Geom_Line) L = new Geom_Line(P,gp_Vec(P,P1));
- Handle(Geom_TrimmedCurve) CT =
- new Geom_TrimmedCurve(L, 0., P.Distance(P1));
- Sprintf(name,"%s%d","ext_",i);
- char* temp = name; // portage WNT
- DrawTrSurf::Set(temp, CT);
- di << name << " ";
- }
- else
- {
- Sprintf(name,"%s%d","ext_",i);
- char* temp = name; // portage WNT
- DrawTrSurf::Set(temp, P1);
- di << name << " ";
- UU = proj.Parameter(i);
- di << " Le point est sur la courbe." << "\n";
- di << " Son parametre est U = " << UU << "\n";
- }
+ showProjSolution(di, i, P, P1, UU, UU, Standard_False);
}
}
static Standard_Integer surfapp(Draw_Interpretor& di, Standard_Integer n, const char** a)
{
- if ( n < 5 ) return 1;
+ if (n < 5) return 1;
- Standard_Integer i,j;
+ Standard_Integer i, j;
Standard_Integer Nu = Draw::Atoi(a[2]);
Standard_Integer Nv = Draw::Atoi(a[3]);
- TColgp_Array2OfPnt Points (1, Nu, 1, Nv);
+ TColgp_Array2OfPnt Points(1, Nu, 1, Nv);
+ Standard_Boolean IsPeriodic = Standard_False;
+ Standard_Boolean RemoveLast = Standard_False;
- if ( n == 5) {
+ if (n >= 5 && n <= 6) {
Handle(Geom_Surface) Surf = DrawTrSurf::GetSurface(a[4]);
- if ( Surf.IsNull()) return 1;
+ if (Surf.IsNull()) return 1;
Standard_Real U, V, U1, V1, U2, V2;
- Surf->Bounds( U1, U2, V1, V2);
- for ( j = 1; j <= Nv; j++) {
- V = V1 + (j-1) * (V2-V1) / (Nv-1);
- for ( i = 1; i <= Nu; i++) {
- U = U1 + (i-1) * (U2-U1) / (Nu-1);
- Points(i,j) = Surf->Value(U,V);
+ Surf->Bounds(U1, U2, V1, V2);
+ for (j = 1; j <= Nv; j++) {
+ V = V1 + (j - 1) * (V2 - V1) / (Nv - 1);
+ for (i = 1; i <= Nu; i++) {
+ U = U1 + (i - 1) * (U2 - U1) / (Nu - 1);
+ Points(i, j) = Surf->Value(U, V);
+ }
+ }
+ if (n == 6)
+ {
+ Standard_Integer ip = Draw::Atoi(a[5]);
+ if (ip > 0) IsPeriodic = Standard_True;
+ }
+ if (IsPeriodic)
+ {
+ for (j = 1; j <= Nv; j++)
+ {
+ Standard_Real d = Points(1, j).Distance(Points(Nu, j));
+ if (d <= Precision::Confusion())
+ {
+ RemoveLast = Standard_True;
+ break;
+ }
}
- }
+ }
}
- else if ( n >= 16) {
+ else if (n >= 16) {
Standard_Integer Count = 4;
- for ( j = 1; j <= Nv; j++) {
- for ( i = 1; i <= Nu; i++) {
- if ( Count > n) return 1;
- Points(i,j) = gp_Pnt(Draw::Atof(a[Count]),Draw::Atof(a[Count+1]),Draw::Atof(a[Count+2]));
- Count += 3;
+ for (j = 1; j <= Nv; j++) {
+ for (i = 1; i <= Nu; i++) {
+ if (Count > n) return 1;
+ Points(i, j) = gp_Pnt(Draw::Atof(a[Count]), Draw::Atof(a[Count + 1]), Draw::Atof(a[Count + 2]));
+ Count += 3;
}
}
}
char name[100];
Standard_Integer Count = 1;
- for ( j = 1; j <= Nv; j++) {
- for ( i = 1; i <= Nu; i++) {
- Sprintf(name,"point_%d",Count++);
+ for (j = 1; j <= Nv; j++) {
+ for (i = 1; i <= Nu; i++) {
+ Sprintf(name, "point_%d", Count++);
char* temp = name; // portage WNT
- DrawTrSurf::Set(temp,Points(i,j));
+ DrawTrSurf::Set(temp, Points(i, j));
}
- }
+ }
+
+ GeomAPI_PointsToBSplineSurface anApprox;
+ if (RemoveLast)
+ {
+ TColgp_Array2OfPnt Points1(1, Nu - 1, 1, Nv);
+ for (j = 1; j <= Nv; j++)
+ {
+ for (i = 1; i <= Nu - 1; i++) {
+ Points1(i, j) = Points(i, j);
+ }
+ }
+ anApprox.Init(Points1, Approx_ChordLength, 3, 8, GeomAbs_C2, 1.e-3, IsPeriodic);
+ }
+ else
+ {
+ anApprox.Init(Points, Approx_ChordLength, 3, 8, GeomAbs_C2, 1.e-3, IsPeriodic);
+ }
+
+ if (anApprox.IsDone())
+ {
+ Handle(Geom_BSplineSurface) S = anApprox.Surface();
+ DrawTrSurf::Set(a[1], S);
+ di << a[1];
+ }
+
+ return 0;
+}
+
+//=======================================================================
+//function : surfint
+//purpose :
+//=======================================================================
+
+static Standard_Integer surfint(Draw_Interpretor& di, Standard_Integer n, const char** a)
+{
+ if (n < 5) return 1;
+
+ Handle(Geom_Surface) Surf = DrawTrSurf::GetSurface(a[2]);
+ if (Surf.IsNull()) return 1;
+ Standard_Integer i, j;
+ Standard_Integer Nu = Draw::Atoi(a[3]);
+ Standard_Integer Nv = Draw::Atoi(a[4]);
+ TColgp_Array2OfPnt Points(1, Nu, 1, Nv);
+
+ Standard_Real U, V, U1, V1, U2, V2;
+ Surf->Bounds(U1, U2, V1, V2);
+ for (j = 1; j <= Nv; j++) {
+ V = V1 + (j - 1) * (V2 - V1) / (Nv - 1);
+ for (i = 1; i <= Nu; i++) {
+ U = U1 + (i - 1) * (U2 - U1) / (Nu - 1);
+ Points(i, j) = Surf->Value(U, V);
+ }
+ }
+
+ char name[100];
+ Standard_Integer Count = 1;
+ for (j = 1; j <= Nv; j++) {
+ for (i = 1; i <= Nu; i++) {
+ Sprintf(name, "point_%d", Count++);
+ char* temp = name; // portage WNT
+ DrawTrSurf::Set(temp, Points(i, j));
+ }
+ }
+
+ Standard_Boolean IsPeriodic = Standard_False;
+ if (n > 5)
+ {
+ Standard_Integer ip = Draw::Atoi(a[5]);
+ if (ip > 0) IsPeriodic = Standard_True;
+ }
+ Standard_Boolean RemoveLast = Standard_False;
+ if (IsPeriodic)
+ {
+ for (j = 1; j <= Nv; j++)
+ {
+ Standard_Real d = Points(1, j).Distance(Points(Nu, j));
+ if (d <= Precision::Confusion())
+ {
+ RemoveLast = Standard_True;
+ break;
+ }
+ }
+ }
+ const Approx_ParametrizationType ParType = Approx_ChordLength;
+ GeomAPI_PointsToBSplineSurface anApprox;
+ if (RemoveLast)
+ {
+ TColgp_Array2OfPnt Points1(1, Nu-1, 1, Nv);
+ for (j = 1; j <= Nv; j++)
+ {
+ for (i = 1; i <= Nu-1; i++) {
+ Points1(i, j) = Points(i, j);
+ }
+ }
+ anApprox.Interpolate(Points1, ParType, IsPeriodic);
+ }
+ else
+ {
+ anApprox.Interpolate(Points, ParType, IsPeriodic);
+ }
+ if (anApprox.IsDone())
+ {
+ Handle(Geom_BSplineSurface) S = anApprox.Surface();
+ DrawTrSurf::Set(a[1], S);
+ di << a[1];
+ }
+ else
+ {
+ di << "Interpolation not done \n";
+ }
- Handle(Geom_BSplineSurface) S = GeomAPI_PointsToBSplineSurface(Points);
- DrawTrSurf::Set(a[1],S);
- di << a[1];
return 0;
}
return 1;
}
- Handle(Geom_Curve) GC1, GC2;
+ Handle(Geom_Curve) GC1, GC2;
Handle(Geom_Surface) GS1, GS2;
- Standard_Boolean C1 = Standard_False;
- Standard_Boolean C2 = Standard_False;
- Standard_Boolean S1 = Standard_False;
- Standard_Boolean S2 = Standard_False;
+ Standard_Boolean isInfinitySolutions = Standard_False;
+ Standard_Real aMinDist = RealLast();
Standard_Real U1f, U1l, U2f, U2l, V1f = 0., V1l = 0., V2f = 0., V2l = 0.;
GS1 = DrawTrSurf::GetSurface(a[1]);
if ( GS1.IsNull())
return 1;
- S1 = Standard_True;
+
GS1->Bounds(U1f,U1l,V1f,V1l);
}
else {
- C1 = Standard_True;
U1f = GC1->FirstParameter();
U1l = GC1->LastParameter();
}
GS2 = DrawTrSurf::GetSurface(a[2]);
if ( GS2.IsNull())
return 1;
- S2 = Standard_True;
GS2->Bounds(U2f,U2l,V2f,V2l);
}
else {
- C2 = Standard_True;
U2f = GC2->FirstParameter();
U2l = GC2->LastParameter();
}
NCollection_Vector<gp_Pnt> aPnts1, aPnts2;
NCollection_Vector<Standard_Real> aPrms[4];
- if (C1 && C2)
+ if (!GC1.IsNull() && !GC2.IsNull())
{
GeomAPI_ExtremaCurveCurve Ex(GC1, GC2, U1f, U1l, U2f, U2l);
- if (!Ex.Extrema().IsParallel())
+
+ // Since GeomAPI cannot provide access to flag directly.
+ isInfinitySolutions = Ex.Extrema().IsParallel();
+ if (isInfinitySolutions)
+ {
+ aMinDist = Ex.LowerDistance();
+ }
+ else
{
for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
{
aPrms[2].Append(aU2);
}
}
- else
- {
- di << "Infinite number of extremas, distance = " << Ex.LowerDistance() << "\n";
- }
}
- else if (C1 && S2)
+ else if (!GC1.IsNull() && !GS2.IsNull())
{
GeomAPI_ExtremaCurveSurface Ex(GC1, GS2, U1f, U1l, U2f, U2l, V2f, V2l);
- for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+
+ isInfinitySolutions = Ex.Extrema().IsParallel();
+ if (isInfinitySolutions)
{
- gp_Pnt aP1, aP2;
- Ex.Points(aJ, aP1, aP2);
- aPnts1.Append(aP1);
- aPnts2.Append(aP2);
-
- Standard_Real aU1, aU2, aV2;
- Ex.Parameters(aJ, aU1, aU2, aV2);
- aPrms[0].Append(aU1);
- aPrms[2].Append(aU2);
- aPrms[3].Append(aV2);
+ aMinDist = Ex.LowerDistance();
+ }
+ else
+ {
+ for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+ {
+ gp_Pnt aP1, aP2;
+ Ex.Points(aJ, aP1, aP2);
+ aPnts1.Append(aP1);
+ aPnts2.Append(aP2);
+
+ Standard_Real aU1, aU2, aV2;
+ Ex.Parameters(aJ, aU1, aU2, aV2);
+ aPrms[0].Append(aU1);
+ aPrms[2].Append(aU2);
+ aPrms[3].Append(aV2);
+ }
}
}
- else if (S1 && C2)
+ else if (!GS1.IsNull() && !GC2.IsNull())
{
GeomAPI_ExtremaCurveSurface Ex(GC2, GS1, U2f, U2l, U1f, U1l, V1f, V1l);
- for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+
+ isInfinitySolutions = Ex.Extrema().IsParallel();
+ if (isInfinitySolutions)
{
- gp_Pnt aP2, aP1;
- Ex.Points(aJ, aP2, aP1);
- aPnts1.Append(aP1);
- aPnts2.Append(aP2);
-
- Standard_Real aU1, aV1, aU2;
- Ex.Parameters(aJ, aU2, aU1, aV1);
- aPrms[0].Append(aU1);
- aPrms[1].Append(aV1);
- aPrms[2].Append(aU2);
+ aMinDist = Ex.LowerDistance();
+ }
+ else
+ {
+ for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+ {
+ gp_Pnt aP2, aP1;
+ Ex.Points(aJ, aP2, aP1);
+ aPnts1.Append(aP1);
+ aPnts2.Append(aP2);
+
+ Standard_Real aU1, aV1, aU2;
+ Ex.Parameters(aJ, aU2, aU1, aV1);
+ aPrms[0].Append(aU1);
+ aPrms[1].Append(aV1);
+ aPrms[2].Append(aU2);
+ }
}
}
- else if (S1 && S2)
+ else if (!GS1.IsNull() && !GS2.IsNull())
{
- GeomAPI_ExtremaSurfaceSurface Ex(
- GS1, GS2, U1f, U1l, V1f, V1l, U2f, U2l, V2f, V2l);
- for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+ GeomAPI_ExtremaSurfaceSurface Ex(GS1, GS2, U1f, U1l, V1f, V1l, U2f, U2l, V2f, V2l);
+ // Since GeomAPI cannot provide access to flag directly.
+ isInfinitySolutions = Ex.Extrema().IsParallel();
+ if (isInfinitySolutions)
+ {
+ aMinDist = Ex.LowerDistance();
+ }
+ else
{
- gp_Pnt aP1, aP2;
- Ex.Points(aJ, aP1, aP2);
- aPnts1.Append(aP1);
- aPnts2.Append(aP2);
-
- Standard_Real aU1, aV1, aU2, aV2;
- Ex.Parameters(aJ, aU1, aV1, aU2, aV2);
- aPrms[0].Append(aU1);
- aPrms[1].Append(aV1);
- aPrms[2].Append(aU2);
- aPrms[3].Append(aV2);
+ for (Standard_Integer aJ = 1; aJ <= Ex.NbExtrema(); ++aJ)
+ {
+ gp_Pnt aP1, aP2;
+ Ex.Points(aJ, aP1, aP2);
+ aPnts1.Append(aP1);
+ aPnts2.Append(aP2);
+
+ Standard_Real aU1, aV1, aU2, aV2;
+ Ex.Parameters(aJ, aU1, aV1, aU2, aV2);
+ aPrms[0].Append(aU1);
+ aPrms[1].Append(aV1);
+ aPrms[2].Append(aU2);
+ aPrms[3].Append(aV2);
+ }
}
}
// Output points.
const Standard_Integer aPntCount = aPnts1.Size();
- if (aPntCount == 0)
+ if (aPntCount == 0 || isInfinitySolutions)
{
- di << "No solutions!\n";
+ // Infinity solutions flag may be set with 0 number of
+ // solutions in analytic extrema Curve/Curve.
+ if (isInfinitySolutions)
+ di << "Infinite number of extremas, distance = " << aMinDist << "\n";
+ else
+ di << "No solutions!\n";
}
for (Standard_Integer aJ = 1; aJ <= aPntCount; aJ++)
{
di << "Extrema is point : " << P1.X() << " " << P1.Y() << " " << P1.Z() << "\n";
}
else {
- di << "Extrema is segment of line" << "\n";
+ di << "Extrema is segment of line\n";
Handle(Geom_Line) L = new Geom_Line(P1,gp_Vec(P1,P2));
Handle(Geom_TrimmedCurve) CT =
new Geom_TrimmedCurve(L, 0., P1.Distance(P2));
}
else {
- di << "Curves are infinite and parallel" << "\n";
+ di << "Curves are infinite and parallel\n";
}
di << "Minimal distance : " << Ex.TotalLowerDistance() << "\n";
done = Standard_True;
- theCommands.Add("proj", "proj curve/surf x y z [extrema algo: g(grad)/t(tree)]",__FILE__, proj);
+ theCommands.Add("proj", "proj curve/surf x y z [{extrema algo: g(grad)/t(tree)}|{u v}]\n"
+ "\t\tOptional parameters are relevant to surf only.\n"
+ "\t\tIf initial {u v} are given then local extrema is called",__FILE__, proj);
theCommands.Add("appro", "appro result nbpoint [curve]",__FILE__, appro);
theCommands.Add("surfapp","surfapp result nbupoint nbvpoint x y z ....",
__FILE__,
surfapp);
+
+ theCommands.Add("surfint", "surfint result surf nbupoint nbvpoint [uperiodic]",
+ __FILE__,
+ surfint);
+
theCommands.Add("grilapp",
"grilapp result nbupoint nbvpoint X0 dX Y0 dY z11 z12 .. z1nu .... ",
__FILE__,grilapp);