1
0
mirror of https://git.dev.opencascade.org/repos/occt.git synced 2025-08-09 13:22:24 +03:00

0027112: GeomAPI_Interpolate produces wrong result

Support of the approximation of natural boundary condition (f''= 0)
This commit is contained in:
aml
2016-02-07 09:31:58 +03:00
parent fa6d1712fd
commit 94e01d856b
12 changed files with 847 additions and 451 deletions

View File

@@ -41,6 +41,8 @@ Geom_Geometry.hxx
Geom_HSequenceOfBSplineSurface.hxx
Geom_Hyperbola.cxx
Geom_Hyperbola.hxx
Geom_Interpolate.cxx
Geom_Interpolate.hxx
Geom_Line.cxx
Geom_Line.hxx
Geom_OffsetCurve.cxx

View File

@@ -0,0 +1,73 @@
// Copyright (c) 1994-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.
#include <Geom_Interpolate.hxx>
//=======================================================================
//function : computeLagrangeTangent2d
//purpose : compute tangent of 3-point template in middle point
// using Lagrange polynomial.
//=======================================================================
gp_Vec2d Geom_Interpolate::ComputeLagrangeTangent2d(const TColgp_Array1OfPnt2d& thePointsArray,
const TColStd_Array1OfReal& theParametersArray,
const Standard_Integer theIdx)
{
gp_Vec2d aCurTangent(0.0, 0.0);
if (theIdx == thePointsArray.Lower() ||
theIdx == thePointsArray.Upper())
return aCurTangent;
// Compute tangent in the second point.
Standard_Real t0 = theParametersArray(theIdx - 1);
Standard_Real t1 = theParametersArray(theIdx);
Standard_Real t2 = theParametersArray(theIdx + 1);
gp_Pnt2d p0 = thePointsArray(theIdx - 1);
gp_Pnt2d p1 = thePointsArray(theIdx);
gp_Pnt2d p2 = thePointsArray(theIdx + 1);
aCurTangent.SetXY(p0.XY() * (t1-t2)/((t0-t1)*(t0-t2))
+ p1.XY() * (2*t1-t0-t2)/((t1-t0)*(t1-t2))
+ p2.XY() * (t1-t0)/((t2-t0)*(t2-t1)));
return aCurTangent;
}
//=======================================================================
//function : computeLagrangeTangent
//purpose : compute tangent of 3-point template in middle point
// using Lagrange polynomial.
//=======================================================================
gp_Vec Geom_Interpolate::ComputeLagrangeTangent(const TColgp_Array1OfPnt& thePointsArray,
const TColStd_Array1OfReal& theParametersArray,
const Standard_Integer theIdx)
{
gp_Vec aCurTangent(0.0, 0.0, 0.0);
if (theIdx == thePointsArray.Lower() ||
theIdx == thePointsArray.Upper())
return aCurTangent;
// Compute tangent in the second point.
Standard_Real t0 = theParametersArray(theIdx - 1);
Standard_Real t1 = theParametersArray(theIdx);
Standard_Real t2 = theParametersArray(theIdx + 1);
gp_Pnt p0 = thePointsArray(theIdx - 1);
gp_Pnt p1 = thePointsArray(theIdx);
gp_Pnt p2 = thePointsArray(theIdx + 1);
aCurTangent.SetXYZ(p0.XYZ() * (t1-t2)/((t0-t1)*(t0-t2))
+ p1.XYZ() * (2*t1-t0-t2)/((t1-t0)*(t1-t2))
+ p2.XYZ() * (t1-t0)/((t2-t0)*(t2-t1)));
return aCurTangent;
}

View File

@@ -0,0 +1,48 @@
// Copyright (c) 1994-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 _Geom_Interpolate_HeaderFile
#define _Geom_Interpolate_HeaderFile
#include <gp_Vec2d.hxx>
#include <gp_Vec.hxx>
#include <TColgp_Array1OfPnt2d.hxx>
#include <TColgp_Array1OfPnt.hxx>
#include <TColStd_Array1OfReal.hxx>
//! This enumeration represent boundary condition type for interpolation algorithm.
enum GeomInterpolate_BCType
{
GeomInterpolate_CLAMPED,
GeomInterpolate_NATURAL
};
//! This class contains static method which are used in interpolation algorithms:
//! GeomAPI_Interpolate
//! Geom2dAPI_Interpolate
class Geom_Interpolate
{
public:
Standard_EXPORT static gp_Vec2d ComputeLagrangeTangent2d(const TColgp_Array1OfPnt2d &thePointsArray,
const TColStd_Array1OfReal &theParametersArray,
const Standard_Integer theIdx);
Standard_EXPORT static gp_Vec ComputeLagrangeTangent(const TColgp_Array1OfPnt &thePointsArray,
const TColStd_Array1OfReal &theParametersArray,
const Standard_Integer theIdx);
};
#endif

View File

