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

0024988: Wrong result done by projection algorithm

Browse Source 
Wrong border 1.0e-9 jump has deleted. Added periodicity information when projecting to surface.
Period "jump" bug fixes.

AppCont_LeastSquare conversion to non cdl class.
AppCont_Function + AppCont_FunctionTool combined in one class providing the same functionality and converted to non cdl.
Testcase modification.

Test cases for issue CR24988

Fixed incorrect comparison.
This commit is contained in:
aml
2014-12-04 15:04:22 +03:00
committed by bugmaster
parent e8feb725a4
commit 368cdde60e
37 changed files with 1164 additions and 1704 deletions
Download Patch File Download Diff File Expand all files Collapse all files

27
src/AppCont/AppCont.cdl
View File

@@ -29,45 +29,22 @@ package AppCont
---Level : Advanced.
-- All methods of all classes will be advanced.
uses AppParCurves, Geom, math, StdFail, TCollection, TColStd, gp,
TColgp, Standard
is
-------------------------------
--- Algorithms:
-------------------------------
generic class LeastSquare;
imported LeastSquare;
------------------------------------------------------
--- Necessary class for approximation a function f(t):
------------------------------------------------------
deferred class Function;
class FunctionTool;
---------------------------------------------------------
--- Necessary class for approximation a 2d function f(t):
---------------------------------------------------------
deferred class Function2d;
class FunctionTool2d;
class FitFunction instantiates LeastSquare from AppCont
(Function from AppCont, FunctionTool from AppCont);
class FitFunction2d instantiates LeastSquare from AppCont
(Function2d from AppCont, FunctionTool2d from AppCont);
imported Function;
end AppCont;

50
src/AppCont/AppCont_Function.cdl
View File

@@ -1,50 +0,0 @@
-- Created on: 1993-09-01
-- Created by: Laurent PAINNOT
-- Copyright (c) 1993-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.
deferred class Function from AppCont
---Purpose: deferred class describing a continous 3d function f(u)
-- This class must be provided by the user to use the
-- approximation algorithm FittingCurve.
uses Pnt from gp,
Vec from gp
is
Delete(me:out) is virtual;
---C++: alias "Standard_EXPORT virtual ~AppCont_Function(){Delete() ; }"
FirstParameter(me) returns Real
---Purpose: returns the first parameter of the function.
is deferred;
LastParameter(me) returns Real
---Purpose: returns the last parameter of the function.
is deferred;
Value(me; U: Real) returns Pnt
---Purpose: returns the point at parameter <U>.
is deferred;
D1(me; U: Real; P: in out Pnt; V: in out Vec) returns Boolean
---Purpose: returns the point and the derivative values at
-- the parameter <U>.
is deferred;
end Function;

18
src/AppCont/AppCont_Function.cxx
View File

@@ -1,18 +0,0 @@
// 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.
#include <AppCont_Function.ixx>
void AppCont_Function::Delete()
{}

91
src/AppCont/AppCont_Function.hxx Normal file
View File

@@ -0,0 +1,91 @@
// 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 AppCont_Function_HeaderFile
#define AppCont_Function_HeaderFile
#include <gp_Pnt.hxx>
#include <gp_Pnt2d.hxx>
#include <gp_Vec.hxx>
#include <gp_Vec2d.hxx>
#include <NCollection_Array1.hxx>
#include <Standard_Integer.hxx>
//! Class describing a continous 3d and/or function f(u).
//! This class must be provided by the user to use the approximation algorithm FittingCurve.
class AppCont_Function
{
public:
Standard_EXPORT AppCont_Function()
{
myNbPnt = -1;
myNbPnt2d = -1;
}
//! Get number of 3d and 2d points returned by "Value" and "D1" functions.
Standard_EXPORT void GetNumberOfPoints(Standard_Integer& theNbPnt,
Standard_Integer& theNbPnt2d) const
{
theNbPnt = myNbPnt;
theNbPnt2d = myNbPnt2d;
}
//! Get number of 3d points returned by "Value" and "D1" functions.
Standard_EXPORT Standard_Integer GetNbOf3dPoints() const
{
return myNbPnt;
}
//! Get number of 2d points returned by "Value" and "D1" functions.
Standard_EXPORT Standard_Integer GetNbOf2dPoints() const
{
return myNbPnt2d;
}
Standard_EXPORT virtual ~AppCont_Function() {}
//! Returns the first parameter of the function.
Standard_EXPORT virtual Standard_Real FirstParameter() const = 0;
//! Returns the last parameter of the function.
Standard_EXPORT virtual Standard_Real LastParameter() const = 0;
//! Returns the point at parameter <theU>.
Standard_EXPORT virtual Standard_Boolean Value(const Standard_Real theU,
NCollection_Array1<gp_Pnt2d>& thePnt2d,
NCollection_Array1<gp_Pnt>& thePnt) const = 0;
//! Returns the derivative at parameter <theU>.
Standard_EXPORT virtual Standard_Boolean D1(const Standard_Real theU,
NCollection_Array1<gp_Vec2d>& theVec2d,
NCollection_Array1<gp_Vec>& theVec) const = 0;
//! Return information about peridicity in output paramateters space.
//! @param theDimIdx Defines index in output parameters space. 1 <= theDimIdx <= 3 * myNbPnt + 2 * myNbPnt2d.
Standard_EXPORT virtual void PeriodInformation(const Standard_Integer /*theDimIdx*/,
Standard_Boolean& IsPeriodic,
Standard_Real& thePeriod) const
{
IsPeriodic = Standard_False;
thePeriod = 0.0;
};
protected:
Standard_Integer myNbPnt;
Standard_Integer myNbPnt2d;
};
#endif

50
src/AppCont/AppCont_Function2d.cdl
View File

@@ -1,50 +0,0 @@
-- Created on: 1993-09-01
-- Created by: Laurent PAINNOT
-- Copyright (c) 1993-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.
deferred class Function2d from AppCont
---Purpose: deferred class describing a continous 2d function f(u)
-- This class must be provided by the user to use the
-- approximation algorithm FittingCurve2d.
uses Pnt2d from gp,
Vec2d from gp
is
Delete(me:out) is virtual;
---C++: alias "Standard_EXPORT virtual ~AppCont_Function2d(){Delete() ; }"
FirstParameter(me) returns Real
---Purpose: returns the first parameter of the function.
is deferred;
LastParameter(me) returns Real
---Purpose: returns the last parameter of the function.
is deferred;
Value(me; U: Real) returns Pnt2d
---Purpose: returns the point at parameter <U>.
is deferred;
D1(me; U: Real; P: in out Pnt2d; V: in out Vec2d) returns Boolean
---Purpose: returns the point and the derivative values at
-- the parameter <U>.
is deferred;
end Function2d;

18
src/AppCont/AppCont_Function2d.cxx
View File

@@ -1,18 +0,0 @@
// 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.
#include <AppCont_Function2d.ixx>
void AppCont_Function2d::Delete()
{}

86
src/AppCont/AppCont_FunctionTool.cdl
View File

