[occt.git] / src / Contap / Contap_ContAna.cxx

// Created on: 1993-03-04
// Created by: Jacques GOUSSARD
// Copyright (c) 1993-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 <Contap_ContAna.ixx>

#include <gp_XYZ.hxx>
#include <gp.hxx>

#define Tolpetit 1.e-8

Contap_ContAna::Contap_ContAna (): done(Standard_False) {}

void Contap_ContAna::Perform (const gp_Sphere& S,
                              const gp_Dir& D)
{
  done  = Standard_False;
  typL  = GeomAbs_Circle;
  pt1   = S.Location();
  dir1  = D;
  if (Abs(D.Dot(S.XAxis().Direction())) < 0.9999999999999) {
    dir2 = D.Crossed(S.XAxis().Direction());
  }
  else {
    dir2 = D.Crossed(S.YAxis().Direction());
  }
  prm   = S.Radius();
  nbSol = 1;
  done  = Standard_True;
}

void Contap_ContAna::Perform (const gp_Sphere& S,
                              const gp_Dir& D,
                              const Standard_Real Angle)
{
  done  = Standard_False;
  typL  = GeomAbs_Circle;

  dir1  = D;
  if (Abs(D.Dot(S.XAxis().Direction())) < 0.9999999999999) {
    dir2 = D.Crossed(S.XAxis().Direction());
  }
  else {
    dir2 = D.Crossed(S.YAxis().Direction());
  }
  Standard_Real alpha = (S.Direct() ? Angle : -Angle); 
  pt1.SetXYZ(S.Location().XYZ() - S.Radius()*sin(alpha)*D.XYZ()) ;
  prm   = S.Radius()*cos(alpha);
  nbSol = 1;
  done  = Standard_True;
}