@@ -127,114 +127,150 @@ static void BuildParameters(const Standard_Boolean PeriodicFlag,
ParametersPtr->Value(ii) + distance) ;
}
}
//=======================================================================
//function : BuildPeriodicTangents
//purpose :
//=======================================================================
static void BuildPeriodicTangent(
const TColgp_Array1OfPnt2d& PointsArray,
TColgp_Array1OfVec2d& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray)
static void BuildPeriodicTangent(const TColgp_Array1OfPnt2d& PointsArray,
TColgp_Array1OfVec2d& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray,
const GeomInterpolate_BCType theBCType)
{
Standard_Integer
ii,
degree ;
Standard_Real *point_array,
*parameter_array,
eval_result[2][2] ;
gp_Vec2d a_vector ;
if (PointsArray.Length() < 3) {
Standard_ConstructionError::Raise();
}
if (!TangentFlags.Value(1)) {
degree = 3 ;
if (PointsArray.Length() == 3) {
degree = 2 ;
}
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ;
parameter_array =
(Standard_Real *) &ParametersArray.Value(1) ;
TangentFlags.SetValue(1,Standard_True) ;
PLib::EvalLagrange(ParametersArray.Value(1),
1,
degree,
2,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 2 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(1,a_vector) ;
gp_Vec2d aCurTangent;
Standard_Integer degree = 3;
if (PointsArray.Length() < 3)
Standard_ConstructionError::Raise();
if (TangentFlags.Value(1))
return; // yet computed.
else
TangentFlags.SetValue(1,Standard_True);
if (theBCType == GeomInterpolate_NATURAL)
{
// Compute tangent in second point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent2d(PointsArray, ParametersArray, PointsArray.Lower() + 1);
// Compute tangent in the first point to satisfy
// f''(0) = 0 condition approximation:
// This is approximation due to usage of the tangent approximation in the second point
// instead of the equation below directly in solver.
//
// Formula from the Golovanov book "Geometrical modeling", 2011, pp. 18:
//
// q0 = 1.5 * (p1 - p0) / (t1 - t0) - 0.5 q1
//
// q1=p0*(t1-t2)/((t0-t1)(t0-t2)) + p1 (2t1-t0-t2)/((t1-t0)(t1-t2))+ p2(t1-t0)/((t2-t0)(t2-t1))
//
// Where:
// t - parameter
// p - point
// q - tangent vectors
gp_Vec2d aVec(PointsArray(1), PointsArray(2));
aCurTangent = 1.5 * aVec / (ParametersArray(2) - ParametersArray(1)) - 0.5 * aCurTangent;
}
else if (theBCType == GeomInterpolate_CLAMPED)
{
Standard_Real *point_array, *parameter_array, eval_result[2][2];
if (PointsArray.Length() < 3)
Standard_ConstructionError::Raise();
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower());
parameter_array = (Standard_Real *) &ParametersArray.Value(1);
PLib::EvalLagrange(ParametersArray.Value(1), 1, degree, 2, point_array[0],
parameter_array[0], eval_result[0][0]);
for (Standard_Integer ii = 1 ; ii <= 2 ; ii++)
aCurTangent.SetCoord(ii,eval_result[1][ii-1]);
}
TangentsArray.SetValue(1,aCurTangent);
}
//=======================================================================
//function : BuildTangents
//purpose :
//=======================================================================
static void BuildTangents(const TColgp_Array1OfPnt2d& PointsArray,
TColgp_Array1OfVec2d& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray)
TColgp_Array1OfVec2d& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray,
const GeomInterpolate_BCType theBCType)
{
Standard_Integer ii,
degree ;
Standard_Real *point_array,
*parameter_array,
eval_result[2][2] ;
gp_Vec2d a_vector ;
degree = 3 ;
gp_Vec2d aCurTangent,
aFirstTangent, // Resulting tangents.
aLastTangent; // Resulting tangents.
Standard_Integer ii, degree = 3;
if (PointsArray.Length() == 3)
degree = 2;
if ( PointsArray.Length() < 3) {
Standard_ConstructionError::Raise();
}
if (PointsArray.Length() == 3) {
degree = 2 ;
}
if (!TangentFlags.Value(1)) {
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ;
parameter_array =
(Standard_Real *) &ParametersArray.Value(1) ;
TangentFlags.SetValue(1,Standard_True) ;
PLib::EvalLagrange(ParametersArray.Value(1),
1,
degree,
2,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 2 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(1,a_vector) ;
}
if (! TangentFlags.Value(TangentFlags.Upper())) {
point_array =
(Standard_Real *) &PointsArray.Value(PointsArray.Upper() - degree) ;
TangentFlags.SetValue(TangentFlags.Upper(),Standard_True) ;
parameter_array =
(Standard_Real *)&ParametersArray.Value(ParametersArray.Upper() - degree) ;
PLib::EvalLagrange(ParametersArray.Value(ParametersArray.Upper()),
1,
degree,
2,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 2 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(TangentsArray.Upper(),a_vector) ;
}
}
if (TangentFlags.Value(1) &&
TangentFlags.Value(TangentFlags.Upper()))
{
return;
}
if (theBCType == GeomInterpolate_NATURAL)
{
// The formulas for tangent computation is stored in BuildPeriodicTangent.
if (!TangentFlags.Value(1))
{
// Compute tangent in second point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent2d(PointsArray, ParametersArray, PointsArray.Lower() + 1);
gp_Vec2d aVec(PointsArray(1), PointsArray(2));
aFirstTangent = 1.5 * aVec / (ParametersArray(2) - ParametersArray(1)) - 0.5 * aCurTangent;
} // if (!TangentFlags.Value(1))
if (!TangentFlags.Value(TangentFlags.Upper()))
{
// Compute tangent in last but one point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent2d(PointsArray, ParametersArray, PointsArray.Upper() - 1);
Standard_Integer aLastIdx = PointsArray.Upper();
gp_Vec2d aVec(PointsArray(aLastIdx - 1), PointsArray(aLastIdx));
aLastTangent = 1.5 * aVec / (ParametersArray(aLastIdx) - ParametersArray(aLastIdx - 1)) - 0.5 * aCurTangent;
}
}
else if (theBCType == GeomInterpolate_CLAMPED)
{
Standard_Real *point_array, *parameter_array, eval_result[2][2];
if (!TangentFlags.Value(1))
{
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower());
parameter_array = (Standard_Real *) &ParametersArray.Value(1);
TangentFlags.SetValue(1,Standard_True);
PLib::EvalLagrange(ParametersArray.Value(1), 1, degree, 3, point_array[0], parameter_array[0], eval_result[0][0]) ;
for (ii = 1 ; ii <= 3 ; ii++)
aFirstTangent.SetCoord(ii,eval_result[1][ii-1]);
} // if (!TangentFlags.Value(1))
if (!TangentFlags.Value(TangentFlags.Upper()))
{
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Upper() - degree);
TangentFlags.SetValue(TangentFlags.Upper(),Standard_True);
parameter_array = (Standard_Real *) &ParametersArray.Value(ParametersArray.Upper() - degree);
PLib::EvalLagrange(ParametersArray.Value(ParametersArray.Upper()), 1, degree, 3,point_array[0], parameter_array[0], eval_result[0][0]);
for (ii = 1 ; ii <= 3 ; ii++)
aLastTangent.SetCoord(ii,eval_result[1][ii-1]);
}
}
if (!TangentFlags.Value(1))
{
TangentFlags.SetValue(1, Standard_True);
TangentsArray.SetValue(1, aFirstTangent);
}
if (!TangentFlags.Value(TangentFlags.Upper()))
{
TangentFlags.SetValue(TangentFlags.Upper(), Standard_True);
TangentsArray.SetValue(TangentsArray.Upper(), aLastTangent);
}
}
//=======================================================================
//function : BuildTangents
//purpose : scale the given tangent so that they have the length of
@@ -310,41 +346,28 @@ static void ScaleTangents(const TColgp_Array1OfPnt2d& PointsArray,
//purpose :
//=======================================================================
Geom2dAPI_Interpolate::Geom2dAPI_Interpolate
(const Handle(TColgp_HArray1OfPnt2d)& PointsPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance) :
myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False)
Geom2dAPI_Interpolate::Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d)& PointsPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType)
: myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False),
myBCType(theBCType)
{
Standard_Integer ii ;
Standard_Boolean result =
CheckPoints(PointsPtr->Array1(),
Tolerance) ;
myTangents =
new TColgp_HArray1OfVec2d(myPoints->Lower(),
myPoints->Upper()) ;
myTangentFlags =
new TColStd_HArray1OfBoolean(myPoints->Lower(),
myPoints->Upper()) ;
Standard_Boolean result = CheckPoints(PointsPtr->Array1(), Tolerance);
myTangents = new TColgp_HArray1OfVec2d(myPoints->Lower(), myPoints->Upper());
myTangentFlags = new TColStd_HArray1OfBoolean(myPoints->Lower(), myPoints->Upper());
if (!result) {
Standard_ConstructionError::Raise();
}
BuildParameters(PeriodicFlag,
PointsPtr->Array1(),
myParameters) ;
if (!result)
Standard_ConstructionError::Raise();
for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
myTangentFlags->SetValue(ii,Standard_False) ;
}
BuildParameters(PeriodicFlag, PointsPtr->Array1(), myParameters);
for (Standard_Integer ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++)
myTangentFlags->SetValue(ii,Standard_False);
}
//=======================================================================
@@ -352,51 +375,40 @@ myTangentRequest(Standard_False)
//purpose :
//=======================================================================
Geom2dAPI_Interpolate::Geom2dAPI_Interpolate
(const Handle(TColgp_HArray1OfPnt2d)& PointsPtr,
const Handle(TColStd_HArray1OfReal)& ParametersPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance) :
myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myParameters(ParametersPtr),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False)
Geom2dAPI_Interpolate::Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d) &PointsPtr,
const Handle(TColStd_HArray1OfReal) &ParametersPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType)
: myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myParameters(ParametersPtr),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False),
myBCType(theBCType)
{
Standard_Integer ii ;
Standard_Boolean result =
CheckPoints(PointsPtr->Array1(),
Tolerance) ;
Standard_Boolean result = CheckPoints(PointsPtr->Array1(), Tolerance);
if (PeriodicFlag) {
if ((PointsPtr->Length()) + 1 != ParametersPtr->Length()) {
Standard_ConstructionError::Raise();
}
}
myTangents =
new TColgp_HArray1OfVec2d(myPoints->Lower(),
myPoints->Upper()) ;
myTangentFlags =
new TColStd_HArray1OfBoolean(myPoints->Lower(),
myPoints->Upper()) ;
if (!result) {
Standard_ConstructionError::Raise();
}
result =
CheckParameters(ParametersPtr->Array1()) ;
if (!result) {
Standard_ConstructionError::Raise();
}
for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
myTangentFlags->SetValue(ii,Standard_False) ;
}
if (PeriodicFlag)
{
if ((PointsPtr->Length()) + 1 != ParametersPtr->Length())
{
Standard_ConstructionError::Raise();
}
}
myTangents = new TColgp_HArray1OfVec2d(myPoints->Lower(), myPoints->Upper());
myTangentFlags = new TColStd_HArray1OfBoolean(myPoints->Lower(), myPoints->Upper());
if (!result)
Standard_ConstructionError::Raise();
result = CheckParameters(ParametersPtr->Array1());
if (!result)
Standard_ConstructionError::Raise();
for (Standard_Integer ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++)
myTangentFlags->SetValue(ii,Standard_False);
}
//=======================================================================
//function : Load
@@ -566,14 +578,15 @@ void Geom2dAPI_Interpolate::PerformPeriodic()
mults.SetValue(num_distinct_knots ,half_order) ;
if (num_points >= 3) {
//
// only enter here if there are more than 3 points otherwise
// it means we have already the tangent
//
//
// only enter here if there are more than 3 points otherwise
// it means we have already the tangent
//
BuildPeriodicTangent(myPoints->Array1(),
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1()) ;
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1(),
myBCType);
}
contact_order_array.SetValue(2,1) ;
parameters.SetValue(1,myParameters->Value(1)) ;
@@ -787,9 +800,10 @@ void Geom2dAPI_Interpolate::PerformNonPeriodic()
// if those where not given in advance
//
BuildTangents(myPoints->Array1(),
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1()) ;
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1(),
myBCType);
}
contact_order_array.SetValue(2,1) ;
parameters.SetValue(1,myParameters->Value(1)) ;

