[occt.git] / src / BRepFill / BRepFill_LocationLaw.cxx

// Created on: 1998-01-14
// Created by: Philippe MANGIN
// Copyright (c) 1998-1999 Matra Datavision
// Copyright (c) 1999-2012 OPEN CASCADE SAS
//
// The content of this file is subject to the Open CASCADE Technology Public
// License Version 6.5 (the "License"). You may not use the content of this file
// except in compliance with the License. Please obtain a copy of the License
// at http://www.opencascade.org and read it completely before using this file.
//
// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
//
// The Original Code and all software distributed under the License is
// distributed on an "AS IS" basis, without warranty of any kind, and the
// Initial Developer hereby disclaims all such warranties, including without
// limitation, any warranties of merchantability, fitness for a particular
// purpose or non-infringement. Please see the License for the specific terms
// and conditions governing the rights and limitations under the License.


#include <BRepFill_LocationLaw.ixx>

#include <BRepTools_WireExplorer.hxx>
#include <BRep_Tool.hxx>
#include <BRep_Builder.hxx>
#include <BRepAdaptor_Curve.hxx>
#include <Adaptor3d_HCurve.hxx>
#include <TopoDS.hxx>
#include <TopoDS_Edge.hxx>
#include <TopExp.hxx>
#include <TopLoc_Location.hxx>

#include <GeomFill_LocationLaw.hxx>
#include <gp_Vec.hxx>
#include <gp_Mat.hxx>
#include <gp_XYZ.hxx>
#include <gp_Trsf.hxx>
#include <GCPnts_AbscissaPoint.hxx>
#include <TColgp_Array1OfPnt2d.hxx>
#include <TColgp_Array1OfVec2d.hxx>
#include <TColStd_SequenceOfInteger.hxx>
#include <Precision.hxx>


//=======================================================================
//function : Norm
//purpose  : Norm of a Matrix
//=======================================================================

static Standard_Real Norm(const gp_Mat& M) {
  Standard_Real R, Norme;
  gp_XYZ Coord;
  Coord = M.Row(1);
  Norme = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
  Coord = M.Row(2);
  R = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
  if (R>Norme) Norme = R;
  Coord = M.Row(3);
  R = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
  if (R>Norme) Norme = R;

  return Norme;
}

//=======================================================================
//function : ToG0
//purpose  : Calculate tranformation T such as T.M2 = M1
//=======================================================================

static void ToG0(const gp_Mat& M1, const gp_Mat& M2, gp_Mat& T) {
  T =  M2.Inverted();
  T *= M1;
}
			  
//=======================================================================
//function : BRepFill_LocationLaw
//purpose  : 
//=======================================================================

void BRepFill_LocationLaw::Init(const TopoDS_Wire& Path)
 
{
  Standard_Integer NbEdge;
  BRepTools_WireExplorer wexp;
// Class BRep_Tool without fields and without Constructor :
//  BRep_Tool B;
  TopoDS_Edge E;

  myPath = Path;
  myTol = 1.e-4;

  for (NbEdge=0, wexp.Init(myPath); 
       wexp.More(); wexp.Next()) 
//    if (! B.Degenerated(wexp.Current())) NbEdge++;
    if (! BRep_Tool::Degenerated(wexp.Current())) NbEdge++;
  

  myLaws = new (GeomFill_HArray1OfLocationLaw)(1, NbEdge);
  myLength = new (TColStd_HArray1OfReal) (1, NbEdge+1);
  myLength->Init(-1.);
  myLength->SetValue(1, 0.);
  myEdges = new (TopTools_HArray1OfShape) (1, NbEdge);
  myDisc.Nullify();
  TangentIsMain();
}

//=======================================================================
//function : GetStatus
//purpose  : 
//=======================================================================
 GeomFill_PipeError BRepFill_LocationLaw::GetStatus() const
{
  Standard_Integer ii, N = myLaws->Length();
  GeomFill_PipeError Status =  GeomFill_PipeOk;
  for (ii=1; ii<=N && (Status == GeomFill_PipeOk); ii++) {
    Status = myLaws->Value(ii)->ErrorStatus();
  }
  return  Status;
}

//=======================================================================
//function : TangentIsMain
//purpose  : 
//=======================================================================
void BRepFill_LocationLaw::TangentIsMain() 
{
  myType = 1;
}

//=======================================================================
//function : NormalIsMain
//purpose  : 
//=======================================================================
void BRepFill_LocationLaw::NormalIsMain() 
{
  myType = 2;
}

//=======================================================================
//function : BiNormalIsMain
//purpose  : 
//=======================================================================
void BRepFill_LocationLaw::BiNormalIsMain() 
{
  myType = 3;
}

//=======================================================================
//function : TransformInCompatibleLaw
//purpose  : Set in continuity of laws
//=======================================================================
 void BRepFill_LocationLaw::TransformInCompatibleLaw(const Standard_Real TolAngular)
{

  Standard_Real First, Last, Angle;
  Standard_Integer ipath;
  gp_Mat Trsf, M1, M2;
  gp_Vec V, T1, T2, N1, N2;
  gp_XYZ OZ(0, 0, 1);

  myLaws->Value(1)->GetDomain(First, Last);

  for (ipath=2; ipath<=myLaws->Length(); ipath++) {
    myLaws->Value(ipath-1)->D0(Last, M1, V);
    myLaws->Value(ipath)->GetDomain(First, Last);
    myLaws->Value(ipath)->D0(First, M2, V);
    T1.SetXYZ(M1.Column(3));
    T2.SetXYZ(M2.Column(3));
    N1.SetXYZ(M1.Column(1));
    N2.SetXYZ(M2.Column(1));
    if (T1.IsParallel(T2, TolAngular ) && 
	!T1.IsOpposite(T2, TolAngular)) { // Correction G0
      ToG0(M1, M2, Trsf);
    }
    else {
      Standard_Real alpha;
      gp_Vec cross(T1);
      cross.Cross(T2);
      alpha = T2.AngleWithRef(T1, cross);
      gp_Ax1 axe(gp::Origin(), cross.XYZ());
      N2.Rotate(axe, alpha); 

#if DEB
      if (N2.Dot(T1) > 1.e-9) {
	cout << "Inprecision in TransformInCompatibleLaw" << endl;
        cout << "--- T1.R(N2) = " << N2.Dot(T1) << endl;
	gp_Vec tt;
	tt = T1;
	tt.Rotate(axe, alpha);
	cout << "--- T1.R(T2) = " << tt.Dot(T1) << endl;
	cout << "--- R(N2).R(T2) = " << N2.Dot(tt) << endl;
      }      
#endif
      Angle = N2.AngleWithRef(N1, T1);
      Trsf.SetRotation(OZ, Angle);
    }
    myLaws->Value(ipath)->SetTrsf(Trsf);
  }
} 