void Contap_ContAna::Perform (const gp_Sphere& S,
                              const gp_Pnt& Eye)
{
  done = Standard_False;

  Standard_Real radius = S.Radius();
  Standard_Real dist = Eye.Distance(S.Location());
  if (dist <= radius) {
    nbSol = 0;
  }
  else {
    prm = radius*sqrt(1.-radius*radius/(dist*dist));
    if (prm < Tolpetit) {
      nbSol = 0;
    }
    else {
      gp_XYZ locxyz(S.Location().XYZ());
      dir1.SetXYZ(Eye.XYZ()-locxyz);
      pt1.SetXYZ(locxyz + (radius*radius/dist)*dir1.XYZ());
      if (Abs(dir1.Dot(S.XAxis().Direction())) < 0.9999999999999) {
	dir2 = dir1.Crossed(S.XAxis().Direction());
      }
      else {
	dir2 = dir1.Crossed(S.YAxis().Direction());
      }
      nbSol = 1;
      typL = GeomAbs_Circle;
    }
  }
  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cylinder& C,
                              const gp_Dir& D)
{
  done = Standard_False;

  gp_XYZ normale(C.Position().Direction().XYZ());
  normale.Cross(D.XYZ());
  if (normale.Modulus() <= 1e-15) {
    nbSol = 0;
  }
  else {
    normale.Normalize();
    typL = GeomAbs_Line;
    dir1 = C.Position().Direction();
    dir2 = dir1;
    pt1.SetXYZ(C.Location().XYZ() + C.Radius()*normale);
    pt2.SetXYZ(C.Location().XYZ() - C.Radius()*normale);
    nbSol = 2;
  }

  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cylinder& C,
                              const gp_Dir& D,
                              const Standard_Real Angle)
{
  done = Standard_False;

  Standard_Real Coefcos = D.Dot(C.Position().XDirection());
  Standard_Real Coefsin = D.Dot(C.Position().YDirection());
  Standard_Real Coefcst = cos(M_PI*0.5 + Angle);

  Standard_Real norm1 = Coefcos*Coefcos + Coefsin*Coefsin;
  Standard_Real norm2 = sqrt(norm1);

  if (Abs(Coefcst) < norm2) {
    typL = GeomAbs_Line;
    nbSol = 2;
    dir1 = dir2 = C.Position().Direction();

    if (!C.Direct()) { // The normal is inverted.
      Coefcos = -Coefcos;
      Coefsin = -Coefsin;
    }

    // Necessary to solve Coefcos*cos(t) + Coefsin*sin(t) = Coefcst
    // and the origins of solution are in the reference of the 
    // cylinder in (R*cost0, R*sint0,0) and (R*cost1,R*sint1,0)
    // By setting cos(phi) = Coefcos/Sqrt(Coefcos**2 + Coefsin**2) and
    //           sin(phi) = Coefsin/Sqrt(Coefcos**2 + Coefsin**2)
    // and by using trigonometric relations the values of cosinus 
    // and sinus to the solutions are obtained.

    prm = Sqrt(norm1 - Coefcst*Coefcst);
    Standard_Real cost0,sint0,cost1,sint1;

    cost0 = (Coefcos*Coefcst - Coefsin*prm)/norm1;
    cost1 = (Coefcos*Coefcst + Coefsin*prm)/norm1;

    sint0 = ( Coefcos*prm + Coefsin*Coefcst)/norm1;
    sint1 = (-Coefcos*prm + Coefsin*Coefcst)/norm1;

    gp_XYZ Xdir(C.Position().XDirection().XYZ());
    gp_XYZ Ydir(C.Position().YDirection().XYZ());
    gp_XYZ dirxyz;

    dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir);
    dirxyz.Multiply(C.Radius());
    pt1.SetXYZ(C.Location().XYZ().Added(dirxyz));

    dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir);
    dirxyz.Multiply(C.Radius());
    pt2.SetXYZ(C.Location().XYZ().Added(dirxyz));
  }
  else {
    nbSol = 0;
  }

  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cylinder& C,
                              const gp_Pnt& Eye)
{
  done = Standard_False;

  Standard_Real radius = C.Radius();
  gp_Lin theaxis(C.Axis());
  Standard_Real dist = theaxis.Distance(Eye);
  if (dist <= radius) {
    nbSol = 0;
  }
  else {
    typL = GeomAbs_Line;
    prm = radius*sqrt(1.-radius*radius/(dist*dist));
    dir1 = C.Axis().Direction();
    dir2 = dir1;
    gp_XYZ axeye(theaxis.Normal(Eye).Direction().XYZ()); // orientate the axis to the outside
    gp_XYZ normale((theaxis.Direction().Crossed(axeye)).XYZ());
//      normale.Normalize();
    pt1.SetXYZ(C.Location().XYZ() + (radius*radius/dist)*axeye);
    pt2.SetXYZ(pt1.XYZ() - prm*normale);
    pt1.SetXYZ(pt1.XYZ() + prm*normale);
    nbSol = 2;
  }
  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cone& C,
                              const gp_Dir& D)
{
  done = Standard_False;

  Standard_Real Tgtalpha = Tan(C.SemiAngle());

  Standard_Real Coefcos = D.Dot(C.Position().XDirection());
  Standard_Real Coefsin = D.Dot(C.Position().YDirection());
  Standard_Real Coefcst = D.Dot(C.Axis().Direction())*Tgtalpha;

  Standard_Real norm1 = Coefcos*Coefcos + Coefsin*Coefsin;
  Standard_Real norm2 = Sqrt(norm1);
//  if (Abs(Abs(Coefcst)-norm2) <= Tolpetit) { // tol angulaire 1.e-8
//    typL = GeomAbs_Line;
//    nbSol = 1;
//    pt1 = C.Apex();
//    dir1 = D;
//  }
//  else if (Abs(Coefcst) < norm2) {

  if (Abs(Coefcst) < norm2) {
    typL = GeomAbs_Line;
    nbSol = 2;
    pt1 = C.Apex();
    pt2 = pt1;
    // Necessary to solve Coefcos*cos(t) + Coefsin*sin(t) = Coefcst
    // and director vectors of solutions are
    // cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
    // By setting cos(phi) = Coefcos/Sqrt(Coefcos**2 + Coefsin**2) and
    //           sin(phi) = Coefsin/Sqrt(Coefcos**2 + Coefsin**2)
    // and by using trigonometric relations the values of cosinus 
    // and sinus to the solutions are obtained.

    prm = Sqrt(norm1 - Coefcst*Coefcst);
    Standard_Real cost0,sint0,cost1,sint1;

    cost0 = (Coefcos*Coefcst - Coefsin*prm)/norm1;
    cost1 = (Coefcos*Coefcst + Coefsin*prm)/norm1;

    sint0 = ( Coefcos*prm + Coefsin*Coefcst)/norm1;
    sint1 = (-Coefcos*prm + Coefsin*Coefcst)/norm1;
    
    gp_XYZ Xdir(C.Position().XDirection().XYZ());
    gp_XYZ Ydir(C.Position().YDirection().XYZ());
    gp_XYZ Zdir(C.Axis().Direction().XYZ());
    gp_XYZ dirxyz;
    dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,1./Tgtalpha,Zdir);
    dir1.SetXYZ(dirxyz);
    pt1.SetXYZ(pt1.XYZ()+dirxyz);
    dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,1./Tgtalpha,Zdir);
    dir2.SetXYZ(dirxyz);
    pt2.SetXYZ(pt2.XYZ()+dirxyz);
  }
  else {
    nbSol = 0;
  }
  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cone& C,
                              const gp_Dir& D,
                              const Standard_Real Angle)
{
  done = Standard_False;
  nbSol = 0;

  Standard_Real Ang = C.SemiAngle();
  Standard_Real Cosa = cos(Ang);
  Standard_Real Sina = sin(Ang);

  Standard_Real Coefcos = D.Dot(C.Position().XDirection());
  Standard_Real Coefsin = D.Dot(C.Position().YDirection());

  Standard_Real Coefcst1 = cos(M_PI*0.5 + Angle);

  Standard_Real norm1 = Coefcos*Coefcos + Coefsin*Coefsin;
  Standard_Real norm2 = Sqrt(norm1);

  Standard_Real Coefnz = D.Dot(C.Axis().Direction())*Sina;
  Standard_Real Coefcst = (Coefcst1 + Coefnz)/Cosa;

  if (Abs(Coefcst) < norm2) {
    typL = GeomAbs_Line;
    nbSol+= 2;
    pt1 = C.Apex();
    pt2 = pt1;

    // It is requiredto solve Coefcos*cos(t) + Coefsin*sin(t) = Coefcst
    // and the director vectors of solutions are
    // cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
    // By setting cos(phi) = Coefcos/Sqrt(Coefcos**2 + Coefsin**2) and
    //           sin(phi) = Coefsin/Sqrt(Coefcos**2 + Coefsin**2)
    // and by using trigonometric relations the values of cosinus 
    // and sinus to the solutions are obtained.

    prm = Sqrt(norm1 - Coefcst*Coefcst);
    Standard_Real cost0,sint0,cost1,sint1;

    cost0 = (Coefcos*Coefcst - Coefsin*prm)/norm1;
    cost1 = (Coefcos*Coefcst + Coefsin*prm)/norm1;

    sint0 = ( Coefcos*prm + Coefsin*Coefcst)/norm1;
    sint1 = (-Coefcos*prm + Coefsin*Coefcst)/norm1;

    gp_XYZ Xdir(C.Position().XDirection().XYZ());
    gp_XYZ Ydir(C.Position().YDirection().XYZ());
    gp_XYZ Zdir(C.Axis().Direction().XYZ());
    if (!C.Direct()) {
      Zdir.Reverse();
    }
    gp_XYZ dirxyz;
    dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,Cosa/Sina,Zdir);
    dir1.SetXYZ(dirxyz);
    pt1.SetXYZ(pt1.XYZ()+dirxyz);
    dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,Cosa/Sina,Zdir);
    dir2.SetXYZ(dirxyz);
    pt2.SetXYZ(pt2.XYZ()+dirxyz);
  }

  Coefcst = (Coefcst1 - Coefnz)/Cosa;

  if (Abs(Coefcst) < norm2) {
    typL = GeomAbs_Line;
    nbSol+= 2;
    pt3 = C.Apex();
    pt4 = pt3;

    prm = Sqrt(norm1 - Coefcst*Coefcst);
    Standard_Real cost0,sint0,cost1,sint1;

    cost0 = (Coefcos*Coefcst - Coefsin*prm)/norm1;
    cost1 = (Coefcos*Coefcst + Coefsin*prm)/norm1;

    sint0 = ( Coefcos*prm + Coefsin*Coefcst)/norm1;
    sint1 = (-Coefcos*prm + Coefsin*Coefcst)/norm1;

    gp_XYZ Xdir(C.Position().XDirection().XYZ());
    gp_XYZ Ydir(C.Position().YDirection().XYZ());
    gp_XYZ Zdir(C.Axis().Direction().XYZ());
    if (!C.Direct()) {
      Zdir.Reverse();
    }
    gp_XYZ dirxyz;
    dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,-Cosa/Sina,Zdir);
    dir3.SetXYZ(dirxyz);
    pt3.SetXYZ(pt3.XYZ()+dirxyz);
    dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,-Cosa/Sina,Zdir);
    dir4.SetXYZ(dirxyz);
    pt4.SetXYZ(pt4.XYZ()+dirxyz);
    if (nbSol == 2) {
      pt1 = pt3;
      pt2 = pt4;
      dir1 = dir3;
      dir2 = dir4;
    }
  }

  done = Standard_True;
}