View File

@@ -17,6 +17,8 @@
#ifndef _Geom2dAPI_Interpolate_HeaderFile
#define _Geom2dAPI_Interpolate_HeaderFile
#include <Geom_Interpolate.hxx>
#include <Standard.hxx>
#include <Standard_DefineAlloc.hxx>
#include <Standard_Handle.hxx>
@@ -31,8 +33,6 @@
class Geom2d_BSplineCurve;
class StdFail_NotDone;
class Standard_ConstructionError;
class gp_Vec2d;
//! This class is used to interpolate a BsplineCurve
//! passing through an array of points, with a C2
@@ -54,16 +54,25 @@ public:
//! Tolerance is to check if the points are not too close to one an other
//! It is also used to check if the tangent vector is not too small.
//! There should be at least 2 points
//! There should be at least 2 points.
//! if PeriodicFlag is True then the curve will be periodic.
Standard_EXPORT Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d)& Points, const Standard_Boolean PeriodicFlag, const Standard_Real Tolerance);
//! Last parameter sets the boundary condition type.
Standard_EXPORT Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d)& Points,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType = GeomInterpolate_CLAMPED);
//! if PeriodicFlag is True then the curve will be periodic
//! Warning:
//! There should be as many parameters as there are points
//! except if PeriodicFlag is True : then there should be one more
//! parameter to close the curve
Standard_EXPORT Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d)& Points, const Handle(TColStd_HArray1OfReal)& Parameters, const Standard_Boolean PeriodicFlag, const Standard_Real Tolerance);
//! There should be as many parameters as there are points.
//! except if PeriodicFlag is True : then there should be one more.
//! parameter to close the curve.
//! Last parameter sets the boundary condition type.
Standard_EXPORT Geom2dAPI_Interpolate(const Handle(TColgp_HArray1OfPnt2d)& Points,
const Handle(TColStd_HArray1OfReal)& Parameters,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType = GeomInterpolate_CLAMPED);
//! Assigns this constrained BSpline curve to be
//! tangential to vectors InitialTangent and FinalTangent
@@ -100,14 +109,17 @@ Standard_EXPORT operator Handle(Geom2d_BSplineCurve)() const;
//! Note: in this case, the result is given by the function Curve.
Standard_EXPORT Standard_Boolean IsDone() const;
//! Get boundary condition type.
GeomInterpolate_BCType GetBoundaryCondition()
{
return myBCType;
}
protected:
//! Set boundary condition type.
void GetBoundaryCondition(const GeomInterpolate_BCType theBCType)
{
myBCType = theBCType;
}
private:
@@ -128,6 +140,7 @@ private:
Handle(TColStd_HArray1OfReal) myParameters;
Standard_Boolean myPeriodic;
Standard_Boolean myTangentRequest;
GeomInterpolate_BCType myBCType;
};