//=======================================================================
//function : TransformInG0Law
//purpose  : Set in continuity of laws
//=======================================================================
 void BRepFill_LocationLaw::TransformInG0Law()
{

  Standard_Real First, Last;
  Standard_Integer ipath;
  gp_Mat  M1, M2, aux;//,Trsf
  gp_Vec V;
  myLaws->Value(1)->GetDomain(First, Last);
  for (ipath=2; ipath<=myLaws->Length(); ipath++) {
    myLaws->Value(ipath-1)->D0(Last, M1, V);
    myLaws->Value(ipath)->GetDomain(First, Last);
    myLaws->Value(ipath)->D0(First, M2, V);
    ToG0(M1, M2, aux);
    myLaws->Value(ipath)->SetTrsf(aux);
  }
  
  // Is the law periodical ?
  if  (myPath.Closed()) {
    myLaws->Value(myLaws->Length())->D0(Last, M1, V);
    myLaws->Value(1)->GetDomain(First, Last);
    myLaws->Value(1)->D0(First, M2, V);
  }
}

//=======================================================================
//function : DeleteTransform
//purpose  : Remove the setting in continuity of law. 
//=======================================================================
 void BRepFill_LocationLaw::DeleteTransform()
{
  gp_Mat Id;
  Id.SetIdentity();
  for (Standard_Integer ii=1; ii<=myEdges->Length(); ii++) {
    myLaws->ChangeValue(ii)->SetTrsf(Id);
  }
  myDisc.Nullify();
}

//=======================================================================
//function : NbHoles
//purpose  : Find "Holes"
//=======================================================================
 Standard_Integer BRepFill_LocationLaw::NbHoles(const Standard_Real Tol) 
{
  if (myDisc.IsNull()) {
    TColStd_SequenceOfInteger Seq;
    Standard_Integer ii, NbDisc;
    for (ii=2, NbDisc=-1; ii<=myLaws->Length()+1; ii++) {
      if (IsG1(ii-1, Tol, 1.e-12) == -1) {
	Seq.Append(ii);
      }
    }
    NbDisc = Seq.Length();
    if ( NbDisc > 0) {
      myDisc = new (TColStd_HArray1OfInteger)(1, NbDisc);
      for (ii=1; ii<=NbDisc; ii++)
	 myDisc->SetValue(ii, Seq(ii));
    }   
  }
  if (myDisc.IsNull()) return 0;
  return myDisc->Length();

}

//=======================================================================
//function : Holes
//purpose  : 
//=======================================================================
 void BRepFill_LocationLaw::Holes(TColStd_Array1OfInteger& Disc) const
{
  if (!myDisc.IsNull()) {
    for (Standard_Integer ii=1; ii<=myDisc->Length(); ii++)
      Disc(ii) = myDisc->Value(ii);
  }
}

//=======================================================================
//function : NbLaw
//purpose  : 
//=======================================================================
 Standard_Integer BRepFill_LocationLaw::NbLaw() const
{
  return myLaws->Length();
}

//=======================================================================
//function : Law
//purpose  : 
//=======================================================================
const Handle(GeomFill_LocationLaw)& 
BRepFill_LocationLaw::Law(const Standard_Integer Index) const
{
  return myLaws->Value(Index);
}

//=======================================================================
//function : Wire
//purpose  : 
//=======================================================================
const TopoDS_Wire& BRepFill_LocationLaw::Wire() const
{
  return myPath;
}

//=======================================================================
//function : Edge
//purpose  : 
//=======================================================================
const TopoDS_Edge& BRepFill_LocationLaw::Edge(const Standard_Integer Index) const
{
  return TopoDS::Edge(myEdges->Value(Index));
}

//=======================================================================
//function : Vertex
//purpose  : 
//=======================================================================
 TopoDS_Vertex BRepFill_LocationLaw::Vertex(const Standard_Integer Index) const
{
  TopoDS_Edge E;
  TopoDS_Vertex V;
  if (Index <= myEdges->Length()) {
    E = TopoDS::Edge(myEdges->Value(Index));
    if (E.Orientation() == TopAbs_REVERSED) 
       V = TopExp::LastVertex(E);
    else V = TopExp::FirstVertex(E);
  }
  else if (Index == myEdges->Length()+1) {
    E = TopoDS::Edge(myEdges->Value(Index-1));
    if (E.Orientation() == TopAbs_REVERSED) 
      V = TopExp::FirstVertex(E);
    else V = TopExp::LastVertex(E);
  }
  return V;
}