void Contap_ContAna::Perform (const gp_Cone& C,
                              const gp_Pnt& Eye)
{
  done = Standard_False;

  Standard_Real Tgtalpha = Tan(C.SemiAngle());

  gp_XYZ apexeye(Eye.XYZ());
  apexeye.Subtract(C.Apex().XYZ());

  Standard_Real Coefcos = apexeye.Dot(C.Position().XDirection().XYZ());
  Standard_Real Coefsin = apexeye.Dot(C.Position().YDirection().XYZ());
  Standard_Real Coefcst = apexeye.Dot(C.Axis().Direction().XYZ())*Tgtalpha;

  Standard_Real norm1 = Coefcos*Coefcos + Coefsin*Coefsin;
  Standard_Real norm2 = Sqrt(Coefcos*Coefcos + Coefsin*Coefsin);
//  if (Abs(Abs(Coefcst)-norm2) <= Tolpetit) { // tol angulaire 1.e-8
//    typL = GeomAbs_Line;
//    nbSol = 1;
//    pt1 = C.Apex();
//    dir1.SetXYZ(apexeye);
//  }
//  else if (Abs(Coefcst) < norm2) {

  if (Abs(Coefcst) < norm2) {
    typL = GeomAbs_Line;
    nbSol = 2;
    pt1 = C.Apex();
    pt2 = pt1;
    // It is required to solve Coefcos*cos(t) + Coefsin*sin(t) = Coefcst
    // and the director vectors of solutions are
    // cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
    // By setting cos(phi) = Coefcos/Sqrt(Coefcos**2 + Coefsin**2) and
    //           sin(phi) = Coefsin/Sqrt(Coefcos**2 + Coefsin**2)
    // and by using trigonometric relations the values of cosinus 
    // and sinus to the solutions are obtained.

    prm = Sqrt(norm1 - Coefcst*Coefcst);
    Standard_Real cost0,sint0,cost1,sint1;

    cost0 = (Coefcos*Coefcst - Coefsin*prm)/norm1;
    cost1 = (Coefcos*Coefcst + Coefsin*prm)/norm1;

    sint0 = ( Coefcos*prm + Coefsin*Coefcst)/norm1;
    sint1 = (-Coefcos*prm + Coefsin*Coefcst)/norm1;

    gp_XYZ Xdir(C.Position().XDirection().XYZ());
    gp_XYZ Ydir(C.Position().YDirection().XYZ());
    gp_XYZ Zdir(C.Axis().Direction().XYZ());
    gp_XYZ dirxyz;
    dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,1./Tgtalpha,Zdir);
    dir1.SetXYZ(dirxyz);
    pt1.SetXYZ(pt1.XYZ()+dirxyz);
    dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,1./Tgtalpha,Zdir);
    dir2.SetXYZ(dirxyz);
    pt2.SetXYZ(pt2.XYZ()+dirxyz);
  }
  else {
    nbSol = 0;
  }
  done = Standard_True;
}