@@ -1,86 +0,0 @@
-- Created on: 1993-09-01
-- Created by: Laurent PAINNOT
-- Copyright (c) 1993-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.
class FunctionTool from AppCont
---Purpose: This class is the inteface between the Function
-- class and the tool asked by LeastSquare.
uses Function from AppCont,
Pnt from gp,
Pnt2d from gp,
Vec from gp,
Vec2d from gp,
Array1OfPnt from TColgp,
Array1OfPnt2d from TColgp,
Array1OfVec from TColgp,
Array1OfVec2d from TColgp
is
FirstParameter(myclass; C: Function from AppCont) returns Real;
---Purpose: returns the first parameter of the Function.
LastParameter(myclass; C: Function from AppCont) returns Real;
---Purpose: returns the last parameter of the Function.
NbP2d(myclass; C: Function from AppCont) returns Integer;
---Purpose: Returns 0.
NbP3d(myclass; C: Function from AppCont) returns Integer;
---Purpose: Returns 1. (the approximation will be done only for one
-- function.
Value(myclass; C: Function from AppCont; U: Real; tabPt: out Array1OfPnt);
---Purpose: <tabP> is an array of only 1 element, the point value at
-- the parameter <U>.
D1(myclass; C: Function from AppCont; U: Real; tabV: out Array1OfVec)
returns Boolean;
---Purpose: <tabV> is an array of only 1 element, the derivative
-- value at the parameter <U>.
----------------------------------------------------------
-- the following methods won t be called by the algorithms
-- but the description must exist in the tool.
----------------------------------------------------------
Value(myclass; C: Function from AppCont;U: Real;
tabPt2d: out Array1OfPnt2d);
Value(myclass; C: Function from AppCont; U: Real;
tabPt: out Array1OfPnt;
tabPt2d: out Array1OfPnt2d);
D1(myclass;C: Function from AppCont;U: Real;
tabV2d: out Array1OfVec2d)
returns Boolean;
D1(myclass; C: Function from AppCont; U: Real;
tabV: out Array1OfVec;
tabV2d: out Array1OfVec2d)
returns Boolean;
end FunctionTool;

109
src/AppCont/AppCont_FunctionTool.cxx
View File

@@ -1,109 +0,0 @@
// 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.
#include <AppCont_FunctionTool.ixx>
#include <AppCont_Function.hxx>
#include <TColgp_Array1OfPnt.hxx>
#include <TColgp_Array1OfVec.hxx>
#include <gp_Pnt.hxx>
#include <gp_Vec.hxx>
Standard_Real AppCont_FunctionTool::FirstParameter
(const AppCont_Function& F)
{
return F.FirstParameter();
}
Standard_Real AppCont_FunctionTool::LastParameter
(const AppCont_Function& F)
{
return F.LastParameter();
}
Standard_Integer AppCont_FunctionTool::NbP2d
(const AppCont_Function&)
{
return (0);
}
Standard_Integer AppCont_FunctionTool::NbP3d
(const AppCont_Function&)
{
return (1);
}
void AppCont_FunctionTool::Value(const AppCont_Function& F,
const Standard_Real U,
TColgp_Array1OfPnt& tabPt)
{
tabPt(tabPt.Lower()) = F.Value(U);
}
Standard_Boolean AppCont_FunctionTool::D1
(const AppCont_Function& F,
const Standard_Real U,
TColgp_Array1OfVec& tabV)
{
gp_Pnt P;
gp_Vec V;
Standard_Boolean Ok = F.D1(U, P, V);
tabV(tabV.Lower()) = V;
return Ok;
}
void AppCont_FunctionTool::Value(const AppCont_Function&,
const Standard_Real,
TColgp_Array1OfPnt2d&)
{
}
void AppCont_FunctionTool::Value(const AppCont_Function&,
const Standard_Real,
TColgp_Array1OfPnt&,
TColgp_Array1OfPnt2d&)
{
}
Standard_Boolean AppCont_FunctionTool::D1
(const AppCont_Function&,
const Standard_Real,
TColgp_Array1OfVec2d&)
{
return (Standard_True);
}
Standard_Boolean AppCont_FunctionTool::D1
(const AppCont_Function&,
const Standard_Real,
TColgp_Array1OfVec&,
TColgp_Array1OfVec2d&)
{
return (Standard_True);
}

86
src/AppCont/AppCont_FunctionTool2d.cdl
View File

@@ -1,86 +0,0 @@
-- Created on: 1993-09-01
-- Created by: Laurent PAINNOT
-- Copyright (c) 1993-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.
class FunctionTool2d from AppCont
---Purpose: This class is the inteface between the Function2d
-- class and the tool asked by LeastSquare.
uses Function2d from AppCont,
Pnt from gp,
Pnt2d from gp,
Vec from gp,
Vec2d from gp,
Array1OfPnt from TColgp,
Array1OfPnt2d from TColgp,
Array1OfVec from TColgp,
Array1OfVec2d from TColgp
is
FirstParameter(myclass; C: Function2d from AppCont) returns Real;
---Purpose: returns the first parameter of the Function.
LastParameter(myclass; C: Function2d from AppCont) returns Real;
---Purpose: returns the last parameter of the Function.
NbP2d(myclass; C: Function2d from AppCont) returns Integer;
---Purpose: Returns 1. (the approximation will be done only for one
-- function.
NbP3d(myclass; C: Function2d from AppCont) returns Integer;
---Purpose: Returns 0.
Value(myclass; C: Function2d from AppCont;
U: Real; tabPt: out Array1OfPnt2d);
---Purpose: <tabP> is an array of only 1 element, the point value at
-- the parameter <U>.
D1(myclass; C: Function2d from AppCont; U: Real; tabV: out Array1OfVec2d)
returns Boolean;
---Purpose: <tabV> is an array of only 1 element, the derivative
-- value at the parameter <U>.
----------------------------------------------------------
-- the following methods won t be called by the algorithms
-- but the description must exist in the tool.
----------------------------------------------------------
Value(myclass; C: Function2d from AppCont;U: Real;
tabPt2d: out Array1OfPnt);
Value(myclass; C: Function2d from AppCont; U: Real;
tabPt: out Array1OfPnt;
tabPt2d: out Array1OfPnt2d);
D1(myclass;C: Function2d from AppCont;U: Real;
tabV2d: out Array1OfVec)
returns Boolean;
D1(myclass; C: Function2d from AppCont; U: Real;
tabV: out Array1OfVec;
tabV2d: out Array1OfVec2d)
returns Boolean;
end FunctionTool2d;

109
src/AppCont/AppCont_FunctionTool2d.cxx
View File

@@ -1,109 +0,0 @@
// 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.
#include <AppCont_FunctionTool2d.ixx>
#include <AppCont_Function2d.hxx>
#include <TColgp_Array1OfPnt2d.hxx>
#include <TColgp_Array1OfVec2d.hxx>
#include <gp_Pnt2d.hxx>
#include <gp_Vec2d.hxx>
Standard_Real AppCont_FunctionTool2d::FirstParameter
(const AppCont_Function2d& F)
{
return F.FirstParameter();
}
Standard_Real AppCont_FunctionTool2d::LastParameter
(const AppCont_Function2d& F)
{
return F.LastParameter();
}
Standard_Integer AppCont_FunctionTool2d::NbP2d
(const AppCont_Function2d&)
{
return (1);
}
Standard_Integer AppCont_FunctionTool2d::NbP3d
(const AppCont_Function2d&)
{
return (0);
}
void AppCont_FunctionTool2d::Value(const AppCont_Function2d& F,
const Standard_Real U,
TColgp_Array1OfPnt2d& tabPt)
{
tabPt(tabPt.Lower()) = F.Value(U);
}
Standard_Boolean AppCont_FunctionTool2d::D1
(const AppCont_Function2d& F,
const Standard_Real U,
TColgp_Array1OfVec2d& tabV)
{
gp_Pnt2d P;
gp_Vec2d V;
Standard_Boolean Ok = F.D1(U, P, V);
tabV(tabV.Lower()) = V;
return Ok;
}
void AppCont_FunctionTool2d::Value(const AppCont_Function2d&,
const Standard_Real,
TColgp_Array1OfPnt&)
{
}
void AppCont_FunctionTool2d::Value(const AppCont_Function2d&,
const Standard_Real,
TColgp_Array1OfPnt&,
TColgp_Array1OfPnt2d&)
{
}
Standard_Boolean AppCont_FunctionTool2d::D1
(const AppCont_Function2d&,
const Standard_Real,
TColgp_Array1OfVec&)
{
return (Standard_False);
}
Standard_Boolean AppCont_FunctionTool2d::D1
(const AppCont_Function2d&,
const Standard_Real,
TColgp_Array1OfVec&,
TColgp_Array1OfVec2d&)
{
return (Standard_False);
}

108
src/AppCont/AppCont_LeastSquare.cdl
View File

@@ -1,108 +0,0 @@
-- Created on: 1993-04-22
-- Created by: Laurent PAINNOT
-- Copyright (c) 1993-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.
generic class LeastSquare from AppCont(MultiLine as any;
LineTool as any)
---as TheToolLine(MultiLine)
---Purpose: Makes an approximation of a continous Line described by
-- the tool TheLineTool.
-- Minimizing the difference between the approximate result
-- Curve and a continous MultiLine
uses Matrix from math,
Vector from math,
Constraint from AppParCurves,
MultiCurve from AppParCurves
raises NotDone from StdFail,
OutOfRange from Standard,
DimensionError from Standard
is
Create(SSP: MultiLine; U0, U1: Real; FirstCons, LastCons: Constraint;
Deg: Integer; NbPoints: Integer = 24)
---Purpose: given a continous MultiLine, this algorithm computes
-- the approximation into Bezier curves.
-- NbPoints points are taken on the initial MultiLine for
-- minimizing the surface between the MultiLine and the
-- Bezier curves doing the approximation.
-- The first point will be the point of parameter U0 with
-- a constraint FirstCons.
-- The last point will be the point of parameter U1 with
-- a constraint LastCons.
returns LeastSquare from AppCont
raises DimensionError from Standard;
IsDone(me)
---Purpose: returns True if all has been correctly done.
returns Boolean
is static;
Value(me: in out)
---Purpose: returns the result of the approximation, i.e. a
-- MultiCurve.
-- An exception is raised if NotDone.
---C++: return const &
returns MultiCurve from AppParCurves
raises NotDone from StdFail
is static;
NbBColumns(me; SSP: MultiLine)
---Purpose: is internally used by the constuctor.
returns Integer
is static protected;
Error(me; F: in out Real; MaxE3d, MaxE2d: in out Real)
---Purpose: F is the sum of the square errors at each of the
-- NbPoints of the MultiLine.
-- MaxE3d is the maximum 3d value of these errors.
-- MaxE2d is the maximum 2d value of these errors.
-- An exception is raised if NotDone.
raises NotDone from StdFail
is static;
fields
Done: Boolean;
SCU: MultiCurve from AppParCurves;
Degre: Integer;
Nbdiscret: Integer;
nbP: Integer;
nbP2d: Integer;
Points: Matrix;
Poles: Matrix;
myParam: Vector;
VB: Matrix;
end LeastSquare from AppCont;

566
src/AppCont/AppCont_LeastSquare.cxx Normal file
View File

@@ -0,0 +1,566 @@
// Created on: 1995-03-14
// Created by: Modelistation
// 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_OutOfRange
#define No_Standard_RangeError
#endif
#include <AppCont_LeastSquare.hxx>
#include <math.hxx>
#include <AppParCurves_MultiPoint.hxx>
#include <AppCont_ContMatrices.hxx>
#include <PLib.hxx>
//=======================================================================
//function : AppCont_LeastSquare
//purpose :
//=======================================================================
void AppCont_LeastSquare::FixSingleBorderPoint(const AppCont_Function& theSSP,
const Standard_Real theU,
const Standard_Real theU0,
const Standard_Real theU1,
NCollection_Array1<gp_Pnt2d>& theFix2d,
NCollection_Array1<gp_Pnt>& theFix)
{
Standard_Real aMaxIter = 15.0;
Standard_Integer j, i2;
NCollection_Array1<gp_Pnt> aTabP(1, Max (myNbP, 1)), aPrevP(1, Max (myNbP, 1));
NCollection_Array1<gp_Pnt2d> aTabP2d(1, Max (myNbP2d, 1)), aPrevP2d(1, Max (myNbP2d, 1));
Standard_Real aMult = ((theU - theU0) > (theU1 - theU)) ? 1.0: -1.0;
Standard_Real aStartParam = (theU0 + theU1) / 2.0,
aCurrParam, aPrevDist = 1.0, aCurrDist = 1.0;
for (Standard_Real anIter = 1.0; anIter < aMaxIter; anIter += 1.0)
{
aCurrParam = aStartParam + aMult * (1 - pow(10, -anIter)) * (theU1 - theU0) / 2.0;
theSSP.Value(aCurrParam, aTabP2d, aTabP);
// from second iteration
if (anIter > 1.5)
{
aCurrDist = 0.0;
i2 = 1;
for (j = 1; j <= myNbP; j++)
{
aCurrDist += aTabP(j).Distance(aPrevP(j));
i2 += 3;
}
for (j = 1; j <= myNbP2d; j++)
{
aCurrDist += aTabP2d(j).Distance(aPrevP2d(j));
i2 += 2;
}
// from the third iteration
if (anIter > 2.5 && aCurrDist / aPrevDist > 10.0)
break;
}
aPrevP = aTabP;
aPrevP2d = aTabP2d;
aPrevDist = aCurrDist;
}
theFix2d = aPrevP2d;
theFix = aPrevP;
}
//=======================================================================
//function : AppCont_LeastSquare
//purpose :
//=======================================================================
AppCont_LeastSquare::AppCont_LeastSquare(const AppCont_Function& SSP,
const Standard_Real U0,
const Standard_Real U1,
const AppParCurves_Constraint FirstCons,
const AppParCurves_Constraint LastCons,
const Standard_Integer Deg,
const Standard_Integer myNbPoints)
: mySCU(Deg+1),
myPoints(1, myNbPoints, 1, 3 * SSP.GetNbOf3dPoints() + 2 * SSP.GetNbOf2dPoints()),
myPoles(1, Deg + 1, 1, 3 * SSP.GetNbOf3dPoints() + 2 * SSP.GetNbOf2dPoints(), 0.0),
myParam(1, myNbPoints),
myVB(1, Deg+1, 1, myNbPoints),
myPerInfo(1, 3 * SSP.GetNbOf3dPoints() + 2 * SSP.GetNbOf2dPoints() )
{
myDone = Standard_False;
myDegre = Deg;
math_Matrix InvM(1, Deg+1, 1, Deg + 1);
Standard_Integer i, j, k, c, i2;
Standard_Integer classe = Deg + 1, cl1 = Deg;
Standard_Real U, dU, Coeff, Coeff2;
Standard_Real IBij, IBPij;
Standard_Integer FirstP = 1, LastP = myNbPoints;
Standard_Integer nbcol = 3 * SSP.GetNbOf3dPoints() + 2 * SSP.GetNbOf2dPoints();
math_Matrix B(1, classe, 1, nbcol, 0.0);
Standard_Integer bdeb = 1, bfin = classe;
AppParCurves_Constraint myFirstC = FirstCons, myLastC = LastCons;
SSP.GetNumberOfPoints(myNbP, myNbP2d);
Standard_Integer i2plus1, i2plus2;
myNbdiscret = myNbPoints;
NCollection_Array1<gp_Pnt> aTabP(1, Max (myNbP, 1));
NCollection_Array1<gp_Pnt2d> aTabP2d(1, Max (myNbP2d, 1));
NCollection_Array1<gp_Vec> aTabV(1, Max (myNbP, 1));
NCollection_Array1<gp_Vec2d> aTabV2d(1, Max (myNbP2d, 1));
for(Standard_Integer aDimIdx = 1; aDimIdx <= myNbP * 3 + myNbP2d * 2; aDimIdx++)
{
SSP.PeriodInformation(aDimIdx,
myPerInfo(aDimIdx).isPeriodic,
myPerInfo(aDimIdx).myPeriod);
}
Standard_Boolean Ok;
if (myFirstC == AppParCurves_TangencyPoint)
{
Ok = SSP.D1(U0, aTabV2d, aTabV);
if (!Ok) myFirstC = AppParCurves_PassPoint;
}
if (myLastC == AppParCurves_TangencyPoint)
{
Ok = SSP.D1(U1, aTabV2d, aTabV);
if (!Ok) myLastC = AppParCurves_PassPoint;
}
// Compute control points params on which approximation will be built.
math_Vector GaussP(1, myNbPoints), GaussW(1, myNbPoints);
math::GaussPoints(myNbPoints, GaussP);
math::GaussWeights(myNbPoints, GaussW);
math_Vector TheWeights(1, myNbPoints), VBParam(1, myNbPoints);
dU = 0.5*(U1-U0);
for (i = FirstP; i <= LastP; i++)
{
U = 0.5 * (U1 + U0) + dU * GaussP(i);
if (i <= (myNbPoints+1)/2)
{
myParam(LastP - i + 1) = U;
VBParam(LastP - i + 1) = 0.5 * (1 + GaussP(i));
TheWeights(LastP - i + 1) = 0.5 * GaussW(i);
}
else
{
VBParam(i - (myNbPoints + 1) / 2) = 0.5*(1 + GaussP(i));
myParam(i - (myNbPoints + 1) / 2) = U;
TheWeights(i - (myNbPoints+ 1) / 2) = 0.5 * GaussW(i);
}
}
// Compute control points.
for (i = FirstP; i <= LastP; i++)
{
U = myParam(i);
SSP.Value(U, aTabP2d, aTabP);
i2 = 1;
for (j = 1; j <= myNbP; j++)
{
(aTabP(j)).Coord(myPoints(i, i2), myPoints(i, i2+1), myPoints(i, i2+2));
i2 += 3;
}
for (j = 1; j <= myNbP2d; j++)
{
(aTabP2d(j)).Coord(myPoints(i, i2), myPoints(i, i2+1));
i2 += 2;
}
}
// Fix possible "period jump".
Standard_Integer aMaxDim = 3 * myNbP + 2 * myNbP2d;
for(Standard_Integer aDimIdx = 1; aDimIdx <= aMaxDim; aDimIdx++)
{
if (myPerInfo(aDimIdx).isPeriodic &&
Abs (myPoints(1, aDimIdx) - myPoints(2, aDimIdx)) > myPerInfo(aDimIdx).myPeriod / 2.01 &&
Abs (myPoints(2, aDimIdx) - myPoints(3, aDimIdx)) < myPerInfo(aDimIdx).myPeriod / 2.01)
{
Standard_Real aPeriodMult = (myPoints(1, aDimIdx) < myPoints(2, aDimIdx)) ? 1.0 : -1.0;
Standard_Real aNewParam = myPoints(1, aDimIdx) + aPeriodMult * myPerInfo(aDimIdx).myPeriod;
myPoints(1, aDimIdx) = aNewParam;
}
}
for (Standard_Integer aPntIdx = 1; aPntIdx < myNbPoints; aPntIdx++)
{
for(Standard_Integer aDimIdx = 1; aDimIdx <= aMaxDim; aDimIdx++)
{
if (myPerInfo(aDimIdx).isPeriodic &&
Abs ( myPoints(aPntIdx, aDimIdx) - myPoints(aPntIdx + 1, aDimIdx) ) > myPerInfo(aDimIdx).myPeriod / 2.01)
{
Standard_Real aPeriodMult = (myPoints(aPntIdx, aDimIdx) > myPoints(aPntIdx + 1, aDimIdx)) ? 1.0 : -1.0;
Standard_Real aNewParam = myPoints(aPntIdx + 1, aDimIdx) + aPeriodMult * myPerInfo(aDimIdx).myPeriod;
myPoints(aPntIdx + 1, aDimIdx) = aNewParam;
}
}
}
VBernstein(classe, myNbPoints, myVB);
// Traitement du second membre:
NCollection_Array1<Standard_Real> tmppoints(1, nbcol);
for (c = 1; c <= classe; c++)
{
tmppoints.Init(0.0);
for (i = 1; i <= myNbPoints; i++)
{
Coeff = TheWeights(i) * myVB(c, i);
for (j = 1; j <= nbcol; j++)
{
tmppoints(j) += myPoints(i, j)*Coeff;
}
}
for (k = 1; k <= nbcol; k++)
{
B(c, k) += tmppoints(k);
}
}
if (myFirstC == AppParCurves_NoConstraint &&
myLastC == AppParCurves_NoConstraint) {
math_Matrix InvM(1, classe, 1, classe);
InvMMatrix(classe, InvM);
// Calcul direct des poles:
for (i = 1; i <= classe; i++) {
for (j = 1; j <= classe; j++) {
IBij = InvM(i, j);
for (k = 1; k <= nbcol; k++) {
myPoles(i, k) += IBij * B(j, k);
}
}
}
}
else
{
math_Matrix M(1, classe, 1, classe);
MMatrix(classe, M);
NCollection_Array1<gp_Pnt2d> aFixP2d(1, Max (myNbP2d, 1));
NCollection_Array1<gp_Pnt> aFixP(1, Max (myNbP, 1));
if (myFirstC == AppParCurves_PassPoint ||
myFirstC == AppParCurves_TangencyPoint)
{
SSP.Value(U0, aTabP2d, aTabP);
FixSingleBorderPoint(SSP, U0, U0, U1, aFixP2d, aFixP);
i2 = 1;
for (k = 1; k<= myNbP; k++)
{
if (aFixP(k).Distance(aTabP(k)) > 0.1)
(aFixP(k)).Coord(myPoles(1, i2), myPoles(1, i2 + 1), myPoles(1, i2 + 2));
else
(aTabP(k)).Coord(myPoles(1, i2), myPoles(1, i2 + 1), myPoles(1, i2 + 2));
i2 += 3;
}
for (k = 1; k<= myNbP2d; k++)
{
if (aFixP2d(k).Distance(aTabP2d(k)) > 0.1)
(aFixP2d(k)).Coord(myPoles(1, i2), myPoles(1, i2 + 1));
else
(aTabP2d(k)).Coord(myPoles(1, i2), myPoles(1, i2 + 1));
i2 += 2;
}
for (Standard_Integer aDimIdx = 1; aDimIdx <= aMaxDim; aDimIdx++)
{
if (myPerInfo(aDimIdx).isPeriodic &&
Abs ( myPoles(1, aDimIdx) - myPoints(1, aDimIdx) ) > myPerInfo(aDimIdx).myPeriod / 2.01 )
{
Standard_Real aMult = myPoles(1, aDimIdx) < myPoints(1, aDimIdx)? 1.0: -1.0;
myPoles(1,aDimIdx) += aMult * myPerInfo(aDimIdx).myPeriod;
}
}
}
if (myLastC == AppParCurves_PassPoint ||
myLastC == AppParCurves_TangencyPoint)
{
SSP.Value(U1, aTabP2d, aTabP);
FixSingleBorderPoint(SSP, U1, U0, U1, aFixP2d, aFixP);
i2 = 1;
for (k = 1; k<= myNbP; k++)
{
if (aFixP(k).Distance(aTabP(k)) > 0.1)
(aFixP(k)).Coord(myPoles(classe, i2), myPoles(classe, i2 + 1), myPoles(classe, i2 + 2));
else
(aTabP(k)).Coord(myPoles(classe, i2), myPoles(classe, i2 + 1), myPoles(classe, i2 + 2));
i2 += 3;
}
for (k = 1; k<= myNbP2d; k++)
{
if (aFixP2d(k).Distance(aTabP2d(k)) > 0.1)
(aFixP2d(k)).Coord(myPoles(classe, i2), myPoles(classe, i2 + 1));
else
(aTabP2d(k)).Coord(myPoles(classe, i2), myPoles(classe, i2 + 1));
i2 += 2;
}
for (Standard_Integer aDimIdx = 1; aDimIdx <= 2; aDimIdx++)
{
if (myPerInfo(aDimIdx).isPeriodic &&
Abs ( myPoles(classe, aDimIdx) - myPoints(myNbPoints, aDimIdx) ) > myPerInfo(aDimIdx).myPeriod / 2.01 )
{
Standard_Real aMult = myPoles(classe, aDimIdx) < myPoints(myNbPoints, aDimIdx)? 1.0: -1.0;
myPoles(classe,aDimIdx) += aMult * myPerInfo(aDimIdx).myPeriod;
}
}
}
if (myFirstC == AppParCurves_PassPoint) {
bdeb = 2;
// mise a jour du second membre:
for (i = 1; i <= classe; i++) {
Coeff = M(i, 1);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= myPoles(1, k)*Coeff;
}
}
}
if (myLastC == AppParCurves_PassPoint) {
bfin = cl1;
for (i = 1; i <= classe; i++) {
Coeff = M(i, classe);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= myPoles(classe, k)*Coeff;
}
}
}
if (myFirstC == AppParCurves_TangencyPoint) {
// On fixe le second pole::
bdeb = 3;
SSP.D1(U0, aTabV2d, aTabV);
i2 = 1;
Coeff = (U1-U0)/myDegre;
for (k = 1; k<= myNbP; k++) {
i2plus1 = i2+1; i2plus2 = i2+2;
myPoles(2, i2) = myPoles(1, i2) + aTabV(k).X()*Coeff;
myPoles(2, i2plus1) = myPoles(1, i2plus1) + aTabV(k).Y()*Coeff;
myPoles(2, i2plus2) = myPoles(1, i2plus2) + aTabV(k).Z()*Coeff;
i2 += 3;
}
for (k = 1; k<= myNbP2d; k++) {
i2plus1 = i2+1;
myPoles(2, i2) = myPoles(1, i2) + aTabV2d(k).X()*Coeff;
myPoles(2, i2plus1) = myPoles(1, i2plus1) + aTabV2d(k).Y()*Coeff;
i2 += 2;
}
for (i = 1; i <= classe; i++) {
Coeff = M(i, 1); Coeff2 = M(i, 2);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= myPoles(1, k)*Coeff+myPoles(2, k)*Coeff2;
}
}
}
if (myLastC == AppParCurves_TangencyPoint) {
bfin = classe-2;
SSP.D1(U1, aTabV2d, aTabV);
i2 = 1;
Coeff = (U1-U0)/myDegre;
for (k = 1; k<= myNbP; k++) {
i2plus1 = i2+1; i2plus2 = i2+2;
myPoles(cl1,i2) = myPoles(classe, i2) - aTabV(k).X()*Coeff;
myPoles(cl1,i2plus1) = myPoles(classe, i2plus1) - aTabV(k).Y()*Coeff;
myPoles(cl1,i2plus2) = myPoles(classe, i2plus2) - aTabV(k).Z()*Coeff;
i2 += 3;
}
for (k = 1; k<= myNbP2d; k++) {
i2plus1 = i2+1;
myPoles(cl1,i2) = myPoles(classe, i2) - aTabV2d(k).X()*Coeff;
myPoles(cl1,i2plus1) = myPoles(classe, i2plus1) - aTabV2d(k).Y()*Coeff;
i2 += 2;
}
for (i = 1; i <= classe; i++) {
Coeff = M(i, classe); Coeff2 = M(i, cl1);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= myPoles(classe, k)*Coeff + myPoles(cl1, k)*Coeff2;
}
}
}
if (bdeb <= bfin) {
math_Matrix B2(bdeb, bfin, 1, B.UpperCol(), 0.0);
for (i = bdeb; i <= bfin; i++) {
for (j = 1; j <= classe; j++) {
Coeff = M(i, j);
for (k = 1; k <= nbcol; k++) {
B2(i, k) += B(j, k)*Coeff;
}
}
}
// Resolution:
// ===========
math_Matrix IBP(bdeb, bfin, bdeb, bfin);
// dans IBPMatrix at IBTMatrix ne sont stockees que les resultats pour
// une classe inferieure ou egale a 26 (pour l instant du moins.)
if (bdeb == 2 && bfin == classe-1 && classe <= 26) {
IBPMatrix(classe, IBP);
}
else if (bdeb == 3 && bfin == classe-2 && classe <= 26) {
IBTMatrix(classe, IBP);
}
else {
math_Matrix MP(1, classe, bdeb, bfin);
for (i = 1; i <= classe; i++) {
for (j = bdeb; j <= bfin; j++) {
MP(i, j) = M(i, j);
}
}
math_Matrix IBP1(bdeb, bfin, bdeb, bfin);
IBP1 = MP.Transposed()*MP;
IBP = IBP1.Inverse();
}
myDone = Standard_True;
for (i = bdeb; i <= bfin; i++) {
for (j = bdeb; j <= bfin; j++) {
IBPij = IBP(i, j);;
for (k = 1; k<= nbcol; k++) {
myPoles(i, k) += IBPij * B2(j, k);
}
}
}
}
}
}
//=======================================================================
//function : Value
//purpose :
//=======================================================================
const AppParCurves_MultiCurve& AppCont_LeastSquare::Value()
{
Standard_Integer i, j, j2;
gp_Pnt Pt;
gp_Pnt2d Pt2d;
Standard_Integer ideb = 1, ifin = myDegre+1;
// On met le resultat dans les curves correspondantes
for (i = ideb; i <= ifin; i++) {
j2 = 1;
AppParCurves_MultiPoint MPole(myNbP, myNbP2d);
for (j = 1; j <= myNbP; j++) {
Pt.SetCoord(myPoles(i, j2), myPoles(i, j2+1), myPoles(i,j2+2));
MPole.SetPoint(j, Pt);
j2 += 3;
}
for (j = myNbP+1;j <= myNbP+myNbP2d; j++) {
Pt2d.SetCoord(myPoles(i, j2), myPoles(i, j2+1));
MPole.SetPoint2d(j, Pt2d);
j2 += 2;
}
mySCU.SetValue(i, MPole);
}
return mySCU;
}
//=======================================================================
//function : Error
//purpose :
//=======================================================================
void AppCont_LeastSquare::Error(Standard_Real& F,
Standard_Real& MaxE3d,
Standard_Real& MaxE2d) const
{
Standard_Integer i, j, k, c, i2, classe = myDegre + 1;
Standard_Real Coeff, err3d = 0.0, err2d = 0.0;
Standard_Integer ncol = myPoints.UpperCol() - myPoints.LowerCol() + 1;
math_Matrix MyPoints(1, myNbdiscret, 1, ncol);
MyPoints = myPoints;
MaxE3d = MaxE2d = F = 0.0;
NCollection_Array1<Standard_Real> tmppoles(1, ncol);
for (c = 1; c <= classe; c++)
{
for (k = 1; k <= ncol; k++)
{
tmppoles(k) = myPoles(c, k);
}
for (i = 1; i <= myNbdiscret; i++)
{
Coeff = myVB(c, i);
for (j = 1; j <= ncol; j++)
{
MyPoints(i, j) -= tmppoles(j) * Coeff;
}
}
}
Standard_Real e1, e2, e3;
for (i = 1; i <= myNbdiscret; i++)
{
i2 = 1;
for (j = 1; j<= myNbP; j++) {
e1 = MyPoints(i, i2);
e2 = MyPoints(i, i2+1);
e3 = MyPoints(i, i2+2);
err3d = e1*e1+e2*e2+e3*e3;
MaxE3d = Max(MaxE3d, err3d);
F += err3d;
i2 += 3;
}
for (j = 1; j<= myNbP2d; j++) {
e1 = MyPoints(i, i2);
e2 = MyPoints(i, i2+1);
err2d = e1*e1+e2*e2;
MaxE2d = Max(MaxE2d, err2d);
F += err2d;
i2 += 2;
}
}
MaxE3d = Sqrt(MaxE3d);
MaxE2d = Sqrt(MaxE2d);
}
//=======================================================================
//function : IsDone
//purpose :
//=======================================================================
Standard_Boolean AppCont_LeastSquare::IsDone() const
{
return myDone;
}

