Skip to content

Commit

Permalink
Projections of arcs of circle, orthogonal case
Browse files Browse the repository at this point in the history
Implementation when circle plane orthogonal to projection plane.
  • Loading branch information
shermelin authored and wwmayer committed Sep 24, 2020
1 parent 3e7e2a9 commit 923607b
Showing 1 changed file with 137 additions and 68 deletions.
205 changes: 137 additions & 68 deletions src/Mod/Sketcher/App/SketchObject.cpp
Expand Up @@ -57,6 +57,7 @@
# include <GC_MakeCircle.hxx>
# include <Standard_Version.hxx>
# include <cmath>
# include <string>
# include <vector>
# include <boost_bind_bind.hpp>
//# include <QtGlobal>
Expand Down Expand Up @@ -5751,11 +5752,11 @@ static gp_Vec ProjVecOnPlane_UVN( const gp_Vec& V, const gp_Pln& Pl)
}

// Auxiliary Method: returns vector projection in XYZ space
/*static gp_Vec ProjVecOnPlane_XYZ( const gp_Vec& V, const gp_Pln& Pl)
static gp_Vec ProjVecOnPlane_XYZ( const gp_Vec& V, const gp_Pln& Pl)
{
return V.Dot(Pl.Position().XDirection()) * Pl.Position().XDirection() +
V.Dot(Pl.Position().YDirection()) * Pl.Position().YDirection();
}*/
}

// Auxiliary Method: returns point projection in UV space of plane
static gp_Vec2d ProjPointOnPlane_UV( const gp_Pnt& P, const gp_Pln& Pl)
Expand All @@ -5765,14 +5766,18 @@ static gp_Vec2d ProjPointOnPlane_UV( const gp_Pnt& P, const gp_Pln& Pl)
}

// Auxiliary Method: returns point projection in UVN space of plane
static gp_Vec ProjPointOnPlane_UVN( const gp_Pnt& P, const gp_Pln& Pl)
static gp_Vec ProjPointOnPlane_UVN(const gp_Pnt& P, const gp_Pln& Pl)
{
gp_Vec2d vec2 = ProjPointOnPlane_UV(P, Pl);
return gp_Vec(vec2.X(), vec2.Y(), 0.0);
}



// Auxiliary Method: returns point projection in XYZ space
static gp_Pnt ProjPointOnPlane_XYZ(const gp_Pnt& P, const gp_Pln& Pl)
{
gp_Vec positionUVN = ProjPointOnPlane_UVN(P, Pl);
return gp_Pnt((positionUVN.X() * Pl.Position().XDirection() + positionUVN.Y() * Pl.Position().YDirection() + gp_Vec(Pl.Location().XYZ())).XYZ());
}