//===================================================================
//function : PerformVertex
//purpose  : Calculate a vertex of sweeping from a vertex of section
//           and the index of the edge in the trajectory
//===================================================================
void BRepFill_LocationLaw::PerformVertex(const Standard_Integer Index,
					 const TopoDS_Vertex& Input,
					 const Standard_Real TolMin,
					 TopoDS_Vertex& Output,
					 const Standard_Integer ILoc) const
{
  BRep_Builder B;
  Standard_Boolean IsBary = (ILoc == 0);
  Standard_Real First, Last;
  gp_Pnt P;
  gp_Vec V1, V2;//, V;
  gp_Mat M1, M2;

  if (Index>0 && Index<myLaws->Length()) {
    if (ILoc <=0) {
      myLaws->Value(Index)->GetDomain(First, Last);
      myLaws->Value(Index)->D0(Last, M1, V1);
    }

    if (ILoc >= 0) {
      myLaws->Value(Index+1)->GetDomain(First, Last);
      if (ILoc == 0)
	myLaws->Value(Index+1)->D0(First, M2, V2);
      else
	myLaws->Value(Index+1)->D0(First, M1, V1);
    }
  }

  if (Index == 0 || Index == myLaws->Length()) {
    if (!myPath.Closed() || (IsG1(Index, TolMin) != 1)) {
      IsBary = Standard_False;
      if (Index == 0) {
	myLaws->Value(1)->GetDomain(First, Last);
	myLaws->Value(1)->D0(First, M1, V1);
      }
      else {
	myLaws->Value(myLaws->Length())->GetDomain(First, Last);
	myLaws->Value(myLaws->Length())->D0(Last, M1, V1);
      } 
    }
    else {
      if (ILoc <=0) {
	myLaws->Value(myLaws->Length())->GetDomain(First, Last);
	myLaws->Value(myLaws->Length())->D0(Last, M1, V1);
      }

      if (ILoc >=0) {      
	myLaws->Value(1)->GetDomain(First, Last);
	if (ILoc==0)
	  myLaws->Value(1)->D0(First, M2, V2);
	else
	  myLaws->Value(1)->D0(First, M1, V1);
      }
    }
  }

  P = BRep_Tool::Pnt(Input);

  if (IsBary) { 
    gp_XYZ P1(P.XYZ()), P2(P.XYZ());
    P1 *= M1;
    P1 += V1.XYZ();
    P2 *= M2;
    P2 += V2.XYZ();
  
    P.ChangeCoord().SetLinearForm(0.5, P1, 0.5, P2);
    P1 -= P2;
    Standard_Real Tol =  P1.Modulus()/2;
    Tol += TolMin;
    B.MakeVertex(Output, P, Tol);
  }
  else {
    P.ChangeCoord() *= M1;
    P.ChangeCoord() += V1.XYZ();
    B.MakeVertex(Output, P, TolMin);    
  }
  
}

//=======================================================================
//function : CurvilinearBounds
//purpose  : 
//=======================================================================
void BRepFill_LocationLaw::CurvilinearBounds(const Standard_Integer Index,
					     Standard_Real& First, 
					     Standard_Real& Last) const
{
  First = myLength->Value(Index);
  Last  = myLength->Value(Index+1);
  if (Last<0) { //It is required to carry out the calculation 
    Standard_Integer ii, NbE = myEdges->Length();
    Standard_Real Length, f, l;
    GCPnts_AbscissaPoint AbsC;

    for (ii=1, Length=0.; ii<=NbE; ii++) {
      myLaws->Value(ii)->GetDomain(f, l);
      Length += AbsC.Length(myLaws->Value(ii)->GetCurve()->GetCurve(), myTol);
      myLength->SetValue(ii+1, Length);
    }

    First = myLength->Value(Index);
    Last  = myLength->Value(Index+1); 
  }
}

 Standard_Boolean BRepFill_LocationLaw::IsClosed() const
{
  return myPath.Closed();
}

//=======================================================================
//function : IsG1
//purpose  : Evaluate the continuity of the law by a vertex
//=======================================================================
 Standard_Integer 
 BRepFill_LocationLaw::IsG1(const Standard_Integer Index,
			    const Standard_Real SpatialTolerance,
			    const Standard_Real AngularTolerance) const
{
  gp_Vec V1, DV1, V2, DV2;
  gp_Mat M1, M2, DM1, DM2;
  Standard_Real First, Last, EpsNul = 1.e-12;
  Standard_Real TolEps = SpatialTolerance;
#ifndef DEB
  Standard_Boolean Ok_D1 = Standard_False;
#else
  Standard_Boolean Ok_D1;
#endif
  TopoDS_Vertex V;
  TopoDS_Edge E;
  TColgp_Array1OfPnt2d Bid1 (1,1);
  TColgp_Array1OfVec2d Bid2 (1,1);
  
  if (Index>0 && Index<myLaws->Length()) {
    myLaws->Value(Index)->GetDomain(First, Last);
    Ok_D1 = myLaws->Value(Index)->D1(Last, M1, V1, DM1, DV1, 
				     Bid1, Bid2);
    if (!Ok_D1)  myLaws->Value(Index)->D0(Last, M1, V1);

    myLaws->Value(Index+1)->GetDomain(First, Last);
    if (Ok_D1)
      Ok_D1 = myLaws->Value(Index+1)->D1(First, M2, V2, DM2, DV2, 
					 Bid1, Bid2);
    if (!Ok_D1)  myLaws->Value(Index+1)->D0(First, M2, V2);

    E = TopoDS::Edge(myEdges->Value(Index+1));
  }
  if (Index == 0 || Index == myLaws->Length()) {
    if (!myPath.Closed()) return -1;
    myLaws->Value(myLaws->Length())->GetDomain(First, Last);
    Ok_D1 = myLaws->Value(myLaws->Length())->D1(Last, M1, V1, DM1, DV1, 
						Bid1, Bid2);
    if (!Ok_D1)  myLaws->Value(myLaws->Length())->D0(Last, M1, V1);

    myLaws->Value(1)->GetDomain(First, Last);
    if (Ok_D1)
      myLaws->Value(1)->D1(First, M2, V2, DM2, DV2, 
			   Bid1, Bid2);
    if (!Ok_D1)  myLaws->Value(1)->D0(First, M2, V2);

    E = TopoDS::Edge(myEdges->Value(1));
  }

  if (E.Orientation() == TopAbs_REVERSED)
    V = TopExp::LastVertex(E);
  else 
    V = TopExp::FirstVertex(E);
  
  TolEps += 2*BRep_Tool::Tolerance(V);

  Standard_Boolean isG0 = Standard_True;
  Standard_Boolean isG1 = Standard_True;
  
  if ((V1-V2).Magnitude() > TolEps) isG0 = Standard_False;
  if (Norm(M1-M2) > SpatialTolerance) isG0 = Standard_False;
  
  if (!isG0) return -1;
  if (!Ok_D1) return 0; // No control of the derivative
  
  if ( (DV1.Magnitude()>EpsNul) && (DV2.Magnitude()>EpsNul)
       && (DV1.Angle(DV2) > AngularTolerance) ) isG1 = Standard_False;

  // For the next, the tests are mostly empirical
  Standard_Real Norm1 = Norm(DM1);
  Standard_Real Norm2 = Norm(DM2);
  // It two 2 norms are null, it is good
  if ((Norm1 > EpsNul) || (Norm2 > EpsNul)) {
    // otherwise the normalized matrices are compared
    if ((Norm1 > EpsNul) && (Norm2 > EpsNul)) {
      DM1 /= Norm1;
      DM2 /= Norm2;
      if (Norm(DM1 - DM2) > AngularTolerance) isG1 = Standard_False;
    }
    else isG1 = Standard_False; // 1 Null the other is not   
  }
  
  if (isG1) return 1;
  else return 0;
}