gp_Lin Contap_ContAna::Line (const Standard_Integer Index) const
{
  if (!done) {StdFail_NotDone::Raise();}
  if (typL != GeomAbs_Line || nbSol == 0) {Standard_DomainError::Raise();}
  if (Index <=0 || Index > nbSol) {Standard_OutOfRange::Raise();}
  switch (Index) {
  case 1:
    return gp_Lin(pt1,dir1);
  case 2:
    return gp_Lin(pt2,dir2);
  case 3:
    return gp_Lin(pt3,dir3);
  case 4:
    return gp_Lin(pt4,dir4);
  }
  Standard_OutOfRange::Raise("Program error in Contap_ContAna");
  return gp_Lin();
}
Commit	Line	Data
b311480e	1	// Created on: 1993-03-04
	2	// Created by: Jacques GOUSSARD
	3	// Copyright (c) 1993-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
733a0e55	21
7fd59977	22	#include <Contap_ContAna.ixx>
	23
	24	#include <gp_XYZ.hxx>
	25	#include <gp.hxx>
	26
7fd59977	27	#define Tolpetit 1.e-8
7fd59977	28
7fd59977	29	Contap_ContAna::Contap_ContAna (): done(Standard_False) {}
	30
	31	void Contap_ContAna::Perform (const gp_Sphere& S,
733a0e55	32	const gp_Dir& D)
7fd59977	33	{
	34	done = Standard_False;
	35	typL = GeomAbs_Circle;
	36	pt1 = S.Location();
	37	dir1 = D;
	38	if (Abs(D.Dot(S.XAxis().Direction())) < 0.9999999999999) {
	39	dir2 = D.Crossed(S.XAxis().Direction());
	40	}
	41	else {
	42	dir2 = D.Crossed(S.YAxis().Direction());
	43	}
	44	prm = S.Radius();
	45	nbSol = 1;
	46	done = Standard_True;
	47	}
	48
	49	void Contap_ContAna::Perform (const gp_Sphere& S,
733a0e55 S	50	const gp_Dir& D,
733a0e55 S	51	const Standard_Real Angle)
7fd59977	52	{
	53	done = Standard_False;
	54	typL = GeomAbs_Circle;
	55
	56	dir1 = D;
	57	if (Abs(D.Dot(S.XAxis().Direction())) < 0.9999999999999) {
	58	dir2 = D.Crossed(S.XAxis().Direction());
	59	}
	60	else {
	61	dir2 = D.Crossed(S.YAxis().Direction());
	62	}
	63	Standard_Real alpha = (S.Direct() ? Angle : -Angle);
	64	pt1.SetXYZ(S.Location().XYZ() - S.Radius()sin(alpha)D.XYZ()) ;
	65	prm = S.Radius()*cos(alpha);
	66	nbSol = 1;
	67	done = Standard_True;
	68	}
	69
	70	void Contap_ContAna::Perform (const gp_Sphere& S,
733a0e55	71	const gp_Pnt& Eye)
7fd59977	72	{
	73	done = Standard_False;
	74
	75	Standard_Real radius = S.Radius();
	76	Standard_Real dist = Eye.Distance(S.Location());
	77	if (dist <= radius) {
	78	nbSol = 0;
	79	}
	80	else {
	81	prm = radiussqrt(1.-radiusradius/(dist*dist));
	82	if (prm < Tolpetit) {
	83	nbSol = 0;
	84	}
	85	else {
	86	gp_XYZ locxyz(S.Location().XYZ());
	87	dir1.SetXYZ(Eye.XYZ()-locxyz);
	88	pt1.SetXYZ(locxyz + (radiusradius/dist)dir1.XYZ());
	89	if (Abs(dir1.Dot(S.XAxis().Direction())) < 0.9999999999999) {
	90	dir2 = dir1.Crossed(S.XAxis().Direction());
	91	}
	92	else {
	93	dir2 = dir1.Crossed(S.YAxis().Direction());
	94	}
	95	nbSol = 1;
	96	typL = GeomAbs_Circle;
	97	}
	98	}
	99	done = Standard_True;
	100	}
	101