View File

@@ -94,147 +94,178 @@ static Standard_Boolean CheckParameters(const
//=======================================================================
//function : BuildParameters
//purpose :
//=======================================================================
//=======================================================================
static void BuildParameters(const Standard_Boolean PeriodicFlag,
const TColgp_Array1OfPnt& PointsArray,
Handle(TColStd_HArray1OfReal)& ParametersPtr)
const TColgp_Array1OfPnt& PointsArray,
Handle(TColStd_HArray1OfReal)& ParametersPtr)
{
Standard_Integer ii,
index ;
Standard_Real distance ;
Standard_Integer
num_parameters = PointsArray.Length() ;
if (PeriodicFlag) {
num_parameters += 1 ;
Standard_Integer ii, index;
Standard_Real distance;
Standard_Integer num_parameters = PointsArray.Length();
if (PeriodicFlag)
num_parameters += 1;
ParametersPtr = new TColStd_HArray1OfReal(1, num_parameters);
ParametersPtr->SetValue(1, 0.0);
index = 2;
for (ii = PointsArray.Lower() ; ii < PointsArray.Upper() ; ii++)
{
distance = PointsArray.Value(ii).Distance(PointsArray.Value(ii+1));
ParametersPtr->SetValue(index, ParametersPtr->Value(ii) + distance);
index += 1;
}
ParametersPtr =
new TColStd_HArray1OfReal(1,
num_parameters) ;
ParametersPtr->SetValue(1,0.0e0) ;
index = 2 ;
for (ii = PointsArray.Lower() ; ii < PointsArray.Upper() ; ii++) {
distance =
PointsArray.Value(ii).Distance(PointsArray.Value(ii+1)) ;
ParametersPtr->SetValue(index,
ParametersPtr->Value(ii) + distance) ;
index += 1 ;
}
if (PeriodicFlag) {
distance =
PointsArray.Value(PointsArray.Upper()).
Distance(PointsArray.Value(PointsArray.Lower())) ;
ParametersPtr->SetValue(index,
ParametersPtr->Value(ii) + distance) ;
if (PeriodicFlag)
{
distance = PointsArray.Value(PointsArray.Upper()).Distance(PointsArray.Value(PointsArray.Lower()));
ParametersPtr->SetValue(index, ParametersPtr->Value(ii) + distance);
}
}
//=======================================================================
//function : BuildPeriodicTangents
//purpose :
//purpose : Compute initial condition in periodic case
//=======================================================================
static void BuildPeriodicTangent(
const TColgp_Array1OfPnt& PointsArray,
TColgp_Array1OfVec& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray)
static void BuildPeriodicTangent(const TColgp_Array1OfPnt& PointsArray,
TColgp_Array1OfVec& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray,
GeomInterpolate_BCType theBCType)
{
Standard_Integer
ii,
degree ;
Standard_Real *point_array,
*parameter_array,
eval_result[2][3] ;
gp_Vec a_vector ;
if (PointsArray.Length() < 3) {
Standard_ConstructionError::Raise();
}
if (!TangentFlags.Value(1)) {
degree = 3 ;
if (PointsArray.Length() == 3) {
degree = 2 ;
}
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ;
parameter_array =
(Standard_Real *) &ParametersArray.Value(1) ;
TangentFlags.SetValue(1,Standard_True) ;
PLib::EvalLagrange(ParametersArray.Value(1),
1,
degree,
3,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 3 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(1,a_vector) ;
gp_Vec aCurTangent;
Standard_Integer degree = 3;
if (PointsArray.Length() == 3)
degree = 2;
if (TangentFlags.Value(1))
return; // yet computed.
else
TangentFlags.SetValue(1,Standard_True);
if (theBCType == GeomInterpolate_NATURAL)
{
// Compute tangent in second point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent(PointsArray, ParametersArray, PointsArray.Lower() + 1);
// Compute tangent approximation in the first point to satisfy
// f''(0) = 0 condition approximation:
// This is approximation due to usage of the tangent approximation in the second point
// instead of the equation below directly in solver.
//
// Formula from the Golovanov book "Geometrical modeling", 2011, pp. 18:
//
// q0 = 1.5 * (p1 - p0) / (t1 - t0) - 0.5 q1
//
// q1=p0*(t1-t2)/((t0-t1)(t0-t2)) + p1 (2t1-t0-t2)/((t1-t0)(t1-t2))+ p2(t1-t0)/((t2-t0)(t2-t1))
//
// Where:
// t - parameter
// p - point
// q - tangent vector
gp_Vec aVec(PointsArray(1), PointsArray(2));
aCurTangent = 1.5 * aVec / (ParametersArray(2) - ParametersArray(1)) - 0.5 * aCurTangent;
}
}
else if (theBCType == GeomInterpolate_CLAMPED)
{
Standard_Real *point_array, *parameter_array, eval_result[2][3];
if (PointsArray.Length() < 3)
Standard_ConstructionError::Raise();
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower());
parameter_array = (Standard_Real *) &ParametersArray.Value(1);
PLib::EvalLagrange(ParametersArray.Value(1), 1, degree, 3, point_array[0],
parameter_array[0], eval_result[0][0]);
for (Standard_Integer ii = 1 ; ii <= 3 ; ii++)
aCurTangent.SetCoord(ii,eval_result[1][ii-1]);
}
TangentsArray.SetValue(1,aCurTangent);
}
//=======================================================================
//function : BuildTangents
//purpose :
//=======================================================================
static void BuildTangents(const TColgp_Array1OfPnt& PointsArray,
TColgp_Array1OfVec& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray)
{
Standard_Integer ii,
degree ;
Standard_Real *point_array,
*parameter_array,
eval_result[2][3] ;
gp_Vec a_vector ;
degree = 3 ;
if ( PointsArray.Length() < 3) {
Standard_ConstructionError::Raise();
}
if (PointsArray.Length() == 3) {
degree = 2 ;
}
if (!TangentFlags.Value(1)) {
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower()) ;
parameter_array =
(Standard_Real *) &ParametersArray.Value(1) ;
TangentFlags.SetValue(1,Standard_True) ;
PLib::EvalLagrange(ParametersArray.Value(1),
1,
degree,
3,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 3 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(1,a_vector) ;
}
if (! TangentFlags.Value(TangentFlags.Upper())) {
point_array =
(Standard_Real *) &PointsArray.Value(PointsArray.Upper() - degree) ;
TangentFlags.SetValue(TangentFlags.Upper(),Standard_True) ;
parameter_array =
(Standard_Real *) &ParametersArray.Value(ParametersArray.Upper() - degree) ;
PLib::EvalLagrange(ParametersArray.Value(ParametersArray.Upper()),
1,
degree,
3,
point_array[0],
parameter_array[0],
eval_result[0][0]) ;
for (ii = 1 ; ii <= 3 ; ii++) {
a_vector.SetCoord(ii,eval_result[1][ii-1]) ;
}
TangentsArray.SetValue(TangentsArray.Upper(),a_vector) ;
}
}
static void BuildTangents(const TColgp_Array1OfPnt& PointsArray,
TColgp_Array1OfVec& TangentsArray,
TColStd_Array1OfBoolean& TangentFlags,
const TColStd_Array1OfReal& ParametersArray,
const GeomInterpolate_BCType theBCType)
{
gp_Vec aFirstTangent, // Resulting tangents.
aLastTangent; // Resulting tangents.
Standard_Integer ii, degree = 3;
if (PointsArray.Length() == 3)
degree = 2;
gp_Vec aCurTangent;
if (TangentFlags.Value(1) &&
TangentFlags.Value(TangentFlags.Upper()))
{
return;
}
if (theBCType == GeomInterpolate_NATURAL)
{
// The formulas for tangent computation is stored in BuildPeriodicTangent.
if (!TangentFlags.Value(1))
{
// Compute tangent in second point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent(PointsArray, ParametersArray, PointsArray.Lower() + 1);
gp_Vec aVec(PointsArray(1), PointsArray(2));
aFirstTangent = 1.5 * aVec / (ParametersArray(2) - ParametersArray(1)) - 0.5 * aCurTangent;
} // if (!TangentFlags.Value(1))
if (!TangentFlags.Value(TangentFlags.Upper()))
{
// Compute tangent in last but one point.
aCurTangent = Geom_Interpolate::
ComputeLagrangeTangent(PointsArray, ParametersArray, PointsArray.Upper() - 1);
Standard_Integer aLastIdx = PointsArray.Upper();
gp_Vec aVec(PointsArray(aLastIdx - 1), PointsArray(aLastIdx));
aLastTangent = 1.5 * aVec / (ParametersArray(aLastIdx) - ParametersArray(aLastIdx - 1)) - 0.5 * aCurTangent;
}
}
else if (theBCType == GeomInterpolate_CLAMPED)
{
Standard_Real *point_array, *parameter_array, eval_result[2][3];
if (!TangentFlags.Value(1))
{
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Lower());
parameter_array = (Standard_Real *) &ParametersArray.Value(1);
PLib::EvalLagrange(ParametersArray.Value(1), 1, degree, 3, point_array[0], parameter_array[0], eval_result[0][0]) ;
for (ii = 1 ; ii <= 3 ; ii++)
aFirstTangent.SetCoord(ii,eval_result[1][ii-1]);
} // if (!TangentFlags.Value(1))
if (!TangentFlags.Value(TangentFlags.Upper()))
{
point_array = (Standard_Real *) &PointsArray.Value(PointsArray.Upper() - degree);
parameter_array = (Standard_Real *) &ParametersArray.Value(ParametersArray.Upper() - degree);
PLib::EvalLagrange(ParametersArray.Value(ParametersArray.Upper()), 1, degree, 3,point_array[0], parameter_array[0], eval_result[0][0]);
for (ii = 1 ; ii <= 3 ; ii++)
aLastTangent.SetCoord(ii,eval_result[1][ii-1]);
}
}
if (!TangentFlags.Value(1))
{
TangentFlags.SetValue(1, Standard_True);
TangentsArray.SetValue(1, aFirstTangent);
}
if (!TangentFlags.Value(TangentFlags.Upper()))
{
TangentFlags.SetValue(TangentFlags.Upper(), Standard_True);
TangentsArray.SetValue(TangentsArray.Upper(), aLastTangent);
}
}
//=======================================================================
//function : BuildTangents
//purpose : scale the given tangent so that they have the length of
@@ -310,40 +341,28 @@ static void ScaleTangents(const TColgp_Array1OfPnt& PointsArray,
//purpose :
//=======================================================================
GeomAPI_Interpolate::GeomAPI_Interpolate
(const Handle(TColgp_HArray1OfPnt)& PointsPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance) :
myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False)
GeomAPI_Interpolate::GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& PointsPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType)
: myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False),
myBCType(theBCType)
{
Standard_Integer ii ;
Standard_Boolean result =
CheckPoints(PointsPtr->Array1(),
Tolerance) ;
myTangents =
new TColgp_HArray1OfVec(myPoints->Lower(),
myPoints->Upper()) ;
myTangentFlags =
new TColStd_HArray1OfBoolean(myPoints->Lower(),
myPoints->Upper()) ;
Standard_Boolean result = CheckPoints(PointsPtr->Array1(), Tolerance);
myTangents = new TColgp_HArray1OfVec(myPoints->Lower(), myPoints->Upper());
myTangentFlags = new TColStd_HArray1OfBoolean(myPoints->Lower(), myPoints->Upper());
if (!result) {
Standard_ConstructionError::Raise();
}
BuildParameters(PeriodicFlag,
PointsPtr->Array1(),
myParameters) ;
if (!result)
Standard_ConstructionError::Raise();
for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
myTangentFlags->SetValue(ii,Standard_False) ;
}
BuildParameters(PeriodicFlag, PointsPtr->Array1(), myParameters);
for ( Standard_Integer ii = myPoints->Lower(); ii <= myPoints->Upper(); ii++)
myTangentFlags->SetValue(ii,Standard_False);
}
//=======================================================================
@@ -351,51 +370,39 @@ myTangentRequest(Standard_False)
//purpose :
//=======================================================================
GeomAPI_Interpolate::GeomAPI_Interpolate
(const Handle(TColgp_HArray1OfPnt)& PointsPtr,
const Handle(TColStd_HArray1OfReal)& ParametersPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance) :
myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myParameters(ParametersPtr),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False)
GeomAPI_Interpolate::GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& PointsPtr,
const Handle(TColStd_HArray1OfReal)& ParametersPtr,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType)
: myTolerance(Tolerance),
myPoints(PointsPtr),
myIsDone(Standard_False),
myParameters(ParametersPtr),
myPeriodic(PeriodicFlag),
myTangentRequest(Standard_False),
myBCType(theBCType)
{
Standard_Integer ii ;
Standard_Boolean result =
CheckPoints(PointsPtr->Array1(),
Tolerance) ;
Standard_Boolean result = CheckPoints(PointsPtr->Array1(), Tolerance);
if (PeriodicFlag)
{
if ((PointsPtr->Length()) + 1 != ParametersPtr->Length())
{
Standard_ConstructionError::Raise();
}
}
myTangents = new TColgp_HArray1OfVec(myPoints->Lower(), myPoints->Upper());
myTangentFlags = new TColStd_HArray1OfBoolean(myPoints->Lower(), myPoints->Upper());
if (PeriodicFlag) {
if ((PointsPtr->Length()) + 1 != ParametersPtr->Length()) {
Standard_ConstructionError::Raise();
}
}
myTangents =
new TColgp_HArray1OfVec(myPoints->Lower(),
myPoints->Upper()) ;
myTangentFlags =
new TColStd_HArray1OfBoolean(myPoints->Lower(),
myPoints->Upper()) ;
if (!result) {
Standard_ConstructionError::Raise();
}
result =
CheckParameters(ParametersPtr->Array1()) ;
if (!result) {
Standard_ConstructionError::Raise();
}
for (ii = myPoints->Lower() ; ii <= myPoints->Upper() ; ii++) {
myTangentFlags->SetValue(ii,Standard_False) ;
}
if (!result)
Standard_ConstructionError::Raise();
result = CheckParameters(ParametersPtr->Array1());
if (!result)
Standard_ConstructionError::Raise();
for (Standard_Integer ii = myPoints->Lower(); ii <= myPoints->Upper(); ii++)
myTangentFlags->SetValue(ii,Standard_False);
}
//=======================================================================
//function : Load
@@ -573,14 +580,14 @@ void GeomAPI_Interpolate::PerformPeriodic()
mults.SetValue(num_distinct_knots ,half_order) ;
if (num_points >= 3) {
//
// only enter here if there are more than 3 points otherwise
// it means we have already the tangent
//
// only enter here if there are more than 3 points otherwise
// it means we have already the tangent
BuildPeriodicTangent(myPoints->Array1(),
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1()) ;
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1(),
myBCType);
}
contact_order_array.SetValue(2,1) ;
parameters.SetValue(1,myParameters->Value(1)) ;
@@ -789,14 +796,15 @@ void GeomAPI_Interpolate::PerformNonPeriodic()
// check if the boundary conditions are set
//
if (num_points >= 3) {
//
// cannot build the tangents with degree 3 with only 2 points
// if those where not given in advance
//
//
// cannot build the tangents with degree 3 with only 2 points
// if those where not given in advance
//
BuildTangents(myPoints->Array1(),
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1()) ;
myTangents->ChangeArray1(),
myTangentFlags->ChangeArray1(),
myParameters->Array1(),
myBCType);
}
contact_order_array.SetValue(2,1) ;
parameters.SetValue(1,myParameters->Value(1)) ;