//=======================================================================
//function : Parameter
//purpose  : 
//=======================================================================
 void BRepFill_LocationLaw::Parameter(const Standard_Real Abcissa,
				      Standard_Integer& Index,
				      Standard_Real& U)
{
  Standard_Integer iedge, NbE=myEdges->Length();
  Standard_Boolean Trouve =  Standard_False;

  //Control that the lengths are calculated
  if (myLength->Value(NbE+1) < 0) {
    Standard_Real f, l;
    CurvilinearBounds(NbE, f, l);
  }

  // Find the interval
  for (iedge=1; iedge<=NbE && !Trouve; ) {
    if (myLength->Value(iedge+1) >= Abcissa) {
      Trouve = Standard_True;
    }
    else iedge++;
  }
  
  if (Trouve) {
    Standard_Real f, l;
    const Handle(GeomFill_LocationLaw)& Law =  myLaws->Value(iedge);
    Law->GetDomain(f, l);

    if  (Abcissa == myLength->Value(iedge+1)) {
      U = l;
    }
    else if  (Abcissa == myLength->Value(iedge)) {
      U = f;
    } 
    else {
      GCPnts_AbscissaPoint 
	AbsC(myTol, 
	     myLaws->Value(iedge)->GetCurve()->GetCurve(), 
	     Abcissa-myLength->Value(iedge), f);
      U =  AbsC.Parameter();
    }
    Index = iedge;
  }
  else {
    Index = 0;
  } 
}


//===================================================================
//function : D0
//purpose  : Position of a section, with a given curviline abscissa
//===================================================================
 void BRepFill_LocationLaw::D0(const Standard_Real Abcissa,
			       TopoDS_Shape& W)
{
  Standard_Real u;
  Standard_Integer ind;
  gp_Mat M;
  gp_Vec V;

  Parameter(Abcissa, ind, u);
  if (ind != 0) {
    // Positionement
    myLaws->Value(ind)->D0(u, M, V);
    gp_Trsf fila;
    fila.SetValues(M(1,1), M(1,2), M(1,3), V.X(),
		   M(2,1), M(2,2), M(2,3), V.Y(),
		   M(3,1), M(3,2), M(3,3), V.Z(),
		   1.e-12, 1.e-14);
    TopLoc_Location Loc(fila);
    W.Location(Loc.Multiplied(W.Location()));
  }
  else {
    W.Nullify();
#if DEB    
    cout << "BRepFill_LocationLaw::D0 : Attention position out of limits" 
         << endl;
#endif
  }
}

//=======================================================================
//function : Abscissa
//purpose  : Calculate the abscissa of a point
//=======================================================================
 Standard_Real BRepFill_LocationLaw::Abscissa(const Standard_Integer Index,
					      const Standard_Real Param)
{
  GCPnts_AbscissaPoint AbsC;
  Standard_Real Length = myLength->Value(Index);
  if (Length < 0) {
    Standard_Real bid;
    CurvilinearBounds(Index, bid, Length);
  }
 
  Length += AbsC.Length(myLaws->Value(Index)->GetCurve()->GetCurve(),
			myLaws->Value(Index)->GetCurve()->FirstParameter(),
			Param, myTol);
  return Length;
}
Commit	Line	Data
b311480e	1	// Created on: 1998-01-14
	2	// Created by: Philippe MANGIN
	3	// Copyright (c) 1998-1999 Matra Datavision
	4	// Copyright (c) 1999-2012 OPEN CASCADE SAS
	5	//
	6	// The content of this file is subject to the Open CASCADE Technology Public
	7	// License Version 6.5 (the "License"). You may not use the content of this file
	8	// except in compliance with the License. Please obtain a copy of the License
	9	// at http://www.opencascade.org and read it completely before using this file.
	10	//
	11	// The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
	12	// main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
	13	//
	14	// The Original Code and all software distributed under the License is
	15	// distributed on an "AS IS" basis, without warranty of any kind, and the
	16	// Initial Developer hereby disclaims all such warranties, including without
	17	// limitation, any warranties of merchantability, fitness for a particular
	18	// purpose or non-infringement. Please see the License for the specific terms
	19	// and conditions governing the rights and limitations under the License.
	20
7fd59977	21
	22
	23	#include <BRepFill_LocationLaw.ixx>
	24
	25	#include <BRepTools_WireExplorer.hxx>
	26	#include <BRep_Tool.hxx>
	27	#include <BRep_Builder.hxx>
	28	#include <BRepAdaptor_Curve.hxx>
	29	#include <Adaptor3d_HCurve.hxx>
	30	#include <TopoDS.hxx>
	31	#include <TopoDS_Edge.hxx>
	32	#include <TopExp.hxx>
	33	#include <TopLoc_Location.hxx>
	34
	35	#include <GeomFill_LocationLaw.hxx>
	36	#include <gp_Vec.hxx>
	37	#include <gp_Mat.hxx>
	38	#include <gp_XYZ.hxx>
	39	#include <gp_Trsf.hxx>
	40	#include <GCPnts_AbscissaPoint.hxx>
	41	#include <TColgp_Array1OfPnt2d.hxx>
	42	#include <TColgp_Array1OfVec2d.hxx>
	43	#include <TColStd_SequenceOfInteger.hxx>
	44	#include <Precision.hxx>
	45
	46
	47	//=======================================================================
	48	//function : Norm