504
src/AppCont/AppCont_LeastSquare.gxx
View File

@@ -1,504 +0,0 @@
// Created on: 1995-03-14
// Created by: Modelistation
// 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_OutOfRange
#define No_Standard_RangeError
#endif
#include <math.hxx>
#include <math_Vector.hxx>
#include <math_Matrix.hxx>
#include <TColgp_Array1OfPnt.hxx>
#include <TColgp_Array1OfPnt2d.hxx>
#include <gp_Pnt2d.hxx>
#include <gp_Pnt.hxx>
#include <gp_Vec.hxx>
#include <gp_Vec2d.hxx>
#include <TColgp_Array1OfVec.hxx>
#include <TColgp_Array1OfVec2d.hxx>
#include <AppParCurves_MultiPoint.hxx>
#include <AppCont_ContMatrices.hxx>
#include <PLib.hxx>
//=======================================================================
//function : AppCont_LeastSquare
//purpose :
//=======================================================================
AppCont_LeastSquare::AppCont_LeastSquare
(const MultiLine& SSP,
const Standard_Real U0,
const Standard_Real U1,
const AppParCurves_Constraint FirstCons,
const AppParCurves_Constraint LastCons,
const Standard_Integer Deg,
const Standard_Integer NbPoints):
SCU(Deg+1),
Points(1, NbPoints, 1, NbBColumns(SSP)),
Poles(1, Deg+1, 1, NbBColumns(SSP), 0.0),
myParam(1, NbPoints),
VB(1, Deg+1, 1, NbPoints)
{
Done = Standard_False;
Degre = Deg;
math_Matrix InvM(1, Deg+1, 1, Deg+1);
Standard_Integer i, j, k, c, i2;
Standard_Integer classe = Deg+1, cl1 = Deg;
Standard_Real U, dU, Coeff, Coeff2;
Standard_Real IBij, IBPij;
Standard_Integer FirstP = 1, LastP = NbPoints;
Standard_Integer nbcol = NbBColumns(SSP);
math_Matrix B(1, classe, 1, nbcol, 0.0);
Standard_Integer bdeb = 1, bfin = classe;
AppParCurves_Constraint myFirstC = FirstCons, myLastC = LastCons;
nbP = LineTool::NbP3d(SSP);
nbP2d = LineTool::NbP2d(SSP);
Standard_Integer mynbP = nbP, mynbP2d = nbP2d;
if (nbP == 0) mynbP = 1;
if (nbP2d == 0) mynbP2d = 1;
Standard_Integer i2plus1, i2plus2;
Nbdiscret = NbPoints;
TColgp_Array1OfPnt TabP(1, mynbP);
TColgp_Array1OfVec TabV(1, mynbP);
TColgp_Array1OfPnt2d TabP2d(1, mynbP2d);
TColgp_Array1OfVec2d TabV2d(1, mynbP2d);
Standard_Boolean Ok;
if (myFirstC == AppParCurves_TangencyPoint) {
if (nbP != 0 && nbP2d != 0) Ok=LineTool::D1(SSP, U0, TabV, TabV2d);
else if (nbP != 0) Ok=LineTool::D1(SSP, U0, TabV);
else Ok=LineTool::D1(SSP, U0, TabV2d);
if (!Ok) myFirstC = AppParCurves_PassPoint;
}
if (myLastC == AppParCurves_TangencyPoint) {
if (nbP != 0 && nbP2d != 0) Ok=LineTool::D1(SSP, U1, TabV, TabV2d);
else if (nbP != 0) Ok=LineTool::D1(SSP, U1, TabV);
else Ok=LineTool::D1(SSP, U1, TabV2d);
if (!Ok) myLastC = AppParCurves_PassPoint;
}
math_Vector GaussP(1, NbPoints), GaussW(1, NbPoints);
math::GaussPoints(NbPoints, GaussP);
math::GaussWeights(NbPoints, GaussW);
math_Vector TheWeights(1, NbPoints), VBParam(1, NbPoints);
dU = 0.5*(U1-U0);
// calcul et mise en ordre des parametres et des poids:
for (i = FirstP; i <= LastP; i++) {
U = 0.5*(U1+U0) + dU*GaussP(i);
if (i <= (NbPoints+1)/2) {
myParam(LastP-i+1) = U;
VBParam(LastP-i+1) = 0.5*(1 + GaussP(i));
TheWeights(LastP-i+1) = 0.5*GaussW(i);
}
else {
VBParam(i-(NbPoints+1)/2) = 0.5*(1 + GaussP(i));
myParam(i-(NbPoints+1)/2) = U;
TheWeights(i-(NbPoints+1)/2) = 0.5*GaussW(i);
}
}
for (i = FirstP; i <= LastP; i++) {
U = myParam(i);
if (nbP != 0 && nbP2d != 0) LineTool::Value(SSP, U, TabP, TabP2d);
else if (nbP != 0) LineTool::Value(SSP, U, TabP);
else LineTool::Value(SSP, U, TabP2d);
i2 = 1;
for (j = 1; j <= nbP; j++) {
(TabP(j)).Coord(Points(i, i2), Points(i, i2+1), Points(i, i2+2));
i2 += 3;
}
for (j = 1; j <= nbP2d; j++) {
(TabP2d(j)).Coord(Points(i, i2), Points(i, i2+1));
i2 += 2;
}
}
// Calcul du VB ( Valeur des fonctions de Bernstein):
// for (i = 1; i <= classe; i++) {
// for (j = 1; j <= NbPoints; j++) {
// VB(i,j)=PLib::Binomial(cl1,i-1)*Pow((1-VBParam(j)),classe-i)*
// Pow(VBParam(j),i-1);
// }
// }
VBernstein(classe, NbPoints, VB);
// Traitement du second membre:
Standard_Real *tmppoints, *tmpbis;
tmppoints = new Standard_Real[nbcol];
for (c = 1; c <= classe; c++) {
tmpbis = tmppoints;
for (k = 1; k <= nbcol; k++, tmpbis++) {
*tmpbis = 0.0;
}
for (i = 1; i <= NbPoints; i++) {
Coeff = TheWeights(i)*VB(c, i);
tmpbis = tmppoints;
for (j = 1; j <= nbcol; j++, tmpbis++) {
*tmpbis += Points(i, j)*Coeff;
//B(c, j) += Points(i, j)*Coeff;
}
}
tmpbis = tmppoints;
for (k = 1; k <= nbcol; k++, tmpbis++) {
B(c, k) += *tmpbis;
}
}
delete [] tmppoints;
if (myFirstC == AppParCurves_NoConstraint &&
myLastC == AppParCurves_NoConstraint) {
math_Matrix InvM(1, classe, 1, classe);
InvMMatrix(classe, InvM);
// Calcul direct des poles:
for (i = 1; i <= classe; i++) {
for (j = 1; j <= classe; j++) {
IBij = InvM(i, j);
for (k = 1; k <= nbcol; k++) {
Poles(i, k) += IBij * B(j, k);
}
}
}
}
else {
math_Matrix M(1, classe, 1, classe);
MMatrix(classe, M);
if (myFirstC == AppParCurves_PassPoint ||
myFirstC == AppParCurves_TangencyPoint) {
if (nbP != 0 && nbP2d != 0) LineTool::Value(SSP, U0, TabP, TabP2d);
else if (nbP != 0) LineTool::Value(SSP, U0, TabP);
else LineTool::Value(SSP, U0, TabP2d);
i2 =1;
for (k = 1; k<= nbP; k++) {
(TabP(k)).Coord(Poles(1, i2), Poles(1, i2+1), Poles(1, i2+2));
i2 += 3;
}
for (k = 1; k<= nbP2d; k++) {
(TabP2d(k)).Coord(Poles(1, i2), Poles(1, i2+1));
i2 += 2;
}
}
if (myLastC == AppParCurves_PassPoint ||
myLastC == AppParCurves_TangencyPoint) {
i2 = 1;
if (nbP != 0 && nbP2d != 0) LineTool::Value(SSP, U1, TabP, TabP2d);
else if (nbP != 0) LineTool::Value(SSP, U1, TabP);
else LineTool::Value(SSP, U1, TabP2d);
for (k = 1; k<= nbP; k++) {
(TabP(k)).Coord(Poles(classe,i2),
Poles(classe,i2+1),
Poles(classe,i2+2));
i2 += 3;
}
for (k = 1; k<= nbP2d; k++) {
(TabP2d(k)).Coord(Poles(classe, i2), Poles(classe, i2+1));
i2 += 2;
}
}
if (myFirstC == AppParCurves_PassPoint) {
bdeb = 2;
// mise a jour du second membre:
for (i = 1; i <= classe; i++) {
Coeff = M(i, 1);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= Poles(1, k)*Coeff;
}
}
}
if (myLastC == AppParCurves_PassPoint) {
bfin = cl1;
for (i = 1; i <= classe; i++) {
Coeff = M(i, classe);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= Poles(classe, k)*Coeff;
}
}
}
if (myFirstC == AppParCurves_TangencyPoint) {
// On fixe le second pole::
bdeb = 3;
if (nbP != 0 && nbP2d != 0) LineTool::D1(SSP, U0, TabV, TabV2d);
else if (nbP != 0) LineTool::D1(SSP, U0, TabV);
else LineTool::D1(SSP, U0, TabV2d);
i2 = 1;
Coeff = (U1-U0)/Degre;
for (k = 1; k<= nbP; k++) {
i2plus1 = i2+1; i2plus2 = i2+2;
Poles(2, i2) = Poles(1, i2) + TabV(k).X()*Coeff;
Poles(2, i2plus1) = Poles(1, i2plus1) + TabV(k).Y()*Coeff;
Poles(2, i2plus2) = Poles(1, i2plus2) + TabV(k).Z()*Coeff;
i2 += 3;
}
for (k = 1; k<= nbP2d; k++) {
i2plus1 = i2+1;
Poles(2, i2) = Poles(1, i2) + TabV2d(k).X()*Coeff;
Poles(2, i2plus1) = Poles(1, i2plus1) + TabV2d(k).Y()*Coeff;
i2 += 2;
}
for (i = 1; i <= classe; i++) {
Coeff = M(i, 1); Coeff2 = M(i, 2);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= Poles(1, k)*Coeff+Poles(2, k)*Coeff2;
}
}
}
if (myLastC == AppParCurves_TangencyPoint) {
bfin = classe-2;
if (nbP != 0 && nbP2d != 0) LineTool::D1(SSP, U1, TabV, TabV2d);
else if (nbP != 0) LineTool::D1(SSP, U1, TabV);
else LineTool::D1(SSP, U1, TabV2d);
i2 = 1;
Coeff = (U1-U0)/Degre;
for (k = 1; k<= nbP; k++) {
i2plus1 = i2+1; i2plus2 = i2+2;
Poles(cl1,i2) = Poles(classe, i2) - TabV(k).X()*Coeff;
Poles(cl1,i2plus1) = Poles(classe, i2plus1) - TabV(k).Y()*Coeff;
Poles(cl1,i2plus2) = Poles(classe, i2plus2) - TabV(k).Z()*Coeff;
i2 += 3;
}
for (k = 1; k<= nbP2d; k++) {
i2plus1 = i2+1;
Poles(cl1,i2) = Poles(classe, i2) - TabV2d(k).X()*Coeff;
Poles(cl1,i2plus1) = Poles(classe, i2plus1) - TabV2d(k).Y()*Coeff;
i2 += 2;
}
for (i = 1; i <= classe; i++) {
Coeff = M(i, classe); Coeff2 = M(i, cl1);
for (k = 1; k <= nbcol; k++) {
B(i, k) -= Poles(classe, k)*Coeff + Poles(cl1, k)*Coeff2;
}
}
}
if (bdeb <= bfin) {
math_Matrix B2(bdeb, bfin, 1, B.UpperCol(), 0.0);
for (i = bdeb; i <= bfin; i++) {
for (j = 1; j <= classe; j++) {
Coeff = M(i, j);
for (k = 1; k <= nbcol; k++) {
B2(i, k) += B(j, k)*Coeff;
}
}
}
// Resolution:
// ===========
math_Matrix IBP(bdeb, bfin, bdeb, bfin);
// dans IBPMatrix at IBTMatrix ne sont stockees que les resultats pour
// une classe inferieure ou egale a 26 (pour l instant du moins.)
if (bdeb == 2 && bfin == classe-1 && classe <= 26) {
IBPMatrix(classe, IBP);
}
else if (bdeb == 3 && bfin == classe-2 && classe <= 26) {
IBTMatrix(classe, IBP);
}
else {
math_Matrix MP(1, classe, bdeb, bfin);
for (i = 1; i <= classe; i++) {
for (j = bdeb; j <= bfin; j++) {
MP(i, j) = M(i, j);
}
}
math_Matrix IBP1(bdeb, bfin, bdeb, bfin);
IBP1 = MP.Transposed()*MP;
IBP = IBP1.Inverse();
}
Done = Standard_True;
for (i = bdeb; i <= bfin; i++) {
for (j = bdeb; j <= bfin; j++) {
IBPij = IBP(i, j);;
for (k = 1; k<= nbcol; k++) {
Poles(i, k) += IBPij * B2(j, k);
}
}
}
}
}
}
//=======================================================================
//function : NbBColumns
//purpose :
//=======================================================================
Standard_Integer AppCont_LeastSquare::NbBColumns(
const MultiLine& SSP) const
{
Standard_Integer BCol;
BCol = (LineTool::NbP3d(SSP))*3 +
(LineTool::NbP2d(SSP))*2;
return BCol;
}
//=======================================================================
//function : Value
//purpose :
//=======================================================================
const AppParCurves_MultiCurve& AppCont_LeastSquare::Value()
{
Standard_Integer i, j, j2;
gp_Pnt Pt;
gp_Pnt2d Pt2d;
Standard_Integer ideb = 1, ifin = Degre+1;
// On met le resultat dans les curves correspondantes
for (i = ideb; i <= ifin; i++) {
j2 = 1;
AppParCurves_MultiPoint MPole(nbP, nbP2d);
for (j = 1; j <= nbP; j++) {
Pt.SetCoord(Poles(i, j2), Poles(i, j2+1), Poles(i,j2+2));
MPole.SetPoint(j, Pt);
j2 += 3;
}
for (j = nbP+1;j <= nbP+nbP2d; j++) {
Pt2d.SetCoord(Poles(i, j2), Poles(i, j2+1));
MPole.SetPoint2d(j, Pt2d);
j2 += 2;
}
SCU.SetValue(i, MPole);
}
return SCU;
}
//=======================================================================
//function : Error
//purpose :
//=======================================================================
void AppCont_LeastSquare::Error(Standard_Real& F,
Standard_Real& MaxE3d,
Standard_Real& MaxE2d) const
{
Standard_Integer i, j, k, c, i2, classe = Degre+1;
// Standard_Real Coeff, val = 0.0, err3d = 0.0, err2d =0.0;
Standard_Real Coeff, err3d = 0.0, err2d =0.0;
Standard_Integer ncol = Points.UpperCol()-Points.LowerCol()+1;
math_Matrix MyPoints(1, Nbdiscret, 1, ncol);
MyPoints = Points;
MaxE3d = MaxE2d = F = 0.0;
Standard_Real *tmppoles, *tmpbis;
tmppoles = new Standard_Real[ncol];
for (c = 1; c <= classe; c++) {
tmpbis = tmppoles;
for (k = 1; k <= ncol; k++, tmpbis++) {
*tmpbis = Poles(c, k);
}
for (i = 1; i <= Nbdiscret; i++) {
Coeff = VB(c, i);
tmpbis = tmppoles;
for (j = 1; j <= ncol; j++, tmpbis++) {
MyPoints(i, j) -= (*tmpbis)*Coeff; // Poles(c, j)*Coeff;
}
}
}
delete [] tmppoles;
Standard_Real e1, e2, e3;
for (i = 1; i <= Nbdiscret; i++) {
i2 = 1;
for (j = 1; j<= nbP; j++) {
e1 = MyPoints(i, i2);
e2 = MyPoints(i, i2+1);
e3 = MyPoints(i, i2+2);
err3d = e1*e1+e2*e2+e3*e3;
MaxE3d = Max(MaxE3d, err3d);
F += err3d;
i2 += 3;
}
for (j = 1; j<= nbP2d; j++) {
e1 = MyPoints(i, i2);
e2 = MyPoints(i, i2+1);
err2d = e1*e1+e2*e2;
MaxE2d = Max(MaxE2d, err2d);
F += err2d;
i2 += 2;
}
}
MaxE3d = Sqrt(MaxE3d);
MaxE2d = Sqrt(MaxE2d);
}
//=======================================================================
//function : IsDone
//purpose :
//=======================================================================
Standard_Boolean AppCont_LeastSquare::IsDone() const
{
return Done;
}