View File

@@ -17,6 +17,8 @@
#ifndef _GeomAPI_Interpolate_HeaderFile
#define _GeomAPI_Interpolate_HeaderFile
#include <Geom_Interpolate.hxx>
#include <Standard.hxx>
#include <Standard_DefineAlloc.hxx>
#include <Standard_Handle.hxx>
@@ -30,8 +32,6 @@
#include <TColgp_Array1OfVec.hxx>
#include <Geom_BSplineCurve.hxx>
class gp_Vec;
//! This class is used to interpolate a BsplineCurve
//! passing through an array of points, with a C2
@@ -94,7 +94,11 @@ public:
//! - conditions relating to the respective
//! number of elements in the parallel tables
//! Points and Parameters are not respected.
Standard_EXPORT GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& Points, const Standard_Boolean PeriodicFlag, const Standard_Real Tolerance);
//! Last parameter sets the boundary condition type.
Standard_EXPORT GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& Points,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType = GeomInterpolate_CLAMPED);
//! Initializes an algorithm for constructing a
//! constrained BSpline curve passing through the points of the table
@@ -134,7 +138,12 @@ public:
//! - conditions relating to the respective
//! number of elements in the parallel tables
//! Points and Parameters are not respected.
Standard_EXPORT GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& Points, const Handle(TColStd_HArray1OfReal)& Parameters, const Standard_Boolean PeriodicFlag, const Standard_Real Tolerance);
//! Last parameter sets the boundary condition type.
Standard_EXPORT GeomAPI_Interpolate(const Handle(TColgp_HArray1OfPnt)& Points,
const Handle(TColStd_HArray1OfReal)& Parameters,
const Standard_Boolean PeriodicFlag,
const Standard_Real Tolerance,
const GeomInterpolate_BCType theBCType = GeomInterpolate_CLAMPED);
//! Assigns this constrained BSpline curve to be
//! tangential to vectors InitialTangent and FinalTangent
@@ -173,14 +182,17 @@ Standard_EXPORT operator Handle(Geom_BSplineCurve)() const;
//! Note: in this case, the result is given by the function Curve.
Standard_EXPORT Standard_Boolean IsDone() const;
//! Get boundary condition type.
GeomInterpolate_BCType GetBoundaryCondition()
{
return myBCType;
}
protected:
//! Set boundary condition type.
void GetBoundaryCondition(const GeomInterpolate_BCType theBCType)
{
myBCType = theBCType;
}
private:
@@ -201,6 +213,7 @@ private:
Handle(TColStd_HArray1OfReal) myParameters;
Standard_Boolean myPeriodic;
Standard_Boolean myTangentRequest;
GeomInterpolate_BCType myBCType; // type of boundary conditions
};