0d969553	49	//purpose : Norm of a Matrix
7fd59977	50	//=======================================================================
	51
	52	static Standard_Real Norm(const gp_Mat& M) {
	53	Standard_Real R, Norme;
	54	gp_XYZ Coord;
	55	Coord = M.Row(1);
	56	Norme = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
	57	Coord = M.Row(2);
	58	R = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
	59	if (R>Norme) Norme = R;
	60	Coord = M.Row(3);
	61	R = Abs(Coord.X()) + Abs(Coord.Y())+ Abs(Coord.Z());
	62	if (R>Norme) Norme = R;
	63
	64	return Norme;
	65	}
	66
	67	//=======================================================================
	68	//function : ToG0
0d969553	69	//purpose : Calculate tranformation T such as T.M2 = M1
7fd59977	70	//=======================================================================
	71
	72	static void ToG0(const gp_Mat& M1, const gp_Mat& M2, gp_Mat& T) {
	73	T = M2.Inverted();
	74	T *= M1;
	75	}
	76
	77	//=======================================================================
	78	//function : BRepFill_LocationLaw
	79	//purpose :
	80	//=======================================================================
	81
	82	void BRepFill_LocationLaw::Init(const TopoDS_Wire& Path)
	83
	84	{
	85	Standard_Integer NbEdge;
	86	BRepTools_WireExplorer wexp;
	87	// Class BRep_Tool without fields and without Constructor :
	88	// BRep_Tool B;
	89	TopoDS_Edge E;
	90
	91	myPath = Path;
	92	myTol = 1.e-4;
	93
	94	for (NbEdge=0, wexp.Init(myPath);
	95	wexp.More(); wexp.Next())
	96	// if (! B.Degenerated(wexp.Current())) NbEdge++;
	97	if (! BRep_Tool::Degenerated(wexp.Current())) NbEdge++;
	98
	99
	100	myLaws = new (GeomFill_HArray1OfLocationLaw)(1, NbEdge);
	101	myLength = new (TColStd_HArray1OfReal) (1, NbEdge+1);
	102	myLength->Init(-1.);
	103	myLength->SetValue(1, 0.);
	104	myEdges = new (TopTools_HArray1OfShape) (1, NbEdge);
	105	myDisc.Nullify();
	106	TangentIsMain();
	107	}
	108
	109	//=======================================================================
	110	//function : GetStatus
	111	//purpose :
	112	//=======================================================================
	113	GeomFill_PipeError BRepFill_LocationLaw::GetStatus() const
	114	{
	115	Standard_Integer ii, N = myLaws->Length();
	116	GeomFill_PipeError Status = GeomFill_PipeOk;
	117	for (ii=1; ii<=N && (Status == GeomFill_PipeOk); ii++) {
	118	Status = myLaws->Value(ii)->ErrorStatus();
	119	}
	120	return Status;
	121	}
	122
	123	//=======================================================================
	124	//function : TangentIsMain
	125	//purpose :
	126	//=======================================================================
	127	void BRepFill_LocationLaw::TangentIsMain()
	128	{
	129	myType = 1;
	130	}
	131
	132	//=======================================================================
	133	//function : NormalIsMain
134	//purpose :
135	//=======================================================================
136	void BRepFill_LocationLaw::NormalIsMain()
137	{
138	myType = 2;
139	}
140
141	//=======================================================================
142	//function : BiNormalIsMain
143	//purpose :
144	//=======================================================================
145	void BRepFill_LocationLaw::BiNormalIsMain()
146	{
147	myType = 3;
148	}
149
150	//=======================================================================
151	//function : TransformInCompatibleLaw
0d969553	152	//purpose : Set in continuity of laws
7fd59977	153	//=======================================================================
	154	void BRepFill_LocationLaw::TransformInCompatibleLaw(const Standard_Real TolAngular)
	155	{
	156
	157	Standard_Real First, Last, Angle;
	158	Standard_Integer ipath;
	159	gp_Mat Trsf, M1, M2;
	160	gp_Vec V, T1, T2, N1, N2;
	161	gp_XYZ OZ(0, 0, 1);
	162
	163	myLaws->Value(1)->GetDomain(First, Last);
	164
	165	for (ipath=2; ipath<=myLaws->Length(); ipath++) {
	166	myLaws->Value(ipath-1)->D0(Last, M1, V);
	167	myLaws->Value(ipath)->GetDomain(First, Last);
	168	myLaws->Value(ipath)->D0(First, M2, V);
	169	T1.SetXYZ(M1.Column(3));
	170	T2.SetXYZ(M2.Column(3));
	171	N1.SetXYZ(M1.Column(1));
	172	N2.SetXYZ(M2.Column(1));
	173	if (T1.IsParallel(T2, TolAngular ) &&
	174	!T1.IsOpposite(T2, TolAngular)) { // Correction G0
	175	ToG0(M1, M2, Trsf);
	176	}
	177	else {
	178	Standard_Real alpha;
	179	gp_Vec cross(T1);
	180	cross.Cross(T2);
	181	alpha = T2.AngleWithRef(T1, cross);
	182	gp_Ax1 axe(gp::Origin(), cross.XYZ());
	183	N2.Rotate(axe, alpha);
	184
	185	#if DEB
	186	if (N2.Dot(T1) > 1.e-9) {
0d969553	187	cout << "Inprecision in TransformInCompatibleLaw" << endl;
7fd59977	188	cout << "--- T1.R(N2) = " << N2.Dot(T1) << endl;
	189	gp_Vec tt;
	190	tt = T1;
	191	tt.Rotate(axe, alpha);
	192	cout << "--- T1.R(T2) = " << tt.Dot(T1) << endl;
	193	cout << "--- R(N2).R(T2) = " << N2.Dot(tt) << endl;
	194	}
	195	#endif
	196	Angle = N2.AngleWithRef(N1, T1);
	197	Trsf.SetRotation(OZ, Angle);
	198	}
	199	myLaws->Value(ipath)->SetTrsf(Trsf);
	200	}
	201	}
	202
	203	//=======================================================================
	204	//function : TransformInG0Law