75
src/AppCont/AppCont_LeastSquare.hxx Normal file
View File

@@ -0,0 +1,75 @@
// Created on: 1995-03-14
// Created by: Modelistation
// 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 AppCont_LeastSquare_HeaderFile
#define AppCont_LeastSquare_HeaderFile
#include <AppCont_Function.hxx>
#include <AppParCurves_MultiCurve.hxx>
#include <math_Vector.hxx>
#include <math_Matrix.hxx>
#include <NCollection_Array1.hxx>
#include <AppParCurves_Constraint.hxx>
struct PeriodicityInfo
{
Standard_Boolean isPeriodic;
Standard_Real myPeriod;
};
class AppCont_LeastSquare
{
public:
Standard_EXPORT AppCont_LeastSquare(const AppCont_Function& SSP,
const Standard_Real U0,
const Standard_Real U1,
const AppParCurves_Constraint FirstCons,
const AppParCurves_Constraint LastCons,
const Standard_Integer Deg,
const Standard_Integer NbPoints);
Standard_EXPORT const AppParCurves_MultiCurve& Value();
Standard_EXPORT void Error(Standard_Real& F,
Standard_Real& MaxE3d,
Standard_Real& MaxE2d) const;
Standard_EXPORT Standard_Boolean IsDone() const;
private:
//! Fix border point evaluation.
void FixSingleBorderPoint(const AppCont_Function& theSSP,
const Standard_Real theU,
const Standard_Real theU0,
const Standard_Real theU1,
NCollection_Array1<gp_Pnt2d>& theFix2d,
NCollection_Array1<gp_Pnt>& theFix);
AppParCurves_MultiCurve mySCU;
math_Matrix myPoints;
math_Matrix myPoles;
math_Vector myParam;
math_Matrix myVB;
NCollection_Array1<PeriodicityInfo> myPerInfo;
Standard_Boolean myDone;
Standard_Integer myDegre;
Standard_Integer myNbdiscret, myNbP, myNbP2d;
};
#endif

3
src/AppCont/FILES
View File

@@ -4,3 +4,6 @@ AppCont_ContMatrices_1.cxx
AppCont_ContMatrices_2.cxx
AppCont_ContMatrices_3.cxx
AppCont_ContMatrices_4.cxx
AppCont_Function.hxx
AppCont_LeastSquare.hxx
AppCont_LeastSquare.cxx