View File

@@ -448,11 +448,22 @@ static Standard_Integer lintang (Draw_Interpretor& di,Standard_Integer n, const
static Standard_Integer interpol (Draw_Interpretor& di,Standard_Integer n, const char** a)
//==================================================================================
{
if (n == 1) {
di <<"give a name to your curve !\n";
if (n == 1)
{
di <<"Interpolate: \n"
<< "interpol resCurveName \n"
<< "interpolation via getting points from the axonometric view \n"
<< "\n"
<< "interpol resCurveName fileWithPoints.txt \n"
<< "interpolate curve using data from file \n"
<< "\n"
<< "interpol [3d/2d] [Natural/Clamped] [x y [z]] \n"
<< "interpolate dataset from input \n";
return 0;
}
if (n == 2) {
if (n == 2)
{
// Read points from axonometric viewer.
Standard_Integer id,XX,YY,b, i, j;
di << "Pick points \n";
dout.Select(id, XX, YY, b);
@@ -463,7 +474,7 @@ static Standard_Integer interpol (Draw_Interpretor& di,Standard_Integer n, const
gp_Pnt2d P2d;
Standard_Boolean newcurve;
if (dout.Is3D(id)) {
if (dout.Is3D(id)){
Handle(Draw_Marker3D) mark;
Handle(TColgp_HArray1OfPnt) Points = new TColgp_HArray1OfPnt(1, 1);
P.SetCoord((Standard_Real)XX/zoom,(Standard_Real)YY/zoom, 0.0);
@@ -584,51 +595,171 @@ static Standard_Integer interpol (Draw_Interpretor& di,Standard_Integer n, const
}
}
else if (n == 3) {
// lecture du fichier.
// nbpoints, 2d ou 3d, puis valeurs.
else if (n == 3)
{
// Read data set from file
const char* nomfic = a[2];
ifstream iFile(nomfic, ios::in);
if (!iFile) return 1;
Standard_Integer nbp, i;
if (!iFile)
return 1;
Standard_Integer nbp;
Standard_Real x, y, z;
iFile >> nbp;
char dimen[3];
iFile >> dimen;
if (!strcmp(dimen,"3d")) {
Handle(TColgp_HArray1OfPnt) Point =
new TColgp_HArray1OfPnt(1, nbp);
for (i = 1; i <= nbp; i++) {
std::string aBC; // boundary condition type.
iFile >> dimen >> aBC;
GeomInterpolate_BCType aBCType = GeomInterpolate_CLAMPED;
if (!aBC.compare("Natural"))
aBCType = GeomInterpolate_NATURAL;
else if (!aBC.compare("Clamped"))
aBCType = GeomInterpolate_CLAMPED;
else
{
di << "Incorrect boundary type, \n"
<< "supported boundary conditions are: \"Natural\" or \"Clamped\" \n";
return 1;
}
if (!strcmp(dimen,"3d"))
{
// 3D case.
Handle(TColgp_HArray1OfPnt) Point = new TColgp_HArray1OfPnt(1, nbp);
for (Standard_Integer i = 1; i <= nbp; i++)
{
iFile >> x >> y >> z;
Point->SetValue(i, gp_Pnt(x, y, z));
}
GeomAPI_Interpolate anInterpolator(Point,
Standard_False,
1.0e-5) ;
anInterpolator.Perform() ;
if (anInterpolator.IsDone()) {
Handle(Geom_BSplineCurve) C =
anInterpolator.Curve();
GeomAPI_Interpolate anInterpolator(Point, Standard_False, 1.0e-5, aBCType);
anInterpolator.Perform();
if (anInterpolator.IsDone())
{
Handle(Geom_BSplineCurve) C = anInterpolator.Curve();
DrawTrSurf::Set(a[1], C);
}
}
else if (!strcmp(dimen,"2d")) {
Handle(TColgp_HArray1OfPnt2d) PointPtr =
new TColgp_HArray1OfPnt2d(1, nbp);
for (i = 1; i <= nbp; i++) {
else if (!strcmp(dimen,"2d"))
{
// 2D case.
Handle(TColgp_HArray1OfPnt2d) PointPtr = new TColgp_HArray1OfPnt2d(1, nbp);
for (Standard_Integer i = 1; i <= nbp; i++)
{
iFile >> x >> y;
PointPtr->SetValue(i, gp_Pnt2d(x, y));
}
Geom2dAPI_Interpolate a2dInterpolator(PointPtr,
Standard_False,
1.0e-5);
a2dInterpolator.Perform() ;
if (a2dInterpolator.IsDone()) {
Handle(Geom2d_BSplineCurve) C = a2dInterpolator.Curve() ;
Geom2dAPI_Interpolate a2dInterpolator(PointPtr, Standard_False, 1.0e-5, aBCType);
a2dInterpolator.Perform();
if (a2dInterpolator.IsDone())
{
Handle(Geom2d_BSplineCurve) C = a2dInterpolator.Curve();
DrawTrSurf::Set(a[1], C);
}
}
else
{
di << "Incorrect dimension type, \n"
<< "supported dimensions are: \"3d\" or \"2d\" \n";
return 1;
}
}
else if (n > 4)
{
// Read points from draw command.
// a[0] - method name.
// a[1] - result curve.
// a[2] - 3d/2d.
// a[3] - Natural / Clamped.
// a[4...] - points.
std::string aDim(a[2]);
Standard_Integer aDimInt = 0;
if (!aDim.compare("3d"))
aDimInt = 3;
else if (!aDim.compare("2d"))
aDimInt = 2;
else
{
di << "Incorrect dimension type, \n"
<< "supported dimensions are: \"3d\" or \"2d\" \n";
return 1;
}
std::string aBC(a[3]);
GeomInterpolate_BCType aBCType = GeomInterpolate_CLAMPED;
if (!aBC.compare("Natural"))
aBCType = GeomInterpolate_NATURAL;
else if (!aBC.compare("Clamped"))
aBCType = GeomInterpolate_CLAMPED;
else
{
di << "Incorrect boundary type, \n"
<< "supported boundary conditions are: \"Natural\" or \"Clamped\" \n";
return 1;
}
// Check for dimension / consistence.
Standard_Integer aMod = (n - 4) % aDimInt;
if (aMod)
{
di << "Incorrect number of input points \n"
<< "Number of points should correspond to the dimension \n";
return 1;
}
// Perform interpolation.
Standard_Integer aPntNb = (n - 4) / aDimInt,
anIdx = 4;
Standard_Real aCoords[3];
if (aDimInt == 2)
{
// Read points.
Handle(TColgp_HArray1OfPnt2d) PointPtr = new TColgp_HArray1OfPnt2d(1, aPntNb);
for (Standard_Integer aPntIdx = 1; aPntIdx <= aPntNb; ++aPntIdx)
{
for (Standard_Integer aDimIdx = 1; aDimIdx <= aDimInt; ++aDimIdx)
aCoords[aDimIdx - 1] = Atof(a[anIdx++]);
PointPtr->SetValue(aPntIdx, gp_Pnt2d(aCoords[0], aCoords[1]));
}
// Compute result.
Geom2dAPI_Interpolate anInterpolator(PointPtr, Standard_False, 1.0e-5, aBCType);
anInterpolator.Perform();
if (anInterpolator.IsDone())
{
Handle(Geom2d_BSplineCurve) aC2D = anInterpolator.Curve();
DrawTrSurf::Set(a[1], aC2D);
}
} // if (aDimInt == 2)
else // 3d case
{
// Read points.
Handle(TColgp_HArray1OfPnt) PointPtr = new TColgp_HArray1OfPnt(1, aPntNb);
for (Standard_Integer aPntIdx = 1; aPntIdx <= aPntNb; ++aPntIdx)
{
for (Standard_Integer aDimIdx = 1; aDimIdx <= aDimInt; ++aDimIdx)
aCoords[aDimIdx - 1] = Atof(a[anIdx++]);
PointPtr->SetValue(aPntIdx, gp_Pnt(aCoords[0], aCoords[1], aCoords[2]));
}
// Compute result.
GeomAPI_Interpolate anInterpolator(PointPtr, Standard_False, 1.0e-5, aBCType);
anInterpolator.Perform();
if (anInterpolator.IsDone())
{
Handle(Geom_BSplineCurve) aC3D = anInterpolator.Curve();
DrawTrSurf::Set(a[1], aC3D);
}
} // else // 3d case
} // else if (n > 4)
return 0;
}
@@ -816,7 +947,7 @@ void GeometryTest::ConstraintCommands(Draw_Interpretor& theCommands)
theCommands.Add("interpol",
"interpol cname [fic]",
"run \"interpol\" without parameters for help",
__FILE__,
interpol, g);
theCommands.Add("tanginterpol",

View File

@@ -401,8 +401,6 @@ void GeomliteTest::API2dCommands(Draw_Interpretor& theCommands)
theCommands.Add("2dapprox", "2dapprox result nbpoint [curve] [[x] y [x] y...]",__FILE__,
appro,g);
theCommands.Add("2dinterpole", "2dinterpole result nbpoint [curve] [[x] y [x] y ...]",__FILE__,
appro,g);
g = "GEOMETRY curves and surfaces analysis";

View File

@@ -0,0 +1,33 @@
puts "========"
puts "OCC27112"
puts "========"
puts ""
##############################################
# GeomAPI_Interpolate produces wrong result
##############################################
interpol result 3d Natural \
0020.2145846565971 0043.0749267405586 -0004.42510603802374 \
0023.9807408024132 0014.6827020376834 -0012.4125287876527 \
0024.3667948231688 0013.8677143193155 -0012.641915904261 \
0024.6513678260243 0013.5109991196997 -0012.7423191658556 \
0024.7858704778464 0013.4143863517737 -0012.7695125913718 \
0024.9455105694487 0013.3596568726249 -0012.7849172391444 \
0025.0753806027073 0013.3631194708948 -0012.7839426245161 \
0025.1971400830829 0013.4054685177828 -0012.7720226826757 \
0025.3773723647216 0013.5376332070131 -0012.7348225311948 \
0025.7223738766752 0014.0220793196683 -0012.5984677382982 \
0026.5185459585299 0016.3036441304239 -0011.9563151411137 \
0028.102183898912 0025.7227430797677 -0009.30579278073889 \
0029.7854153434029 0043.0749267405586 -0004.42510603802373
# Result length check.
mkedge anEdge result
checkprops anEdge -l 63.2944
# Visual check.
donly result
smallview
fit
display result
checkview -screenshot -2d -path ${imagedir}/${test_image}.png

View File

@@ -0,0 +1,37 @@
puts "========"
puts "OCC27112"
puts "========"
puts ""
##############################################
# GeomAPI_Interpolate produces wrong result
##############################################
interpol result 3d Natural \
1.0 1.0 0.0 \
2.0 0.0 0.0 \
3.0 1.0 0.0 \
3.0 8.0 0.0 \
4.0 10.0 0.0\
5.0 9.0 0.0 \
4.0 7.0 0.0 \
0.0 5.0 0.0 \
0.0 3.0 0.0 \
3.0 2.0 0.0 \
6.0 3.0 0.0 \
6.0 5.0 0.0 \
4.0 6.0 0.0 \
3.0 6.0 0.0 \
2.0 5.0 0.0 \
2.0 4.0 0.0 \
3.0 3.0 0.0
# Result length check.
mkedge anEdge result
checkprops anEdge -l 38.925
# Visual check.
donly result
top
fit
display result
checkview -screenshot -2d -path ${imagedir}/${test_image}.png

View File

@@ -0,0 +1,26 @@
puts "========"
puts "OCC27112"
puts "========"
puts ""
##############################################
# GeomAPI_Interpolate produces wrong result
##############################################
interpol result 3d Natural \
1.0 1.0 0.0 \
8.0 -2.0 0.0 \
15.0 -6.0 0.0 \
16.0 -8.0 0.0 \
16.5 -7.0 0.0 \
52.0 -1.0 0.0 \
# Result length check.
mkedge anEdge result
checkprops anEdge -l 59.3204
# Visual check.
donly result
top
fit
display result
checkview -screenshot -2d -path ${imagedir}/${test_image}.png