0d969553	205	//purpose : Set in continuity of laws
7fd59977	206	//=======================================================================
	207	void BRepFill_LocationLaw::TransformInG0Law()
	208	{
	209
	210	Standard_Real First, Last;
	211	Standard_Integer ipath;
	212	gp_Mat M1, M2, aux;//,Trsf
	213	gp_Vec V;
	214	myLaws->Value(1)->GetDomain(First, Last);
	215	for (ipath=2; ipath<=myLaws->Length(); ipath++) {
	216	myLaws->Value(ipath-1)->D0(Last, M1, V);
	217	myLaws->Value(ipath)->GetDomain(First, Last);
	218	myLaws->Value(ipath)->D0(First, M2, V);
	219	ToG0(M1, M2, aux);
	220	myLaws->Value(ipath)->SetTrsf(aux);
	221	}
	222
0d969553	223	// Is the law periodical ?
7fd59977	224	if (myPath.Closed()) {
	225	myLaws->Value(myLaws->Length())->D0(Last, M1, V);
	226	myLaws->Value(1)->GetDomain(First, Last);
	227	myLaws->Value(1)->D0(First, M2, V);
	228	}
	229	}
	230
	231	//=======================================================================
	232	//function : DeleteTransform
0d969553	233	//purpose : Remove the setting in continuity of law.
7fd59977	234	//=======================================================================
	235	void BRepFill_LocationLaw::DeleteTransform()
	236	{
	237	gp_Mat Id;
	238	Id.SetIdentity();
	239	for (Standard_Integer ii=1; ii<=myEdges->Length(); ii++) {
	240	myLaws->ChangeValue(ii)->SetTrsf(Id);
	241	}
	242	myDisc.Nullify();
	243	}
	244
	245	//=======================================================================
	246	//function : NbHoles
0d969553	247	//purpose : Find "Holes"
7fd59977	248	//=======================================================================
	249	Standard_Integer BRepFill_LocationLaw::NbHoles(const Standard_Real Tol)
	250	{
	251	if (myDisc.IsNull()) {
	252	TColStd_SequenceOfInteger Seq;
	253	Standard_Integer ii, NbDisc;
	254	for (ii=2, NbDisc=-1; ii<=myLaws->Length()+1; ii++) {
	255	if (IsG1(ii-1, Tol, 1.e-12) == -1) {
	256	Seq.Append(ii);
	257	}
	258	}
	259	NbDisc = Seq.Length();
	260	if ( NbDisc > 0) {
	261	myDisc = new (TColStd_HArray1OfInteger)(1, NbDisc);
	262	for (ii=1; ii<=NbDisc; ii++)
	263	myDisc->SetValue(ii, Seq(ii));
	264	}
	265	}
	266	if (myDisc.IsNull()) return 0;
	267	return myDisc->Length();
	268
	269	}
	270
	271	//=======================================================================
	272	//function : Holes
	273	//purpose :
	274	//=======================================================================
	275	void BRepFill_LocationLaw::Holes(TColStd_Array1OfInteger& Disc) const
	276	{
	277	if (!myDisc.IsNull()) {
	278	for (Standard_Integer ii=1; ii<=myDisc->Length(); ii++)
	279	Disc(ii) = myDisc->Value(ii);
	280	}
	281	}
	282
	283	//=======================================================================
	284	//function : NbLaw
	285	//purpose :
	286	//=======================================================================
	287	Standard_Integer BRepFill_LocationLaw::NbLaw() const
	288	{
	289	return myLaws->Length();
	290	}
	291
	292	//=======================================================================
	293	//function : Law
	294	//purpose :
	295	//=======================================================================
	296	const Handle(GeomFill_LocationLaw)&
	297	BRepFill_LocationLaw::Law(const Standard_Integer Index) const
	298	{
	299	return myLaws->Value(Index);
	300	}
	301
	302	//=======================================================================
	303	//function : Wire
	304	//purpose :
	305	//=======================================================================
	306	const TopoDS_Wire& BRepFill_LocationLaw::Wire() const
	307	{
	308	return myPath;
	309	}
	310
	311	//=======================================================================
312	//function : Edge
313	//purpose :
314	//=======================================================================
315	const TopoDS_Edge& BRepFill_LocationLaw::Edge(const Standard_Integer Index) const
316	{
317	return TopoDS::Edge(myEdges->Value(Index));
318	}
319
320	//=======================================================================
321	//function : Vertex
322	//purpose :
323	//=======================================================================
324	TopoDS_Vertex BRepFill_LocationLaw::Vertex(const Standard_Integer Index) const
325	{
326	TopoDS_Edge E;
327	TopoDS_Vertex V;
328	if (Index <= myEdges->Length()) {
329	E = TopoDS::Edge(myEdges->Value(Index));
330	if (E.Orientation() == TopAbs_REVERSED)
331	V = TopExp::LastVertex(E);
332	else V = TopExp::FirstVertex(E);
333	}
334	else if (Index == myEdges->Length()+1) {
335	E = TopoDS::Edge(myEdges->Value(Index-1));
336	if (E.Orientation() == TopAbs_REVERSED)
337	V = TopExp::FirstVertex(E);
338	else V = TopExp::LastVertex(E);
339	}
340	return V;
341	}
342
0d969553	343	//===================================================================
7fd59977	344	//function : PerformVertex
0d969553 Y	345	//purpose : Calculate a vertex of sweeping from a vertex of section
	346	// and the index of the edge in the trajectory
	347	//===================================================================