7fd59977	102	void Contap_ContAna::Perform (const gp_Cylinder& C,
733a0e55	103	const gp_Dir& D)
7fd59977	104	{
	105	done = Standard_False;
	106
	107	gp_XYZ normale(C.Position().Direction().XYZ());
	108	normale.Cross(D.XYZ());
	109	if (normale.Modulus() <= 1e-15) {
	110	nbSol = 0;
	111	}
	112	else {
	113	normale.Normalize();
	114	typL = GeomAbs_Line;
	115	dir1 = C.Position().Direction();
	116	dir2 = dir1;
	117	pt1.SetXYZ(C.Location().XYZ() + C.Radius()*normale);
	118	pt2.SetXYZ(C.Location().XYZ() - C.Radius()*normale);
	119	nbSol = 2;
	120	}
	121
	122	done = Standard_True;
	123	}
	124
	125	void Contap_ContAna::Perform (const gp_Cylinder& C,
733a0e55 S	126	const gp_Dir& D,
733a0e55 S	127	const Standard_Real Angle)
7fd59977	128	{
	129	done = Standard_False;
	130
	131	Standard_Real Coefcos = D.Dot(C.Position().XDirection());
	132	Standard_Real Coefsin = D.Dot(C.Position().YDirection());
c6541a0c	133	Standard_Real Coefcst = cos(M_PI*0.5 + Angle);
7fd59977	134
	135	Standard_Real norm1 = CoefcosCoefcos + CoefsinCoefsin;
	136	Standard_Real norm2 = sqrt(norm1);
	137
	138	if (Abs(Coefcst) < norm2) {
	139	typL = GeomAbs_Line;
	140	nbSol = 2;
	141	dir1 = dir2 = C.Position().Direction();
	142
0d969553	143	if (!C.Direct()) { // The normal is inverted.
7fd59977	144	Coefcos = -Coefcos;
	145	Coefsin = -Coefsin;
	146	}
	147
0d969553 Y	148	// Necessary to solve Coefcoscos(t) + Coefsinsin(t) = Coefcst
	149	// and the origins of solution are in the reference of the
	150	// cylinder in (Rcost0, Rsint0,0) and (Rcost1,Rsint1,0)
	151	// By setting cos(phi) = Coefcos/Sqrt(Coefcos2 + Coefsin2) and
7fd59977	152	// sin(phi) = Coefsin/Sqrt(Coefcos2 + Coefsin2)
0d969553 Y	153	// and by using trigonometric relations the values of cosinus
0d969553 Y	154	// and sinus to the solutions are obtained.
7fd59977	155
	156	prm = Sqrt(norm1 - Coefcst*Coefcst);
	157	Standard_Real cost0,sint0,cost1,sint1;
	158
	159	cost0 = (CoefcosCoefcst - Coefsinprm)/norm1;
	160	cost1 = (CoefcosCoefcst + Coefsinprm)/norm1;
	161
	162	sint0 = ( Coefcosprm + CoefsinCoefcst)/norm1;
	163	sint1 = (-Coefcosprm + CoefsinCoefcst)/norm1;
	164
	165	gp_XYZ Xdir(C.Position().XDirection().XYZ());
	166	gp_XYZ Ydir(C.Position().YDirection().XYZ());
	167	gp_XYZ dirxyz;
	168
	169	dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir);
	170	dirxyz.Multiply(C.Radius());
	171	pt1.SetXYZ(C.Location().XYZ().Added(dirxyz));
	172
	173	dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir);
	174	dirxyz.Multiply(C.Radius());
	175	pt2.SetXYZ(C.Location().XYZ().Added(dirxyz));
	176	}
	177	else {
	178	nbSol = 0;
	179	}
	180
	181	done = Standard_True;
	182	}
	183
	184	void Contap_ContAna::Perform (const gp_Cylinder& C,
733a0e55	185	const gp_Pnt& Eye)
7fd59977	186	{
	187	done = Standard_False;
	188
	189	Standard_Real radius = C.Radius();
	190	gp_Lin theaxis(C.Axis());
	191	Standard_Real dist = theaxis.Distance(Eye);
	192	if (dist <= radius) {
	193	nbSol = 0;
	194	}
	195	else {
	196	typL = GeomAbs_Line;
	197	prm = radiussqrt(1.-radiusradius/(dist*dist));
	198	dir1 = C.Axis().Direction();
	199	dir2 = dir1;
0d969553	200	gp_XYZ axeye(theaxis.Normal(Eye).Direction().XYZ()); // orientate the axis to the outside
7fd59977	201	gp_XYZ normale((theaxis.Direction().Crossed(axeye)).XYZ());
	202	// normale.Normalize();
	203	pt1.SetXYZ(C.Location().XYZ() + (radiusradius/dist)axeye);
	204	pt2.SetXYZ(pt1.XYZ() - prm*normale);
	205	pt1.SetXYZ(pt1.XYZ() + prm*normale);
	206	nbSol = 2;
	207	}
	208	done = Standard_True;
	209	}
	210