// Auxiliary method
Part::Geometry* projectLine(const BRepAdaptor_Curve& curve, const Handle(Geom_Plane)& gPlane, const Base::Placement& invPlm)
Expand Down Expand Up @@ -6046,69 +6051,133 @@ void SketchObject::rebuildExternalGeometry(void)
}
}
else {
// creates an ellipse or a segment

if (!curve.IsClosed()) {
throw Base::NotImplementedError("Non parallel arcs of circle not yet supported geometry for external geometry");
}
else {
gp_Dir vec1 = sketchPlane.Axis().Direction();
gp_Dir vec2 = curve.Circle().Axis().Direction();
gp_Circ origCircle = curve.Circle();

if (vec1.IsNormal(vec2, Precision::Angular())) { // circle's normal vector in plane:
// projection is a line
// define center by projection
gp_Pnt cnt = origCircle.Location();
GeomAPI_ProjectPointOnSurf proj(cnt,gPlane);
cnt = proj.NearestPoint();
// gp_Dir dirLine(vec1.Crossed(vec2));
gp_Dir dirLine(vec1 ^ vec2);

Part::GeomLineSegment * projectedSegment = new Part::GeomLineSegment();
Geom_Line ligne(cnt, dirLine); // helper object to compute end points
gp_Pnt P1, P2; // end points of the segment, OCC style

ligne.D0(-origCircle.Radius(), P1);
ligne.D0( origCircle.Radius(), P2);

Base::Vector3d p1(P1.X(),P1.Y(),P1.Z()); // ends of segment FCAD style
Base::Vector3d p2(P2.X(),P2.Y(),P2.Z());
invPlm.multVec(p1,p1);
invPlm.multVec(p2,p2);

projectedSegment->setPoints(p1, p2);
projectedSegment->Construction = true;
ExternalGeo.push_back(projectedSegment);
}
else { // general case, full circle
gp_Pnt cnt = origCircle.Location();
GeomAPI_ProjectPointOnSurf proj(cnt,gPlane);
cnt = proj.NearestPoint(); // projection of circle center on sketch plane, 3D space
Base::Vector3d p(cnt.X(),cnt.Y(),cnt.Z()); // converting to FCAD style vector
invPlm.multVec(p,p); // transforming towards sketch's (x,y) coordinates


gp_Vec vecMajorAxis = vec1 ^ vec2; // major axis in 3D space

double minorRadius; // TODO use data type of vectors around...
double cosTheta;
cosTheta = fabs(vec1.Dot(vec2)); // cos of angle between the two planes, assuming vectirs are normalized to 1
minorRadius = origCircle.Radius() * cosTheta;

Base::Vector3d vectorMajorAxis(vecMajorAxis.X(),vecMajorAxis.Y(),vecMajorAxis.Z()); // maj axis into FCAD style vector
invRot.multVec(vectorMajorAxis, vectorMajorAxis); // transforming to sketch's (x,y) coordinates
vecMajorAxis.SetXYZ(gp_XYZ(vectorMajorAxis[0], vectorMajorAxis[1], vectorMajorAxis[2])); // back to OCC

gp_Ax2 refFrameEllipse(gp_Pnt(gp_XYZ(p[0], p[1], p[2])), gp_Vec(0, 0, 1), vecMajorAxis); // NB: force normal of ellipse to be normal of sketch's plane.
Handle(Geom_Ellipse) curve = new Geom_Ellipse(refFrameEllipse, origCircle.Radius(), minorRadius);
Part::GeomEllipse* ellipse = new Part::GeomEllipse();
ellipse->setHandle(curve);
ellipse->Construction = true;

ExternalGeo.push_back(ellipse);
}
}
// creates an ellipse or a segment

gp_Dir vec1 = sketchPlane.Axis().Direction();
gp_Dir vec2 = curve.Circle().Axis().Direction();
gp_Circ origCircle = curve.Circle();

if (vec1.IsNormal(vec2, Precision::Angular())) { // circle's normal vector in plane:
// projection is a line
// define center by projection
gp_Pnt cnt = origCircle.Location();
GeomAPI_ProjectPointOnSurf proj(cnt,gPlane);
cnt = proj.NearestPoint();

gp_Dir dirOrientation = gp_Dir(vec1 ^ vec2);
gp_Dir dirLine(dirOrientation);

Part::GeomLineSegment * projectedSegment = new Part::GeomLineSegment();
Geom_Line ligne(cnt, dirLine); // helper object to compute end points
gp_Pnt P1, P2; // end points of the segment, OCC style

ligne.D0(-origCircle.Radius(), P1);
ligne.D0( origCircle.Radius(), P2);

if (!curve.IsClosed()) { // arc of circle

gp_Pnt pntF = curve.Value(curve.FirstParameter()); // start point of arc of circle
gp_Pnt pntL = curve.Value(curve.LastParameter()); // end point of arc of circle

double alpha = dirOrientation.AngleWithRef(curve.Circle().XAxis().Direction(), curve.Circle().Axis().Direction());

double baseAngle = curve.FirstParameter();

int tours = 0;
double startAngle = baseAngle + alpha;
// bring startAngle back in [-pi/2 , 3pi/2[
while(startAngle < -M_PI/2.0 && tours < 10) {
startAngle = baseAngle + ++tours*2.0*M_PI + alpha;
}
while(startAngle >= 3.0*M_PI/2.0 && tours > -10) {
startAngle = baseAngle + --tours*2.0*M_PI + alpha;
}

double endAngle = curve.LastParameter() + startAngle - baseAngle; // apply same offset to end angle

if (startAngle <= 0.0) {
if (endAngle <= 0.0) {
P1 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
P2 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
} else {
if (endAngle <= fabs(startAngle)) {
// P2 = P2 already defined
P1 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
} else if (endAngle < M_PI) {
// P2 = P2, already defined
P1 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
} else {
// P1 = P1, already defined
// P2 = P2, already defined
}
}
} else if (startAngle < M_PI) {
if (endAngle < M_PI) {
P1 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
P2 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
} else if (endAngle < 2.0 * M_PI - startAngle) {
P2 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
// P1 = P1, already defined
} else if (endAngle < 2.0 * M_PI) {
P2 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
// P1 = P1, already defined
} else {
// P1 = P1, alreday defined
// P2 = P2, already defined
}
} else {
if (endAngle < 2*M_PI) {
P1 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
P2 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
} else if (endAngle < 4*M_PI - startAngle) {
P1 = ProjPointOnPlane_XYZ(pntF, sketchPlane);
// P2 = P2, already defined
} else if (endAngle < 3*M_PI) {
// P1 = P1, already defined
P2 = ProjPointOnPlane_XYZ(pntL, sketchPlane);
} else {
// P1 = P1, already defined
// P2 = P2, already defined
}
}
}

Base::Vector3d p1(P1.X(),P1.Y(),P1.Z()); // ends of segment FCAD style
Base::Vector3d p2(P2.X(),P2.Y(),P2.Z());
invPlm.multVec(p1,p1);
invPlm.multVec(p2,p2);

projectedSegment->setPoints(p1, p2);
projectedSegment->Construction = true;
ExternalGeo.push_back(projectedSegment);
}
else { // general case, full circle
gp_Pnt cnt = origCircle.Location();
GeomAPI_ProjectPointOnSurf proj(cnt,gPlane);
cnt = proj.NearestPoint(); // projection of circle center on sketch plane, 3D space
Base::Vector3d p(cnt.X(),cnt.Y(),cnt.Z()); // converting to FCAD style vector
invPlm.multVec(p,p); // transforming towards sketch's (x,y) coordinates


gp_Vec vecMajorAxis = vec1 ^ vec2; // major axis in 3D space

double minorRadius; // TODO use data type of vectors around...
double cosTheta;
cosTheta = fabs(vec1.Dot(vec2)); // cos of angle between the two planes, assuming vectirs are normalized to 1
minorRadius = origCircle.Radius() * cosTheta;

Base::Vector3d vectorMajorAxis(vecMajorAxis.X(),vecMajorAxis.Y(),vecMajorAxis.Z()); // maj axis into FCAD style vector
invRot.multVec(vectorMajorAxis, vectorMajorAxis); // transforming to sketch's (x,y) coordinates
vecMajorAxis.SetXYZ(gp_XYZ(vectorMajorAxis[0], vectorMajorAxis[1], vectorMajorAxis[2])); // back to OCC

gp_Ax2 refFrameEllipse(gp_Pnt(gp_XYZ(p[0], p[1], p[2])), gp_Vec(0, 0, 1), vecMajorAxis); // NB: force normal of ellipse to be normal of sketch's plane.
Handle(Geom_Ellipse) curve = new Geom_Ellipse(refFrameEllipse, origCircle.Radius(), minorRadius);
Part::GeomEllipse* ellipse = new Part::GeomEllipse();
ellipse->setHandle(curve);
ellipse->Construction = true;

ExternalGeo.push_back(ellipse);
}
}
}
else if (curve.GetType() == GeomAbs_Ellipse) {
Expand Down

0 comments on commit 923607b

Please sign in to comment.