7fd59977	348	void BRepFill_LocationLaw::PerformVertex(const Standard_Integer Index,
	349	const TopoDS_Vertex& Input,
	350	const Standard_Real TolMin,
	351	TopoDS_Vertex& Output,
	352	const Standard_Integer ILoc) const
	353	{
	354	BRep_Builder B;
	355	Standard_Boolean IsBary = (ILoc == 0);
	356	Standard_Real First, Last;
	357	gp_Pnt P;
	358	gp_Vec V1, V2;//, V;
	359	gp_Mat M1, M2;
	360
	361	if (Index>0 && Index<myLaws->Length()) {
	362	if (ILoc <=0) {
	363	myLaws->Value(Index)->GetDomain(First, Last);
	364	myLaws->Value(Index)->D0(Last, M1, V1);
	365	}
	366
	367	if (ILoc >= 0) {
	368	myLaws->Value(Index+1)->GetDomain(First, Last);
	369	if (ILoc == 0)
	370	myLaws->Value(Index+1)->D0(First, M2, V2);
	371	else
	372	myLaws->Value(Index+1)->D0(First, M1, V1);
	373	}
	374	}
	375
	376	if (Index == 0 \|\| Index == myLaws->Length()) {
	377	if (!myPath.Closed() \|\| (IsG1(Index, TolMin) != 1)) {
	378	IsBary = Standard_False;
	379	if (Index == 0) {
	380	myLaws->Value(1)->GetDomain(First, Last);
	381	myLaws->Value(1)->D0(First, M1, V1);
	382	}
	383	else {
	384	myLaws->Value(myLaws->Length())->GetDomain(First, Last);
	385	myLaws->Value(myLaws->Length())->D0(Last, M1, V1);
	386	}
	387	}
	388	else {
	389	if (ILoc <=0) {
	390	myLaws->Value(myLaws->Length())->GetDomain(First, Last);
	391	myLaws->Value(myLaws->Length())->D0(Last, M1, V1);
	392	}
	393
	394	if (ILoc >=0) {
	395	myLaws->Value(1)->GetDomain(First, Last);
	396	if (ILoc==0)
	397	myLaws->Value(1)->D0(First, M2, V2);
	398	else
	399	myLaws->Value(1)->D0(First, M1, V1);
	400	}
	401	}
	402	}
	403
	404	P = BRep_Tool::Pnt(Input);
	405
	406	if (IsBary) {
	407	gp_XYZ P1(P.XYZ()), P2(P.XYZ());
	408	P1 *= M1;
	409	P1 += V1.XYZ();
	410	P2 *= M2;
	411	P2 += V2.XYZ();
412
413	P.ChangeCoord().SetLinearForm(0.5, P1, 0.5, P2);
414	P1 -= P2;
415	Standard_Real Tol = P1.Modulus()/2;
416	Tol += TolMin;
417	B.MakeVertex(Output, P, Tol);
418	}
419	else {
420	P.ChangeCoord() *= M1;
421	P.ChangeCoord() += V1.XYZ();
422	B.MakeVertex(Output, P, TolMin);
423	}
424
425	}
426
427	//=======================================================================
428	//function : CurvilinearBounds
429	//purpose :
430	//=======================================================================
431	void BRepFill_LocationLaw::CurvilinearBounds(const Standard_Integer Index,
432	Standard_Real& First,
433	Standard_Real& Last) const
434	{
435	First = myLength->Value(Index);
436	Last = myLength->Value(Index+1);
0d969553	437	if (Last<0) { //It is required to carry out the calculation
7fd59977	438	Standard_Integer ii, NbE = myEdges->Length();
	439	Standard_Real Length, f, l;
	440	GCPnts_AbscissaPoint AbsC;
	441
	442	for (ii=1, Length=0.; ii<=NbE; ii++) {
	443	myLaws->Value(ii)->GetDomain(f, l);
	444	Length += AbsC.Length(myLaws->Value(ii)->GetCurve()->GetCurve(), myTol);
	445	myLength->SetValue(ii+1, Length);
	446	}
	447
	448	First = myLength->Value(Index);
	449	Last = myLength->Value(Index+1);
	450	}
	451	}
	452
	453	Standard_Boolean BRepFill_LocationLaw::IsClosed() const
	454	{
	455	return myPath.Closed();
	456	}
	457
	458	//=======================================================================
	459	//function : IsG1