7fd59977	211	void Contap_ContAna::Perform (const gp_Cone& C,
733a0e55	212	const gp_Dir& D)
7fd59977	213	{
	214	done = Standard_False;
	215
	216	Standard_Real Tgtalpha = Tan(C.SemiAngle());
	217
	218	Standard_Real Coefcos = D.Dot(C.Position().XDirection());
	219	Standard_Real Coefsin = D.Dot(C.Position().YDirection());
	220	Standard_Real Coefcst = D.Dot(C.Axis().Direction())*Tgtalpha;
	221
	222	Standard_Real norm1 = CoefcosCoefcos + CoefsinCoefsin;
	223	Standard_Real norm2 = Sqrt(norm1);
	224	// if (Abs(Abs(Coefcst)-norm2) <= Tolpetit) { // tol angulaire 1.e-8
	225	// typL = GeomAbs_Line;
	226	// nbSol = 1;
	227	// pt1 = C.Apex();
	228	// dir1 = D;
	229	// }
	230	// else if (Abs(Coefcst) < norm2) {
	231
	232	if (Abs(Coefcst) < norm2) {
	233	typL = GeomAbs_Line;
	234	nbSol = 2;
	235	pt1 = C.Apex();
	236	pt2 = pt1;
0d969553 Y	237	// Necessary to solve Coefcoscos(t) + Coefsinsin(t) = Coefcst
0d969553 Y	238	// and director vectors of solutions are
7fd59977	239	// cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
0d969553	240	// By setting cos(phi) = Coefcos/Sqrt(Coefcos2 + Coefsin2) and
7fd59977	241	// sin(phi) = Coefsin/Sqrt(Coefcos2 + Coefsin2)
0d969553 Y	242	// and by using trigonometric relations the values of cosinus
0d969553 Y	243	// and sinus to the solutions are obtained.
7fd59977	244
	245	prm = Sqrt(norm1 - Coefcst*Coefcst);
	246	Standard_Real cost0,sint0,cost1,sint1;
	247
	248	cost0 = (CoefcosCoefcst - Coefsinprm)/norm1;
	249	cost1 = (CoefcosCoefcst + Coefsinprm)/norm1;
	250
	251	sint0 = ( Coefcosprm + CoefsinCoefcst)/norm1;
	252	sint1 = (-Coefcosprm + CoefsinCoefcst)/norm1;
	253
	254	gp_XYZ Xdir(C.Position().XDirection().XYZ());
	255	gp_XYZ Ydir(C.Position().YDirection().XYZ());
	256	gp_XYZ Zdir(C.Axis().Direction().XYZ());
	257	gp_XYZ dirxyz;
	258	dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,1./Tgtalpha,Zdir);
	259	dir1.SetXYZ(dirxyz);
	260	pt1.SetXYZ(pt1.XYZ()+dirxyz);
	261	dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,1./Tgtalpha,Zdir);
	262	dir2.SetXYZ(dirxyz);
	263	pt2.SetXYZ(pt2.XYZ()+dirxyz);
	264	}
	265	else {
	266	nbSol = 0;
	267	}
	268	done = Standard_True;
	269	}
	270
7fd59977	271	void Contap_ContAna::Perform (const gp_Cone& C,
733a0e55 S	272	const gp_Dir& D,
733a0e55 S	273	const Standard_Real Angle)
7fd59977	274	{
	275	done = Standard_False;
	276	nbSol = 0;
	277
	278	Standard_Real Ang = C.SemiAngle();
	279	Standard_Real Cosa = cos(Ang);
	280	Standard_Real Sina = sin(Ang);
	281
	282	Standard_Real Coefcos = D.Dot(C.Position().XDirection());
	283	Standard_Real Coefsin = D.Dot(C.Position().YDirection());
	284
c6541a0c	285	Standard_Real Coefcst1 = cos(M_PI*0.5 + Angle);
7fd59977	286
	287	Standard_Real norm1 = CoefcosCoefcos + CoefsinCoefsin;
	288	Standard_Real norm2 = Sqrt(norm1);
	289
	290	Standard_Real Coefnz = D.Dot(C.Axis().Direction())*Sina;
	291	Standard_Real Coefcst = (Coefcst1 + Coefnz)/Cosa;
	292
	293	if (Abs(Coefcst) < norm2) {
	294	typL = GeomAbs_Line;
	295	nbSol+= 2;
	296	pt1 = C.Apex();
	297	pt2 = pt1;
	298
0d969553 Y	299	// It is requiredto solve Coefcoscos(t) + Coefsinsin(t) = Coefcst
0d969553 Y	300	// and the director vectors of solutions are
7fd59977	301	// cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
0d969553	302	// By setting cos(phi) = Coefcos/Sqrt(Coefcos2 + Coefsin2) and
7fd59977	303	// sin(phi) = Coefsin/Sqrt(Coefcos2 + Coefsin2)
0d969553 Y	304	// and by using trigonometric relations the values of cosinus
0d969553 Y	305	// and sinus to the solutions are obtained.
7fd59977	306
	307	prm = Sqrt(norm1 - Coefcst*Coefcst);
	308	Standard_Real cost0,sint0,cost1,sint1;
	309
	310	cost0 = (CoefcosCoefcst - Coefsinprm)/norm1;
	311	cost1 = (CoefcosCoefcst + Coefsinprm)/norm1;
	312
	313	sint0 = ( Coefcosprm + CoefsinCoefcst)/norm1;
	314	sint1 = (-Coefcosprm + CoefsinCoefcst)/norm1;
	315
	316	gp_XYZ Xdir(C.Position().XDirection().XYZ());
	317	gp_XYZ Ydir(C.Position().YDirection().XYZ());
	318	gp_XYZ Zdir(C.Axis().Direction().XYZ());
	319	if (!C.Direct()) {
	320	Zdir.Reverse();
	321	}
	322	gp_XYZ dirxyz;
	323	dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,Cosa/Sina,Zdir);
	324	dir1.SetXYZ(dirxyz);
	325	pt1.SetXYZ(pt1.XYZ()+dirxyz);
	326	dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,Cosa/Sina,Zdir);
	327	dir2.SetXYZ(dirxyz);
	328	pt2.SetXYZ(pt2.XYZ()+dirxyz);
	329	}
	330
	331	Coefcst = (Coefcst1 - Coefnz)/Cosa;
	332
	333	if (Abs(Coefcst) < norm2) {
	334	typL = GeomAbs_Line;
	335	nbSol+= 2;
	336	pt3 = C.Apex();
	337	pt4 = pt3;
	338
	339	prm = Sqrt(norm1 - Coefcst*Coefcst);
	340	Standard_Real cost0,sint0,cost1,sint1;
	341
	342	cost0 = (CoefcosCoefcst - Coefsinprm)/norm1;
	343	cost1 = (CoefcosCoefcst + Coefsinprm)/norm1;
	344
	345	sint0 = ( Coefcosprm + CoefsinCoefcst)/norm1;
	346	sint1 = (-Coefcosprm + CoefsinCoefcst)/norm1;
	347
	348	gp_XYZ Xdir(C.Position().XDirection().XYZ());
	349	gp_XYZ Ydir(C.Position().YDirection().XYZ());
	350	gp_XYZ Zdir(C.Axis().Direction().XYZ());
	351	if (!C.Direct()) {
	352	Zdir.Reverse();
	353	}
	354	gp_XYZ dirxyz;
	355	dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,-Cosa/Sina,Zdir);
	356	dir3.SetXYZ(dirxyz);
	357	pt3.SetXYZ(pt3.XYZ()+dirxyz);
	358	dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,-Cosa/Sina,Zdir);
	359	dir4.SetXYZ(dirxyz);
	360	pt4.SetXYZ(pt4.XYZ()+dirxyz);
	361	if (nbSol == 2) {
	362	pt1 = pt3;
	363	pt2 = pt4;
	364	dir1 = dir3;
	365	dir2 = dir4;
	366	}
	367	}
	368
	369	done = Standard_True;