0d969553	460	//purpose : Evaluate the continuity of the law by a vertex
7fd59977	461	//=======================================================================
	462	Standard_Integer
	463	BRepFill_LocationLaw::IsG1(const Standard_Integer Index,
	464	const Standard_Real SpatialTolerance,
	465	const Standard_Real AngularTolerance) const
	466	{
	467	gp_Vec V1, DV1, V2, DV2;
	468	gp_Mat M1, M2, DM1, DM2;
	469	Standard_Real First, Last, EpsNul = 1.e-12;
	470	Standard_Real TolEps = SpatialTolerance;
	471	#ifndef DEB
	472	Standard_Boolean Ok_D1 = Standard_False;
	473	#else
	474	Standard_Boolean Ok_D1;
	475	#endif
	476	TopoDS_Vertex V;
	477	TopoDS_Edge E;
	478	TColgp_Array1OfPnt2d Bid1 (1,1);
	479	TColgp_Array1OfVec2d Bid2 (1,1);
	480
	481	if (Index>0 && Index<myLaws->Length()) {
	482	myLaws->Value(Index)->GetDomain(First, Last);
	483	Ok_D1 = myLaws->Value(Index)->D1(Last, M1, V1, DM1, DV1,
	484	Bid1, Bid2);
	485	if (!Ok_D1) myLaws->Value(Index)->D0(Last, M1, V1);
	486
	487	myLaws->Value(Index+1)->GetDomain(First, Last);
	488	if (Ok_D1)
	489	Ok_D1 = myLaws->Value(Index+1)->D1(First, M2, V2, DM2, DV2,
	490	Bid1, Bid2);
	491	if (!Ok_D1) myLaws->Value(Index+1)->D0(First, M2, V2);
	492
	493	E = TopoDS::Edge(myEdges->Value(Index+1));
	494	}
	495	if (Index == 0 \|\| Index == myLaws->Length()) {
	496	if (!myPath.Closed()) return -1;
	497	myLaws->Value(myLaws->Length())->GetDomain(First, Last);
	498	Ok_D1 = myLaws->Value(myLaws->Length())->D1(Last, M1, V1, DM1, DV1,
	499	Bid1, Bid2);
	500	if (!Ok_D1) myLaws->Value(myLaws->Length())->D0(Last, M1, V1);
	501
	502	myLaws->Value(1)->GetDomain(First, Last);
	503	if (Ok_D1)
	504	myLaws->Value(1)->D1(First, M2, V2, DM2, DV2,
	505	Bid1, Bid2);
	506	if (!Ok_D1) myLaws->Value(1)->D0(First, M2, V2);
	507
	508	E = TopoDS::Edge(myEdges->Value(1));
	509	}
	510
	511	if (E.Orientation() == TopAbs_REVERSED)
	512	V = TopExp::LastVertex(E);
	513	else
	514	V = TopExp::FirstVertex(E);
	515
	516	TolEps += 2*BRep_Tool::Tolerance(V);
	517
	518	Standard_Boolean isG0 = Standard_True;
	519	Standard_Boolean isG1 = Standard_True;
	520
	521	if ((V1-V2).Magnitude() > TolEps) isG0 = Standard_False;
	522	if (Norm(M1-M2) > SpatialTolerance) isG0 = Standard_False;
	523
	524	if (!isG0) return -1;
0d969553	525	if (!Ok_D1) return 0; // No control of the derivative
7fd59977	526
	527	if ( (DV1.Magnitude()>EpsNul) && (DV2.Magnitude()>EpsNul)
	528	&& (DV1.Angle(DV2) > AngularTolerance) ) isG1 = Standard_False;
	529
0d969553	530	// For the next, the tests are mostly empirical
7fd59977	531	Standard_Real Norm1 = Norm(DM1);
7fd59977	532	Standard_Real Norm2 = Norm(DM2);
0d969553	533	// It two 2 norms are null, it is good
7fd59977	534	if ((Norm1 > EpsNul) \|\| (Norm2 > EpsNul)) {
0d969553	535	// otherwise the normalized matrices are compared
7fd59977	536	if ((Norm1 > EpsNul) && (Norm2 > EpsNul)) {
	537	DM1 /= Norm1;
	538	DM2 /= Norm2;
	539	if (Norm(DM1 - DM2) > AngularTolerance) isG1 = Standard_False;
	540	}
0d969553	541	else isG1 = Standard_False; // 1 Null the other is not
7fd59977	542	}
	543
	544	if (isG1) return 1;
	545	else return 0;
	546	}
	547
	548
	549	//=======================================================================
	550	//function : Parameter
	551	//purpose :
	552	//=======================================================================
	553	void BRepFill_LocationLaw::Parameter(const Standard_Real Abcissa,
	554	Standard_Integer& Index,
	555	Standard_Real& U)
	556	{
	557	Standard_Integer iedge, NbE=myEdges->Length();
	558	Standard_Boolean Trouve = Standard_False;
	559
0d969553	560	//Control that the lengths are calculated
7fd59977	561	if (myLength->Value(NbE+1) < 0) {
	562	Standard_Real f, l;
	563	CurvilinearBounds(NbE, f, l);
	564	}
	565
0d969553	566	// Find the interval
7fd59977	567	for (iedge=1; iedge<=NbE && !Trouve; ) {
	568	if (myLength->Value(iedge+1) >= Abcissa) {
	569	Trouve = Standard_True;
	570	}
	571	else iedge++;
	572	}
	573
	574	if (Trouve) {
	575	Standard_Real f, l;
	576	const Handle(GeomFill_LocationLaw)& Law = myLaws->Value(iedge);
	577	Law->GetDomain(f, l);
	578
	579	if (Abcissa == myLength->Value(iedge+1)) {
	580	U = l;
	581	}
	582	else if (Abcissa == myLength->Value(iedge)) {
	583	U = f;
	584	}
	585	else {
	586	GCPnts_AbscissaPoint
	587	AbsC(myTol,
	588	myLaws->Value(iedge)->GetCurve()->GetCurve(),
	589	Abcissa-myLength->Value(iedge), f);
	590	U = AbsC.Parameter();
	591	}
	592	Index = iedge;
	593	}
	594	else {
	595	Index = 0;
	596	}
	597	}
	598
	599
0d969553	600	//===================================================================
7fd59977	601	//function : D0
0d969553 Y	602	//purpose : Position of a section, with a given curviline abscissa
0d969553 Y	603	//===================================================================
7fd59977	604	void BRepFill_LocationLaw::D0(const Standard_Real Abcissa,
	605	TopoDS_Shape& W)
	606	{
	607	Standard_Real u;
	608	Standard_Integer ind;
	609	gp_Mat M;
	610	gp_Vec V;
	611
	612	Parameter(Abcissa, ind, u);
	613	if (ind != 0) {
	614	// Positionement
	615	myLaws->Value(ind)->D0(u, M, V);
	616	gp_Trsf fila;
	617	fila.SetValues(M(1,1), M(1,2), M(1,3), V.X(),
	618	M(2,1), M(2,2), M(2,3), V.Y(),
	619	M(3,1), M(3,2), M(3,3), V.Z(),
	620	1.e-12, 1.e-14);
	621	TopLoc_Location Loc(fila);
	622	W.Location(Loc.Multiplied(W.Location()));
	623	}
	624	else {
	625	W.Nullify();
	626	#if DEB
0d969553	627	cout << "BRepFill_LocationLaw::D0 : Attention position out of limits"
7fd59977	628	<< endl;
	629	#endif
	630	}
	631	}
	632
	633	//=======================================================================
	634	//function : Abscissa
0d969553	635	//purpose : Calculate the abscissa of a point
7fd59977	636	//=======================================================================
	637	Standard_Real BRepFill_LocationLaw::Abscissa(const Standard_Integer Index,
	638	const Standard_Real Param)
	639	{
	640	GCPnts_AbscissaPoint AbsC;
	641	Standard_Real Length = myLength->Value(Index);
	642	if (Length < 0) {
	643	Standard_Real bid;
	644	CurvilinearBounds(Index, bid, Length);
	645	}
	646
	647	Length += AbsC.Length(myLaws->Value(Index)->GetCurve()->GetCurve(),
	648	myLaws->Value(Index)->GetCurve()->FirstParameter(),
	649	Param, myTol);
	650	return Length;
	651	}