7fd59977	370	}
	371
	372	void Contap_ContAna::Perform (const gp_Cone& C,
733a0e55	373	const gp_Pnt& Eye)
7fd59977	374	{
	375	done = Standard_False;
	376
	377	Standard_Real Tgtalpha = Tan(C.SemiAngle());
	378
	379	gp_XYZ apexeye(Eye.XYZ());
	380	apexeye.Subtract(C.Apex().XYZ());
	381
	382	Standard_Real Coefcos = apexeye.Dot(C.Position().XDirection().XYZ());
	383	Standard_Real Coefsin = apexeye.Dot(C.Position().YDirection().XYZ());
	384	Standard_Real Coefcst = apexeye.Dot(C.Axis().Direction().XYZ())*Tgtalpha;
	385
	386	Standard_Real norm1 = CoefcosCoefcos + CoefsinCoefsin;
	387	Standard_Real norm2 = Sqrt(CoefcosCoefcos + CoefsinCoefsin);
	388	// if (Abs(Abs(Coefcst)-norm2) <= Tolpetit) { // tol angulaire 1.e-8
	389	// typL = GeomAbs_Line;
	390	// nbSol = 1;
	391	// pt1 = C.Apex();
	392	// dir1.SetXYZ(apexeye);
	393	// }
	394	// else if (Abs(Coefcst) < norm2) {
	395
	396	if (Abs(Coefcst) < norm2) {
	397	typL = GeomAbs_Line;
	398	nbSol = 2;
	399	pt1 = C.Apex();
	400	pt2 = pt1;
0d969553 Y	401	// It is required to solve Coefcoscos(t) + Coefsinsin(t) = Coefcst
0d969553 Y	402	// and the director vectors of solutions are
7fd59977	403	// cos(t0) * XDirection + sin(t0) * YDirection + ZDirection/Tgtalpha
0d969553	404	// By setting cos(phi) = Coefcos/Sqrt(Coefcos2 + Coefsin2) and
7fd59977	405	// sin(phi) = Coefsin/Sqrt(Coefcos2 + Coefsin2)
0d969553 Y	406	// and by using trigonometric relations the values of cosinus
0d969553 Y	407	// and sinus to the solutions are obtained.
7fd59977	408
	409	prm = Sqrt(norm1 - Coefcst*Coefcst);
	410	Standard_Real cost0,sint0,cost1,sint1;
	411
	412	cost0 = (CoefcosCoefcst - Coefsinprm)/norm1;
	413	cost1 = (CoefcosCoefcst + Coefsinprm)/norm1;
	414
	415	sint0 = ( Coefcosprm + CoefsinCoefcst)/norm1;
	416	sint1 = (-Coefcosprm + CoefsinCoefcst)/norm1;
	417
	418	gp_XYZ Xdir(C.Position().XDirection().XYZ());
	419	gp_XYZ Ydir(C.Position().YDirection().XYZ());
	420	gp_XYZ Zdir(C.Axis().Direction().XYZ());
	421	gp_XYZ dirxyz;
	422	dirxyz.SetLinearForm(cost0,Xdir,sint0,Ydir,1./Tgtalpha,Zdir);
	423	dir1.SetXYZ(dirxyz);
	424	pt1.SetXYZ(pt1.XYZ()+dirxyz);
	425	dirxyz.SetLinearForm(cost1,Xdir,sint1,Ydir,1./Tgtalpha,Zdir);
	426	dir2.SetXYZ(dirxyz);
	427	pt2.SetXYZ(pt2.XYZ()+dirxyz);
	428	}
	429	else {
	430	nbSol = 0;
	431	}
	432	done = Standard_True;
7fd59977	433	}
7fd59977	434
7fd59977	435	gp_Lin Contap_ContAna::Line (const Standard_Integer Index) const
	436	{
	437	if (!done) {StdFail_NotDone::Raise();}
	438	if (typL != GeomAbs_Line \|\| nbSol == 0) {Standard_DomainError::Raise();}
	439	if (Index <=0 \|\| Index > nbSol) {Standard_OutOfRange::Raise();}
	440	switch (Index) {
	441	case 1:
	442	return gp_Lin(pt1,dir1);
	443	case 2:
	444	return gp_Lin(pt2,dir2);
	445	case 3:
	446	return gp_Lin(pt3,dir3);
	447	case 4:
	448	return gp_Lin(pt4,dir4);
	449	}
733a0e55	450	Standard_OutOfRange::Raise("Program error in Contap_ContAna");
7fd59977	451	return gp_Lin();
7fd59977	452	}