diff --git a/src/Mod/Draft/DraftGeomUtils.py b/src/Mod/Draft/DraftGeomUtils.py index 4e4dc68ee..9267d61aa 100755 --- a/src/Mod/Draft/DraftGeomUtils.py +++ b/src/Mod/Draft/DraftGeomUtils.py @@ -51,10 +51,10 @@ def vec(edge): return None def edg(p1,p2): - "edg(Vector,Vector): returns an edge from 2 vectors" - if isinstance(p1,FreeCAD.Vector) and isinstance(p2,FreeCAD.Vector): - if DraftVecUtils.equals(p1,p2): return None - else: return Part.Line(p1,p2).toShape() + "edg(Vector,Vector): returns an edge from 2 vectors" + if isinstance(p1,FreeCAD.Vector) and isinstance(p2,FreeCAD.Vector): + if DraftVecUtils.equals(p1,p2): return None + else: return Part.Line(p1,p2).toShape() def getVerts(shape): "getVerts(shape): returns a list containing vectors of each vertex of the shape" @@ -64,8 +64,8 @@ def getVerts(shape): return p def v1(edge): - "v1(edge): returns the first point of an edge" - return edge.Vertexes[0].Point + "v1(edge): returns the first point of an edge" + return edge.Vertexes[0].Point def isNull(something): '''isNull(object): returns true if the given shape is null or the given placement is null or @@ -84,38 +84,39 @@ def isNull(something): return False def isPtOnEdge(pt,edge) : - '''isPtOnEdge(Vector,edge): Tests if a point is on an edge''' - if isinstance(edge.Curve,Part.Line) : - orig = edge.Vertexes[0].Point - if DraftVecUtils.isNull(pt.sub(orig).cross(vec(edge))) : - return pt.sub(orig).Length <= vec(edge).Length and pt.sub(orig).dot(vec(edge)) >= 0 - else : - return False - - elif isinstance(edge.Curve,Part.Circle) : - center = edge.Curve.Center - axis = edge.Curve.Axis ; axis.normalize() - radius = edge.Curve.Radius - if round(pt.sub(center).dot(axis),precision()) == 0 \ - and round(pt.sub(center).Length - radius,precision()) == 0 : - if len(edge.Vertexes) == 1 : - return True # edge is a complete circle - else : - begin = edge.Vertexes[0].Point - end = edge.Vertexes[-1].Point - if DraftVecUtils.isNull(pt.sub(begin)) or DraftVecUtils.isNull(pt.sub(end)) : - return True - else : - # newArc = Part.Arc(begin,pt,end) - # return DraftVecUtils.isNull(newArc.Center.sub(center)) \ - # and DraftVecUtils.isNull(newArc.Axis-axis) \ - # and round(newArc.Radius-radius,precision()) == 0 - angle1 = DraftVecUtils.angle(begin.sub(center)) - angle2 = DraftVecUtils.angle(end.sub(center)) - anglept = DraftVecUtils.angle(pt.sub(center)) - if (angle1 < anglept) and (anglept < angle2): - return True - return False + '''isPtOnEdge(Vector,edge): Tests if a point is on an edge''' + if isinstance(edge.Curve,Part.Line) : + orig = edge.Vertexes[0].Point + #if DraftVecUtils.isNull(pt.sub(orig).cross(vec(edge))) : <== can give unprecise results + if round(pt.sub(orig).getAngle(vec(edge)),precision()) == 0: + return pt.sub(orig).Length <= vec(edge).Length and pt.sub(orig).dot(vec(edge)) >= 0 + else : + return False + + elif isinstance(edge.Curve,Part.Circle) : + center = edge.Curve.Center + axis = edge.Curve.Axis ; axis.normalize() + radius = edge.Curve.Radius + if round(pt.sub(center).dot(axis),precision()) == 0 \ + and round(pt.sub(center).Length - radius,precision()) == 0 : + if len(edge.Vertexes) == 1 : + return True # edge is a complete circle + else : + begin = edge.Vertexes[0].Point + end = edge.Vertexes[-1].Point + if DraftVecUtils.isNull(pt.sub(begin)) or DraftVecUtils.isNull(pt.sub(end)) : + return True + else : + # newArc = Part.Arc(begin,pt,end) + # return DraftVecUtils.isNull(newArc.Center.sub(center)) \ + # and DraftVecUtils.isNull(newArc.Axis-axis) \ + # and round(newArc.Radius-radius,precision()) == 0 + angle1 = DraftVecUtils.angle(begin.sub(center)) + angle2 = DraftVecUtils.angle(end.sub(center)) + anglept = DraftVecUtils.angle(pt.sub(center)) + if (angle1 < anglept) and (anglept < angle2): + return True + return False def hasCurves(shape): "hasCurve(shape): checks if the given shape has curves" @@ -171,27 +172,27 @@ def findEdge(anEdge,aList): if DraftVecUtils.equals(anEdge.Vertexes[-1].Point,aList[e].Vertexes[-1].Point): return(e) return None - + def findIntersection(edge1,edge2,infinite1=False,infinite2=False,ex1=False,ex2=False) : - '''findIntersection(edge1,edge2,infinite1=False,infinite2=False): - returns a list containing the intersection point(s) of 2 edges. - You can also feed 4 points instead of edge1 and edge2''' + '''findIntersection(edge1,edge2,infinite1=False,infinite2=False): + returns a list containing the intersection point(s) of 2 edges. + You can also feed 4 points instead of edge1 and edge2''' - pt1 = None - - if isinstance(edge1,FreeCAD.Vector) and isinstance(edge2,FreeCAD.Vector): - # we got points directly - pt1 = edge1 - pt2 = edge2 - pt3 = infinite1 - pt4 = infinite2 - infinite1 = ex1 - infinite2 = ex2 + pt1 = None - elif isinstance(edge1.Curve,Part.Line) and isinstance(edge2.Curve,Part.Line) : - # we have 2 edges - pt1, pt2, pt3, pt4 = [edge1.Vertexes[0].Point, + if isinstance(edge1,FreeCAD.Vector) and isinstance(edge2,FreeCAD.Vector): + # we got points directly + pt1 = edge1 + pt2 = edge2 + pt3 = infinite1 + pt4 = infinite2 + infinite1 = ex1 + infinite2 = ex2 + + elif isinstance(edge1.Curve,Part.Line) and isinstance(edge2.Curve,Part.Line) : + # we have 2 straight lines + pt1, pt2, pt3, pt4 = [edge1.Vertexes[0].Point, edge1.Vertexes[1].Point, edge2.Vertexes[0].Point, edge2.Vertexes[1].Point] @@ -202,160 +203,164 @@ def findIntersection(edge1,edge2,infinite1=False,infinite2=False,ex1=False,ex2=F return [pt1] elif (pt2 in [pt3,pt4]): return [pt2] + norm1 = pt2.sub(pt1).cross(pt3.sub(pt1)) + norm2 = pt2.sub(pt4).cross(pt3.sub(pt4)) + norm1.normalize() + norm2.normalize() + if DraftVecUtils.isNull(norm1.cross(norm2)): + vec1 = pt2.sub(pt1) + vec2 = pt4.sub(pt3) + if DraftVecUtils.isNull(vec1) or DraftVecUtils.isNull(vec2): + return [] # One of the line has zero-length + vec1.normalize() + vec2.normalize() + norm3 = vec1.cross(vec2) + if not DraftVecUtils.isNull(norm3) : + k = ((pt3.z-pt1.z)*(vec2.x-vec2.y)+(pt3.y-pt1.y)*(vec2.z-vec2.x)+ \ + (pt3.x-pt1.x)*(vec2.y-vec2.z))/(norm3.x+norm3.y+norm3.z) + vec1.scale(k,k,k) + intp = pt1.add(vec1) - #we have 2 straight lines - if DraftVecUtils.isNull(pt2.sub(pt1).cross(pt3.sub(pt1)).cross(pt2.sub(pt4).cross(pt3.sub(pt4)))): - vec1 = pt2.sub(pt1) ; vec2 = pt4.sub(pt3) - if DraftVecUtils.isNull(vec1) or DraftVecUtils.isNull(vec2): - return [] - vec1.normalize() ; vec2.normalize() - cross = vec1.cross(vec2) - if not DraftVecUtils.isNull(cross) : - k = ((pt3.z-pt1.z)*(vec2.x-vec2.y)+(pt3.y-pt1.y)*(vec2.z-vec2.x)+ \ - (pt3.x-pt1.x)*(vec2.y-vec2.z))/(cross.x+cross.y+cross.z) - vec1.scale(k,k,k) - int = pt1.add(vec1) - - if infinite1 == False and not isPtOnEdge(int,edge1) : - return [] - - if infinite2 == False and not isPtOnEdge(int,edge2) : - return [] - - return [int] - else : - return [] # Lines have same direction - else : - return [] # Lines aren't on same plane - - elif isinstance(edge1.Curve,Part.Circle) and isinstance(edge2.Curve,Part.Line) \ - or isinstance(edge1.Curve,Part.Line) and isinstance(edge2.Curve,Part.Circle) : - - # deals with an arc or circle and a line - - edges = [edge1,edge2] - for edge in edges : - if isinstance(edge.Curve,Part.Line) : - line = edge - else : - arc = edge - - dirVec = vec(line) ; dirVec.normalize() - pt1 = line.Vertexes[0].Point - pt2 = line.Vertexes[1].Point - pt3 = arc.Vertexes[0].Point - pt4 = arc.Vertexes[-1].Point - center = arc.Curve.Center - - # first check for coincident endpoints - if (pt1 in [pt3,pt4]): - return [pt1] - elif (pt2 in [pt3,pt4]): - return [pt2] + if infinite1 == False and not isPtOnEdge(intp,edge1) : + return [] + + if infinite2 == False and not isPtOnEdge(intp,edge2) : + return [] + + return [intp] + else : + return [] # Lines have same direction + else : + return [] # Lines aren't on same plane + + elif isinstance(edge1.Curve,Part.Circle) and isinstance(edge2.Curve,Part.Line) \ + or isinstance(edge1.Curve,Part.Line) and isinstance(edge2.Curve,Part.Circle) : + + # deals with an arc or circle and a line + + edges = [edge1,edge2] + for edge in edges : + if isinstance(edge.Curve,Part.Line) : + line = edge + else : + arc = edge - if DraftVecUtils.isNull(pt1.sub(center).cross(pt2.sub(center)).cross(arc.Curve.Axis)) : - # Line and Arc are on same plane - - dOnLine = center.sub(pt1).dot(dirVec) - onLine = Vector(dirVec) ; onLine.scale(dOnLine,dOnLine,dOnLine) - toLine = pt1.sub(center).add(onLine) - - if toLine.Length < arc.Curve.Radius : - dOnLine = (arc.Curve.Radius**2 - toLine.Length**2)**(0.5) - onLine = Vector(dirVec) ; onLine.scale(dOnLine,dOnLine,dOnLine) - int = [center.add(toLine).add(onLine)] - onLine = Vector(dirVec) ; onLine.scale(-dOnLine,-dOnLine,-dOnLine) - int += [center.add(toLine).add(onLine)] - elif round(toLine.Length-arc.Curve.Radius,precision()) == 0 : - int = [center.add(toLine)] - else : - return [] - - else : # Line isn't on Arc's plane - if dirVec.dot(arc.Curve.Axis) != 0 : - toPlane = Vector(arc.Curve.Axis) ; toPlane.normalize() - d = pt1.dot(toPlane) - dToPlane = center.sub(pt1).dot(toPlane) - toPlane = Vector(pt1) - toPlane.scale(dToPlane/d,dToPlane/d,dToPlane/d) - ptOnPlane = toPlane.add(pt1) - if round(ptOnPlane.sub(center).Length - arc.Curve.Radius,precision()) == 0 : - int = [ptOnPlane] - else : - return [] - else : - return [] - - if infinite1 == False : - for i in range(len(int)-1,-1,-1) : - if not isPtOnEdge(int[i],edge1) : - del int[i] - if infinite2 == False : - for i in range(len(int)-1,-1,-1) : - if not isPtOnEdge(int[i],edge2) : - del int[i] - return int - - elif isinstance(edge1.Curve,Part.Circle) and isinstance(edge2.Curve,Part.Circle) : - - # deals with 2 arcs or circles + dirVec = vec(line) ; dirVec.normalize() + pt1 = line.Vertexes[0].Point + pt2 = line.Vertexes[1].Point + pt3 = arc.Vertexes[0].Point + pt4 = arc.Vertexes[-1].Point + center = arc.Curve.Center - cent1, cent2 = edge1.Curve.Center, edge2.Curve.Center - rad1 , rad2 = edge1.Curve.Radius, edge2.Curve.Radius - axis1, axis2 = edge1.Curve.Axis , edge2.Curve.Axis - c2c = cent2.sub(cent1) - - if DraftVecUtils.isNull(axis1.cross(axis2)) : - if round(c2c.dot(axis1),precision()) == 0 : - # circles are on same plane - dc2c = c2c.Length ; - if not DraftVecUtils.isNull(c2c): c2c.normalize() - if round(rad1+rad2-dc2c,precision()) < 0 \ - or round(rad1-dc2c-rad2,precision()) > 0 or round(rad2-dc2c-rad1,precision()) > 0 : - return [] - else : - norm = c2c.cross(axis1) - if not DraftVecUtils.isNull(norm): norm.normalize() - if DraftVecUtils.isNull(norm): x = 0 - else: x = (dc2c**2 + rad1**2 - rad2**2)/(2*dc2c) - y = abs(rad1**2 - x**2)**(0.5) - c2c.scale(x,x,x) - if round(y,precision()) != 0 : - norm.scale(y,y,y) - int = [cent1.add(c2c).add(norm)] - int += [cent1.add(c2c).sub(norm)] - else : - int = [cent1.add(c2c)] - else : - return [] # circles are on parallel planes - else : - # circles aren't on same plane - axis1.normalize() ; axis2.normalize() - U = axis1.cross(axis2) - V = axis1.cross(U) - dToPlane = c2c.dot(axis2) - d = V.add(cent1).dot(axis2) - V.scale(dToPlane/d,dToPlane/d,dToPlane/d) - PtOn2Planes = V.add(cent1) - planeIntersectionVector = U.add(PtOn2Planes) - intTemp = findIntersection(planeIntersectionVector,edge1,True,True) - int = [] - for pt in intTemp : - if round(pt.sub(cent2).Length-rad2,precision()) == 0 : - int += [pt] - - if infinite1 == False : - for i in range(len(int)-1,-1,-1) : - if not isPtOnEdge(int[i],edge1) : - del int[i] - if infinite2 == False : - for i in range(len(int)-1,-1,-1) : - if not isPtOnEdge(int[i],edge2) : - del int[i] - - return int - else : - print "DraftGeomUtils: Unsupported curve type: (" + str(edge1.Curve) + ", " + str(edge2.Curve) + ")" + # first check for coincident endpoints + if (pt1 in [pt3,pt4]): + return [pt1] + elif (pt2 in [pt3,pt4]): + return [pt2] + + if DraftVecUtils.isNull(pt1.sub(center).cross(pt2.sub(center)).cross(arc.Curve.Axis)) : + # Line and Arc are on same plane + + dOnLine = center.sub(pt1).dot(dirVec) + onLine = Vector(dirVec) ; onLine.scale(dOnLine,dOnLine,dOnLine) + toLine = pt1.sub(center).add(onLine) + + if toLine.Length < arc.Curve.Radius : + dOnLine = (arc.Curve.Radius**2 - toLine.Length**2)**(0.5) + onLine = Vector(dirVec) ; onLine.scale(dOnLine,dOnLine,dOnLine) + int = [center.add(toLine).add(onLine)] + onLine = Vector(dirVec) ; onLine.scale(-dOnLine,-dOnLine,-dOnLine) + int += [center.add(toLine).add(onLine)] + elif round(toLine.Length-arc.Curve.Radius,precision()) == 0 : + int = [center.add(toLine)] + else : + return [] + + else : # Line isn't on Arc's plane + if dirVec.dot(arc.Curve.Axis) != 0 : + toPlane = Vector(arc.Curve.Axis) ; toPlane.normalize() + d = pt1.dot(toPlane) + dToPlane = center.sub(pt1).dot(toPlane) + toPlane = Vector(pt1) + toPlane.scale(dToPlane/d,dToPlane/d,dToPlane/d) + ptOnPlane = toPlane.add(pt1) + if round(ptOnPlane.sub(center).Length - arc.Curve.Radius,precision()) == 0 : + int = [ptOnPlane] + else : + return [] + else : + return [] + + if infinite1 == False : + for i in range(len(int)-1,-1,-1) : + if not isPtOnEdge(int[i],edge1) : + del int[i] + if infinite2 == False : + for i in range(len(int)-1,-1,-1) : + if not isPtOnEdge(int[i],edge2) : + del int[i] + return int + + elif isinstance(edge1.Curve,Part.Circle) and isinstance(edge2.Curve,Part.Circle) : + + # deals with 2 arcs or circles + + cent1, cent2 = edge1.Curve.Center, edge2.Curve.Center + rad1 , rad2 = edge1.Curve.Radius, edge2.Curve.Radius + axis1, axis2 = edge1.Curve.Axis , edge2.Curve.Axis + c2c = cent2.sub(cent1) + + if DraftVecUtils.isNull(axis1.cross(axis2)) : + if round(c2c.dot(axis1),precision()) == 0 : + # circles are on same plane + dc2c = c2c.Length ; + if not DraftVecUtils.isNull(c2c): c2c.normalize() + if round(rad1+rad2-dc2c,precision()) < 0 \ + or round(rad1-dc2c-rad2,precision()) > 0 or round(rad2-dc2c-rad1,precision()) > 0 : + return [] + else : + norm = c2c.cross(axis1) + if not DraftVecUtils.isNull(norm): norm.normalize() + if DraftVecUtils.isNull(norm): x = 0 + else: x = (dc2c**2 + rad1**2 - rad2**2)/(2*dc2c) + y = abs(rad1**2 - x**2)**(0.5) + c2c.scale(x,x,x) + if round(y,precision()) != 0 : + norm.scale(y,y,y) + int = [cent1.add(c2c).add(norm)] + int += [cent1.add(c2c).sub(norm)] + else : + int = [cent1.add(c2c)] + else : + return [] # circles are on parallel planes + else : + # circles aren't on same plane + axis1.normalize() ; axis2.normalize() + U = axis1.cross(axis2) + V = axis1.cross(U) + dToPlane = c2c.dot(axis2) + d = V.add(cent1).dot(axis2) + V.scale(dToPlane/d,dToPlane/d,dToPlane/d) + PtOn2Planes = V.add(cent1) + planeIntersectionVector = U.add(PtOn2Planes) + intTemp = findIntersection(planeIntersectionVector,edge1,True,True) + int = [] + for pt in intTemp : + if round(pt.sub(cent2).Length-rad2,precision()) == 0 : + int += [pt] + + if infinite1 == False : + for i in range(len(int)-1,-1,-1) : + if not isPtOnEdge(int[i],edge1) : + del int[i] + if infinite2 == False : + for i in range(len(int)-1,-1,-1) : + if not isPtOnEdge(int[i],edge2) : + del int[i] + + return int + else : + print "DraftGeomUtils: Unsupported curve type: (" + str(edge1.Curve) + ", " + str(edge2.Curve) + ")" def geom(edge,plac=FreeCAD.Placement()): "returns a Line, ArcOfCircle or Circle geom from the given edge, according to the given placement" @@ -386,62 +391,62 @@ def geom(edge,plac=FreeCAD.Placement()): return edge.Curve def mirror (point, edge): - "finds mirror point relative to an edge" - normPoint = point.add(findDistance(point, edge, False)) - if normPoint: - normPoint_point = Vector.sub(point, normPoint) - normPoint_refl = DraftVecUtils.neg(normPoint_point) - refl = Vector.add(normPoint, normPoint_refl) - return refl - else: - return None + "finds mirror point relative to an edge" + normPoint = point.add(findDistance(point, edge, False)) + if normPoint: + normPoint_point = Vector.sub(point, normPoint) + normPoint_refl = DraftVecUtils.neg(normPoint_point) + refl = Vector.add(normPoint, normPoint_refl) + return refl + else: + return None def findClosest(basepoint,pointslist): - ''' - findClosest(vector,list) - in a list of 3d points, finds the closest point to the base point. - an index from the list is returned. - ''' - if not pointslist: return None - smallest = 100000 - for n in range(len(pointslist)): - new = basepoint.sub(pointslist[n]).Length - if new < smallest: - smallest = new - npoint = n - return npoint + ''' + findClosest(vector,list) + in a list of 3d points, finds the closest point to the base point. + an index from the list is returned. + ''' + if not pointslist: return None + smallest = 100000 + for n in range(len(pointslist)): + new = basepoint.sub(pointslist[n]).Length + if new < smallest: + smallest = new + npoint = n + return npoint def concatenate(shape): - "concatenate(shape) -- turns several faces into one" - edges = getBoundary(shape) - edges = sortEdges(edges) - try: - wire=Part.Wire(edges) - face=Part.Face(wire) - except: - print "DraftGeomUtils: Couldn't join faces into one" - return(shape) - else: - if not wire.isClosed(): return(wire) - else: return(face) + "concatenate(shape) -- turns several faces into one" + edges = getBoundary(shape) + edges = sortEdges(edges) + try: + wire=Part.Wire(edges) + face=Part.Face(wire) + except: + print "DraftGeomUtils: Couldn't join faces into one" + return(shape) + else: + if not wire.isClosed(): return(wire) + else: return(face) def getBoundary(shape): - "getBoundary(shape) -- this function returns the boundary edges of a group of faces" - # make a lookup-table where we get the number of occurrences - # to each edge in the fused face - if isinstance(shape,list): - shape = Part.makeCompound(shape) - lut={} - for f in shape.Faces: - for e in f.Edges: - hc= e.hashCode() - if lut.has_key(hc): lut[hc]=lut[hc]+1 - else: lut[hc]=1 - # filter out the edges shared by more than one sub-face - bound=[] - for e in shape.Edges: - if lut[e.hashCode()] == 1: bound.append(e) - return bound + "getBoundary(shape) -- this function returns the boundary edges of a group of faces" + # make a lookup-table where we get the number of occurrences + # to each edge in the fused face + if isinstance(shape,list): + shape = Part.makeCompound(shape) + lut={} + for f in shape.Faces: + for e in f.Edges: + hc= e.hashCode() + if lut.has_key(hc): lut[hc]=lut[hc]+1 + else: lut[hc]=1 + # filter out the edges shared by more than one sub-face + bound=[] + for e in shape.Edges: + if lut[e.hashCode()] == 1: bound.append(e) + return bound def isLine(bsp): "returns True if the given BSpline curve is a straight line" @@ -653,112 +658,112 @@ def superWire(edgeslist,closed=False): return Part.Wire(newedges) def findMidpoint(edge): - "calculates the midpoint of an edge" - first = edge.Vertexes[0].Point - last = edge.Vertexes[-1].Point - if isinstance(edge.Curve,Part.Circle): - center = edge.Curve.Center - radius = edge.Curve.Radius - if len(edge.Vertexes) == 1: - # Circle - dv = first.sub(center) - dv = DraftVecUtils.neg(dv) - return center.add(dv) - axis = edge.Curve.Axis - chord = last.sub(first) - perp = chord.cross(axis) - perp.normalize() - ray = first.sub(center) - apothem = ray.dot(perp) - sagitta = radius - apothem - startpoint = Vector.add(first, DraftVecUtils.scale(chord,0.5)) - endpoint = DraftVecUtils.scaleTo(perp,sagitta) - return Vector.add(startpoint,endpoint) + "calculates the midpoint of an edge" + first = edge.Vertexes[0].Point + last = edge.Vertexes[-1].Point + if isinstance(edge.Curve,Part.Circle): + center = edge.Curve.Center + radius = edge.Curve.Radius + if len(edge.Vertexes) == 1: + # Circle + dv = first.sub(center) + dv = DraftVecUtils.neg(dv) + return center.add(dv) + axis = edge.Curve.Axis + chord = last.sub(first) + perp = chord.cross(axis) + perp.normalize() + ray = first.sub(center) + apothem = ray.dot(perp) + sagitta = radius - apothem + startpoint = Vector.add(first, DraftVecUtils.scale(chord,0.5)) + endpoint = DraftVecUtils.scaleTo(perp,sagitta) + return Vector.add(startpoint,endpoint) - elif isinstance(edge.Curve,Part.Line): - halfedge = DraftVecUtils.scale(last.sub(first),.5) - return Vector.add(first,halfedge) + elif isinstance(edge.Curve,Part.Line): + halfedge = DraftVecUtils.scale(last.sub(first),.5) + return Vector.add(first,halfedge) - else: - return None + else: + return None def complexity(obj): - ''' - tests given object for shape complexity: - 1: line - 2: arc - 3: circle - 4: open wire with no arc - 5: closed wire - 6: wire with arcs - 7: faces - 8: faces with arcs - ''' - shape = obj.Shape - if shape.Faces: - for e in shape.Edges: - if (isinstance(e.Curve,Part.Circle)): return 8 - return 7 - if shape.Wires: - for e in shape.Edges: - if (isinstance(e.Curve,Part.Circle)): return 6 - for w in shape.Wires: - if w.isClosed(): return 5 - return 4 - if (isinstance(shape.Edges[0].Curve,Part.Circle)): - if len(shape.Vertexes) == 1: - return 3 - return 2 - return 1 + ''' + tests given object for shape complexity: + 1: line + 2: arc + 3: circle + 4: open wire with no arc + 5: closed wire + 6: wire with arcs + 7: faces + 8: faces with arcs + ''' + shape = obj.Shape + if shape.Faces: + for e in shape.Edges: + if (isinstance(e.Curve,Part.Circle)): return 8 + return 7 + if shape.Wires: + for e in shape.Edges: + if (isinstance(e.Curve,Part.Circle)): return 6 + for w in shape.Wires: + if w.isClosed(): return 5 + return 4 + if (isinstance(shape.Edges[0].Curve,Part.Circle)): + if len(shape.Vertexes) == 1: + return 3 + return 2 + return 1 def findPerpendicular(point,edgeslist,force=None): - ''' - findPerpendicular(vector,wire,[force]): - finds the shortest perpendicular distance between a point and an edgeslist. - If force is specified, only the edge[force] will be considered, and it will be - considered infinite. - The function will return a list [vector_from_point_to_closest_edge,edge_index] - or None if no perpendicular vector could be found. - ''' - if not isinstance(edgeslist,list): - try: - edgeslist = edgeslist.Edges - except: - return None - if (force == None): - valid = None - for edge in edgeslist: - dist = findDistance(point,edge,strict=True) - if dist: - if not valid: valid = [dist,edgeslist.index(edge)] - else: - if (dist.Length < valid[0].Length): - valid = [dist,edgeslist.index(edge)] - return valid - else: - edge = edgeslist[force] - dist = findDistance(point,edge) - if dist: return [dist,force] - else: return None + ''' + findPerpendicular(vector,wire,[force]): + finds the shortest perpendicular distance between a point and an edgeslist. + If force is specified, only the edge[force] will be considered, and it will be + considered infinite. + The function will return a list [vector_from_point_to_closest_edge,edge_index] + or None if no perpendicular vector could be found. + ''' + if not isinstance(edgeslist,list): + try: + edgeslist = edgeslist.Edges + except: + return None + if (force == None): + valid = None + for edge in edgeslist: + dist = findDistance(point,edge,strict=True) + if dist: + if not valid: valid = [dist,edgeslist.index(edge)] + else: + if (dist.Length < valid[0].Length): + valid = [dist,edgeslist.index(edge)] + return valid + else: + edge = edgeslist[force] + dist = findDistance(point,edge) + if dist: return [dist,force] + else: return None return None def offset(edge,vector): - ''' - offset(edge,vector) - returns a copy of the edge at a certain (vector) distance - if the edge is an arc, the vector will be added at its first point - and a complete circle will be returned - ''' - if (not isinstance(edge,Part.Shape)) or (not isinstance(vector,FreeCAD.Vector)): - return None - if isinstance(edge.Curve,Part.Line): - v1 = Vector.add(edge.Vertexes[0].Point, vector) - v2 = Vector.add(edge.Vertexes[-1].Point, vector) - return Part.Line(v1,v2).toShape() - else: - rad = edge.Vertexes[0].Point.sub(edge.Curve.Center) - newrad = Vector.add(rad,vector).Length - return Part.Circle(edge.Curve.Center,NORM,newrad).toShape() + ''' + offset(edge,vector) + returns a copy of the edge at a certain (vector) distance + if the edge is an arc, the vector will be added at its first point + and a complete circle will be returned + ''' + if (not isinstance(edge,Part.Shape)) or (not isinstance(vector,FreeCAD.Vector)): + return None + if isinstance(edge.Curve,Part.Line): + v1 = Vector.add(edge.Vertexes[0].Point, vector) + v2 = Vector.add(edge.Vertexes[-1].Point, vector) + return Part.Line(v1,v2).toShape() + else: + rad = edge.Vertexes[0].Point.sub(edge.Curve.Center) + newrad = Vector.add(rad,vector).Length + return Part.Circle(edge.Curve.Center,NORM,newrad).toShape() def isReallyClosed(wire): "checks if a wire is really closed" @@ -788,9 +793,9 @@ def getNormal(shape): n = e1.cross(e2).normalize() break if FreeCAD.GuiUp: - import Draft - vdir = Draft.get3DView().getViewDirection() - if n.getAngle(vdir) < 0.78: n = DraftVecUtils.neg(n) + import Draft + vdir = Draft.get3DView().getViewDirection() + if n.getAngle(vdir) < 0.78: n = DraftVecUtils.neg(n) return n def getRotation(v1,v2=FreeCAD.Vector(0,0,1)): @@ -833,7 +838,7 @@ def offsetWire(wire,dvec,bind=False,occ=False): l=abs(dvec.Length) if not l: return None if wire.Wires: - wire = wire.Wires[0] + wire = wire.Wires[0] else: wire = Part.Wire(edges) try: @@ -910,104 +915,104 @@ def connect(edges,closed=False): return None def findDistance(point,edge,strict=False): - ''' - findDistance(vector,edge,[strict]) - Returns a vector from the point to its - closest point on the edge. If strict is True, the vector will be returned - only if its endpoint lies on the edge. - ''' - if isinstance(point, FreeCAD.Vector): - if isinstance(edge.Curve, Part.Line): - segment = vec(edge) - chord = edge.Vertexes[0].Point.sub(point) - norm = segment.cross(chord) - perp = segment.cross(norm) - dist = DraftVecUtils.project(chord,perp) - if not dist: return None - newpoint = point.add(dist) - if (dist.Length == 0): - return None - if strict: - s1 = newpoint.sub(edge.Vertexes[0].Point) - s2 = newpoint.sub(edge.Vertexes[-1].Point) - if (s1.Length <= segment.Length) and (s2.Length <= segment.Length): - return dist - else: - return None - else: return dist - elif isinstance(edge.Curve, Part.Circle): - ve1 = edge.Vertexes[0].Point - if (len(edge.Vertexes) > 1): - ve2 = edge.Vertexes[-1].Point - else: - ve2 = None - center = edge.Curve.Center - segment = center.sub(point) - ratio = (segment.Length - edge.Curve.Radius) / segment.Length - dist = DraftVecUtils.scale(segment,ratio) - newpoint = Vector.add(point, dist) - if (dist.Length == 0): - return None - if strict and ve2: - ang1 = DraftVecUtils.angle(ve1.sub(center)) - ang2 = DraftVecUtils.angle(ve2.sub(center)) - angpt = DraftVecUtils.angle(newpoint.sub(center)) - if ((angpt <= ang2 and angpt >= ang1) or (angpt <= ang1 and angpt >= ang2)): - return dist - else: - return None - else: - return dist - elif isinstance(edge.Curve,Part.BSplineCurve): - try: - pr = edge.Curve.parameter(point) - np = edge.Curve.value(pr) - dist = np.sub(point) - except: - print "DraftGeomUtils: Unable to get curve parameter for point ",point - return None - else: - return dist - else: - print "DraftGeomUtils: Couldn't project point" - return None - else: - print "DraftGeomUtils: Couldn't project point" - return None + ''' + findDistance(vector,edge,[strict]) - Returns a vector from the point to its + closest point on the edge. If strict is True, the vector will be returned + only if its endpoint lies on the edge. + ''' + if isinstance(point, FreeCAD.Vector): + if isinstance(edge.Curve, Part.Line): + segment = vec(edge) + chord = edge.Vertexes[0].Point.sub(point) + norm = segment.cross(chord) + perp = segment.cross(norm) + dist = DraftVecUtils.project(chord,perp) + if not dist: return None + newpoint = point.add(dist) + if (dist.Length == 0): + return None + if strict: + s1 = newpoint.sub(edge.Vertexes[0].Point) + s2 = newpoint.sub(edge.Vertexes[-1].Point) + if (s1.Length <= segment.Length) and (s2.Length <= segment.Length): + return dist + else: + return None + else: return dist + elif isinstance(edge.Curve, Part.Circle): + ve1 = edge.Vertexes[0].Point + if (len(edge.Vertexes) > 1): + ve2 = edge.Vertexes[-1].Point + else: + ve2 = None + center = edge.Curve.Center + segment = center.sub(point) + ratio = (segment.Length - edge.Curve.Radius) / segment.Length + dist = DraftVecUtils.scale(segment,ratio) + newpoint = Vector.add(point, dist) + if (dist.Length == 0): + return None + if strict and ve2: + ang1 = DraftVecUtils.angle(ve1.sub(center)) + ang2 = DraftVecUtils.angle(ve2.sub(center)) + angpt = DraftVecUtils.angle(newpoint.sub(center)) + if ((angpt <= ang2 and angpt >= ang1) or (angpt <= ang1 and angpt >= ang2)): + return dist + else: + return None + else: + return dist + elif isinstance(edge.Curve,Part.BSplineCurve): + try: + pr = edge.Curve.parameter(point) + np = edge.Curve.value(pr) + dist = np.sub(point) + except: + print "DraftGeomUtils: Unable to get curve parameter for point ",point + return None + else: + return dist + else: + print "DraftGeomUtils: Couldn't project point" + return None + else: + print "DraftGeomUtils: Couldn't project point" + return None def angleBisection(edge1, edge2): - "angleBisection(edge,edge) - Returns an edge that bisects the angle between the 2 edges." - if isinstance(edge1.Curve, Part.Line) and isinstance(edge2.Curve, Part.Line): - p1 = edge1.Vertexes[0].Point - p2 = edge1.Vertexes[-1].Point - p3 = edge2.Vertexes[0].Point - p4 = edge2.Vertexes[-1].Point - int = findIntersection(edge1, edge2, True, True) - if int: - line1Dir = p2.sub(p1) - angleDiff = DraftVecUtils.angle(line1Dir, p4.sub(p3)) - ang = angleDiff * 0.5 - origin = int[0] - line1Dir.normalize() - dir = DraftVecUtils.rotate(line1Dir, ang) - return Part.Line(origin,origin.add(dir)).toShape() - else: - diff = p3.sub(p1) - origin = p1.add(DraftVecUtils.scale(diff, 0.5)) - dir = p2.sub(p1); dir.normalize() - return Part.Line(origin,origin.add(dir)).toShape() - else: - return None + "angleBisection(edge,edge) - Returns an edge that bisects the angle between the 2 edges." + if isinstance(edge1.Curve, Part.Line) and isinstance(edge2.Curve, Part.Line): + p1 = edge1.Vertexes[0].Point + p2 = edge1.Vertexes[-1].Point + p3 = edge2.Vertexes[0].Point + p4 = edge2.Vertexes[-1].Point + int = findIntersection(edge1, edge2, True, True) + if int: + line1Dir = p2.sub(p1) + angleDiff = DraftVecUtils.angle(line1Dir, p4.sub(p3)) + ang = angleDiff * 0.5 + origin = int[0] + line1Dir.normalize() + dir = DraftVecUtils.rotate(line1Dir, ang) + return Part.Line(origin,origin.add(dir)).toShape() + else: + diff = p3.sub(p1) + origin = p1.add(DraftVecUtils.scale(diff, 0.5)) + dir = p2.sub(p1); dir.normalize() + return Part.Line(origin,origin.add(dir)).toShape() + else: + return None def findClosestCircle(point,circles): - "findClosestCircle(Vector, list of circles) -- returns the circle with closest center" - dist = 1000000 - closest = None - for c in circles: - if c.Center.sub(point).Length < dist: - dist = c.Center.sub(point).Length - closest = c - return closest + "findClosestCircle(Vector, list of circles) -- returns the circle with closest center" + dist = 1000000 + closest = None + for c in circles: + if c.Center.sub(point).Length < dist: + dist = c.Center.sub(point).Length + closest = c + return closest def isCoplanar(faces): "checks if all faces in the given list are coplanar" @@ -1023,16 +1028,16 @@ def isCoplanar(faces): return True def isPlanar(shape): - "checks if the given shape is planar" - if len(shape.Vertexes) <= 3: - return True - n = getNormal(shape) - for p in shape.Vertexes[1:]: - pv = p.Point.sub(shape.Vertexes[0].Point) - rv = DraftVecUtils.project(pv,n) - if not DraftVecUtils.isNull(rv): - return False - return True + "checks if the given shape is planar" + if len(shape.Vertexes) <= 3: + return True + n = getNormal(shape) + for p in shape.Vertexes[1:]: + pv = p.Point.sub(shape.Vertexes[0].Point) + rv = DraftVecUtils.project(pv,n) + if not DraftVecUtils.isNull(rv): + return False + return True def findWiresOld(edges): '''finds connected edges in the list, and returns a list of lists containing edges @@ -1328,275 +1333,275 @@ def arcFromSpline(edge): # Fillet code graciously donated by Jacques-Antoine Gaudin def fillet(lEdges,r): - ''' Take a list of two Edges & a float as argument, - Returns a list of sorted edges describing a round corner''' + ''' Take a list of two Edges & a float as argument, + Returns a list of sorted edges describing a round corner''' - def getCurveType(edge,existingCurveType = None): - '''Builds or completes a dictionnary containing edges with keys "Arc" and "Line"''' - if not existingCurveType : - existingCurveType = { 'Line' : [], 'Arc' : [] } - if issubclass(type(edge.Curve),Part.Line) : - existingCurveType['Line'] += [edge] - elif issubclass(type(edge.Curve),Part.Circle) : - existingCurveType['Arc'] += [edge] + def getCurveType(edge,existingCurveType = None): + '''Builds or completes a dictionnary containing edges with keys "Arc" and "Line"''' + if not existingCurveType : + existingCurveType = { 'Line' : [], 'Arc' : [] } + if issubclass(type(edge.Curve),Part.Line) : + existingCurveType['Line'] += [edge] + elif issubclass(type(edge.Curve),Part.Circle) : + existingCurveType['Arc'] += [edge] + else : + raise Exception("Edge's curve must be either Line or Arc") + return existingCurveType + + rndEdges = lEdges[0:2] + rndEdges = sortEdges(rndEdges) + + if len(rndEdges) < 2 : + return rndEdges + + if r <= 0 : + print "DraftGeomUtils.fillet : Error : radius is negative." + return rndEdges + + curveType = getCurveType(rndEdges[0]) + curveType = getCurveType(rndEdges[1],curveType) + + lVertexes = rndEdges[0].Vertexes + [rndEdges[1].Vertexes[-1]] + + if len(curveType['Line']) == 2: + + # Deals with 2-line-edges lists -------------------------------------- + + U1 = lVertexes[0].Point.sub(lVertexes[1].Point) ; U1.normalize() + U2 = lVertexes[2].Point.sub(lVertexes[1].Point) ; U2.normalize() + alpha = U1.getAngle(U2) + + if round(alpha,precision()) == 0 or round(alpha - math.pi,precision()) == 0: # Edges have same direction + print "DraftGeomUtils.fillet : Warning : edges have same direction. Did nothing" + return rndEdges + + dToCenter = r / math.sin(alpha/2.) + dToTangent = (dToCenter**2-r**2)**(0.5) + dirVect = Vector(U1) ; dirVect.scale(dToTangent,dToTangent,dToTangent) + arcPt1 = lVertexes[1].Point.add(dirVect) + + dirVect = U2.add(U1) ; dirVect.normalize() + dirVect.scale(dToCenter-r,dToCenter-r,dToCenter-r) + arcPt2 = lVertexes[1].Point.add(dirVect) + + dirVect = Vector(U2) ; dirVect.scale(dToTangent,dToTangent,dToTangent) + arcPt3 = lVertexes[1].Point.add(dirVect) + + if (dToTangent>lEdges[0].Length) or (dToTangent>lEdges[1].Length) : + print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" + return rndEdges + + rndEdges[1] = Part.Edge(Part.Arc(arcPt1,arcPt2,arcPt3)) + rndEdges[0] = Part.Edge(Part.Line(lVertexes[0].Point,arcPt1)) + rndEdges += [Part.Edge(Part.Line(arcPt3,lVertexes[2].Point))] + + return rndEdges + + elif len(curveType['Arc']) == 1 : + + # Deals with lists containing an arc and a line ---------------------------------- + + if lEdges[0] in curveType['Arc'] : + lineEnd = lVertexes[2] ; arcEnd = lVertexes[0] ; arcFirst = True + else : + lineEnd = lVertexes[0] ; arcEnd = lVertexes[2] ; arcFirst = False + arcCenter = curveType['Arc'][0].Curve.Center + arcRadius = curveType['Arc'][0].Curve.Radius + arcAxis = curveType['Arc'][0].Curve.Axis + arcLength = curveType['Arc'][0].Length + + U1 = lineEnd.Point.sub(lVertexes[1].Point) ; U1.normalize() + toCenter = arcCenter.sub(lVertexes[1].Point) + if arcFirst : # make sure the tangent points towards the arc + T = arcAxis.cross(toCenter) + else : + T = toCenter.cross(arcAxis) + + projCenter = toCenter.dot(U1) + if round(abs(projCenter),precision()) > 0 : + normToLine = U1.cross(T).cross(U1) + else : + normToLine = Vector(toCenter) + normToLine.normalize() + + dCenterToLine = toCenter.dot(normToLine) - r + + if round(projCenter,precision()) > 0 : + newRadius = arcRadius - r + elif round(projCenter,precision()) < 0 or (round(projCenter,precision()) == 0 and U1.dot(T) > 0): + newRadius = arcRadius + r + else : + print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" + return rndEdges + + toNewCent = newRadius**2-dCenterToLine**2 + if toNewCent > 0 : + toNewCent = abs(abs(projCenter) - toNewCent**(0.5)) + else : + print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" + return rndEdges + + U1.scale(toNewCent,toNewCent,toNewCent) + normToLine.scale(r,r,r) + newCent = lVertexes[1].Point.add(U1).add(normToLine) + + arcPt1= lVertexes[1].Point.add(U1) + arcPt2= lVertexes[1].Point.sub(newCent); arcPt2.normalize() + arcPt2.scale(r,r,r) ; arcPt2 = arcPt2.add(newCent) + if newRadius == arcRadius - r : + arcPt3= newCent.sub(arcCenter) + else : + arcPt3= arcCenter.sub(newCent) + arcPt3.normalize() + arcPt3.scale(r,r,r) ; arcPt3 = arcPt3.add(newCent) + arcPt = [arcPt1,arcPt2,arcPt3] + + + # Warning : In the following I used a trick for calling the right element + # in arcPt or V : arcFirst is a boolean so - not arcFirst is -0 or -1 + # list[-1] is the last element of a list and list[0] the first + # this way I don't have to proceed tests to know the position of the arc + + myTrick = not arcFirst + + V = [arcPt3] + V += [arcEnd.Point] + + toCenter.scale(-1,-1,-1) + + delLength = arcRadius * V[0].sub(arcCenter).getAngle(toCenter) + if delLength > arcLength or toNewCent > curveType['Line'][0].Length: + print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" + return rndEdges + + arcAsEdge = arcFrom2Pts(V[-arcFirst],V[-myTrick],arcCenter,arcAxis) + + V = [lineEnd.Point,arcPt1] + lineAsEdge = Part.Edge(Part.Line(V[-arcFirst],V[myTrick])) + + rndEdges[not arcFirst] = arcAsEdge + rndEdges[arcFirst] = lineAsEdge + rndEdges[1:1] = [Part.Edge(Part.Arc(arcPt[- arcFirst],arcPt[1],arcPt[- myTrick]))] + + return rndEdges + + elif len(curveType['Arc']) == 2 : + + # Deals with lists of 2 arc-edges -------------------------------------------- + + arcCenter, arcRadius, arcAxis, arcLength, toCenter, T, newRadius = [], [], [], [], [], [], [] + for i in range(2) : + arcCenter += [curveType['Arc'][i].Curve.Center] + arcRadius += [curveType['Arc'][i].Curve.Radius] + arcAxis += [curveType['Arc'][i].Curve.Axis] + arcLength += [curveType['Arc'][i].Length] + toCenter += [arcCenter[i].sub(lVertexes[1].Point)] + T += [arcAxis[0].cross(toCenter[0])] + T += [toCenter[1].cross(arcAxis[1])] + CentToCent = toCenter[1].sub(toCenter[0]) + dCentToCent = CentToCent.Length + + sameDirection = (arcAxis[0].dot(arcAxis[1]) > 0) + TcrossT = T[0].cross(T[1]) + if sameDirection : + if round(TcrossT.dot(arcAxis[0]),precision()) > 0 : + newRadius += [arcRadius[0]+r] + newRadius += [arcRadius[1]+r] + elif round(TcrossT.dot(arcAxis[0]),precision()) < 0 : + newRadius += [arcRadius[0]-r] + newRadius += [arcRadius[1]-r] + elif T[0].dot(T[1]) > 0 : + newRadius += [arcRadius[0]+r] + newRadius += [arcRadius[1]+r] + else : + print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" + return rndEdges + elif not sameDirection : + if round(TcrossT.dot(arcAxis[0]),precision()) > 0 : + newRadius += [arcRadius[0]+r] + newRadius += [arcRadius[1]-r] + elif round(TcrossT.dot(arcAxis[0]),precision()) < 0 : + newRadius += [arcRadius[0]-r] + newRadius += [arcRadius[1]+r] + elif T[0].dot(T[1]) > 0 : + if arcRadius[0] > arcRadius[1] : + newRadius += [arcRadius[0]-r] + newRadius += [arcRadius[1]+r] + elif arcRadius[1] > arcRadius[0] : + newRadius += [arcRadius[0]+r] + newRadius += [arcRadius[1]-r] else : - raise Exception("Edge's curve must be either Line or Arc") - return existingCurveType - - rndEdges = lEdges[0:2] - rndEdges = sortEdges(rndEdges) - - if len(rndEdges) < 2 : - return rndEdges - - if r <= 0 : - print "DraftGeomUtils.fillet : Error : radius is negative." - return rndEdges - - curveType = getCurveType(rndEdges[0]) - curveType = getCurveType(rndEdges[1],curveType) - - lVertexes = rndEdges[0].Vertexes + [rndEdges[1].Vertexes[-1]] - - if len(curveType['Line']) == 2: - - # Deals with 2-line-edges lists -------------------------------------- - - U1 = lVertexes[0].Point.sub(lVertexes[1].Point) ; U1.normalize() - U2 = lVertexes[2].Point.sub(lVertexes[1].Point) ; U2.normalize() - alpha = U1.getAngle(U2) - - if round(alpha,precision()) == 0 or round(alpha - math.pi,precision()) == 0: # Edges have same direction - print "DraftGeomUtils.fillet : Warning : edges have same direction. Did nothing" - return rndEdges - - dToCenter = r / math.sin(alpha/2.) - dToTangent = (dToCenter**2-r**2)**(0.5) - dirVect = Vector(U1) ; dirVect.scale(dToTangent,dToTangent,dToTangent) - arcPt1 = lVertexes[1].Point.add(dirVect) - - dirVect = U2.add(U1) ; dirVect.normalize() - dirVect.scale(dToCenter-r,dToCenter-r,dToCenter-r) - arcPt2 = lVertexes[1].Point.add(dirVect) - - dirVect = Vector(U2) ; dirVect.scale(dToTangent,dToTangent,dToTangent) - arcPt3 = lVertexes[1].Point.add(dirVect) - - if (dToTangent>lEdges[0].Length) or (dToTangent>lEdges[1].Length) : - print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" - return rndEdges - - rndEdges[1] = Part.Edge(Part.Arc(arcPt1,arcPt2,arcPt3)) - rndEdges[0] = Part.Edge(Part.Line(lVertexes[0].Point,arcPt1)) - rndEdges += [Part.Edge(Part.Line(arcPt3,lVertexes[2].Point))] - - return rndEdges - - elif len(curveType['Arc']) == 1 : - - # Deals with lists containing an arc and a line ---------------------------------- - - if lEdges[0] in curveType['Arc'] : - lineEnd = lVertexes[2] ; arcEnd = lVertexes[0] ; arcFirst = True - else : - lineEnd = lVertexes[0] ; arcEnd = lVertexes[2] ; arcFirst = False - arcCenter = curveType['Arc'][0].Curve.Center - arcRadius = curveType['Arc'][0].Curve.Radius - arcAxis = curveType['Arc'][0].Curve.Axis - arcLength = curveType['Arc'][0].Length - - U1 = lineEnd.Point.sub(lVertexes[1].Point) ; U1.normalize() - toCenter = arcCenter.sub(lVertexes[1].Point) - if arcFirst : # make sure the tangent points towards the arc - T = arcAxis.cross(toCenter) - else : - T = toCenter.cross(arcAxis) - - projCenter = toCenter.dot(U1) - if round(abs(projCenter),precision()) > 0 : - normToLine = U1.cross(T).cross(U1) - else : - normToLine = Vector(toCenter) - normToLine.normalize() - - dCenterToLine = toCenter.dot(normToLine) - r - - if round(projCenter,precision()) > 0 : - newRadius = arcRadius - r - elif round(projCenter,precision()) < 0 or (round(projCenter,precision()) == 0 and U1.dot(T) > 0): - newRadius = arcRadius + r - else : - print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" - return rndEdges - - toNewCent = newRadius**2-dCenterToLine**2 - if toNewCent > 0 : - toNewCent = abs(abs(projCenter) - toNewCent**(0.5)) - else : - print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" - return rndEdges - - U1.scale(toNewCent,toNewCent,toNewCent) - normToLine.scale(r,r,r) - newCent = lVertexes[1].Point.add(U1).add(normToLine) - - arcPt1= lVertexes[1].Point.add(U1) - arcPt2= lVertexes[1].Point.sub(newCent); arcPt2.normalize() - arcPt2.scale(r,r,r) ; arcPt2 = arcPt2.add(newCent) - if newRadius == arcRadius - r : - arcPt3= newCent.sub(arcCenter) - else : - arcPt3= arcCenter.sub(newCent) - arcPt3.normalize() - arcPt3.scale(r,r,r) ; arcPt3 = arcPt3.add(newCent) - arcPt = [arcPt1,arcPt2,arcPt3] - - - # Warning : In the following I used a trick for calling the right element - # in arcPt or V : arcFirst is a boolean so - not arcFirst is -0 or -1 - # list[-1] is the last element of a list and list[0] the first - # this way I don't have to proceed tests to know the position of the arc - - myTrick = not arcFirst - - V = [arcPt3] - V += [arcEnd.Point] - - toCenter.scale(-1,-1,-1) - - delLength = arcRadius * V[0].sub(arcCenter).getAngle(toCenter) - if delLength > arcLength or toNewCent > curveType['Line'][0].Length: - print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" - return rndEdges - - arcAsEdge = arcFrom2Pts(V[-arcFirst],V[-myTrick],arcCenter,arcAxis) - - V = [lineEnd.Point,arcPt1] - lineAsEdge = Part.Edge(Part.Line(V[-arcFirst],V[myTrick])) - - rndEdges[not arcFirst] = arcAsEdge - rndEdges[arcFirst] = lineAsEdge - rndEdges[1:1] = [Part.Edge(Part.Arc(arcPt[- arcFirst],arcPt[1],arcPt[- myTrick]))] - - return rndEdges - - elif len(curveType['Arc']) == 2 : - - # Deals with lists of 2 arc-edges -------------------------------------------- - - arcCenter, arcRadius, arcAxis, arcLength, toCenter, T, newRadius = [], [], [], [], [], [], [] - for i in range(2) : - arcCenter += [curveType['Arc'][i].Curve.Center] - arcRadius += [curveType['Arc'][i].Curve.Radius] - arcAxis += [curveType['Arc'][i].Curve.Axis] - arcLength += [curveType['Arc'][i].Length] - toCenter += [arcCenter[i].sub(lVertexes[1].Point)] - T += [arcAxis[0].cross(toCenter[0])] - T += [toCenter[1].cross(arcAxis[1])] - CentToCent = toCenter[1].sub(toCenter[0]) - dCentToCent = CentToCent.Length - - sameDirection = (arcAxis[0].dot(arcAxis[1]) > 0) - TcrossT = T[0].cross(T[1]) - if sameDirection : - if round(TcrossT.dot(arcAxis[0]),precision()) > 0 : - newRadius += [arcRadius[0]+r] - newRadius += [arcRadius[1]+r] - elif round(TcrossT.dot(arcAxis[0]),precision()) < 0 : - newRadius += [arcRadius[0]-r] - newRadius += [arcRadius[1]-r] - elif T[0].dot(T[1]) > 0 : - newRadius += [arcRadius[0]+r] - newRadius += [arcRadius[1]+r] - else : - print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" - return rndEdges - elif not sameDirection : - if round(TcrossT.dot(arcAxis[0]),precision()) > 0 : - newRadius += [arcRadius[0]+r] - newRadius += [arcRadius[1]-r] - elif round(TcrossT.dot(arcAxis[0]),precision()) < 0 : - newRadius += [arcRadius[0]-r] - newRadius += [arcRadius[1]+r] - elif T[0].dot(T[1]) > 0 : - if arcRadius[0] > arcRadius[1] : - newRadius += [arcRadius[0]-r] - newRadius += [arcRadius[1]+r] - elif arcRadius[1] > arcRadius[0] : - newRadius += [arcRadius[0]+r] - newRadius += [arcRadius[1]-r] - else : - print "DraftGeomUtils.fillet : Warning : arcs are coincident. Did nothing" - return rndEdges - else : - print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" - return rndEdges - - if newRadius[0]+newRadius[1] < dCentToCent or \ - newRadius[0]-newRadius[1] > dCentToCent or \ - newRadius[1]-newRadius[0] > dCentToCent : - print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" - return rndEdges - - x = (dCentToCent**2+newRadius[0]**2-newRadius[1]**2)/(2*dCentToCent) - y = (newRadius[0]**2-x**2)**(0.5) - - CentToCent.normalize() ; toCenter[0].normalize() ; toCenter[1].normalize() - if abs(toCenter[0].dot(toCenter[1])) != 1 : - normVect = CentToCent.cross(CentToCent.cross(toCenter[0])) - else : - normVect = T[0] - normVect.normalize() - CentToCent.scale(x,x,x) ; normVect.scale(y,y,y) - newCent = arcCenter[0].add(CentToCent.add(normVect)) - CentToNewCent = [newCent.sub(arcCenter[0]),newCent.sub(arcCenter[1])] - for i in range(2) : - CentToNewCent[i].normalize() - if newRadius[i] == arcRadius[i]+r : - CentToNewCent[i].scale(-r,-r,-r) - else : - CentToNewCent[i].scale(r,r,r) - toThirdPt = lVertexes[1].Point.sub(newCent) ; toThirdPt.normalize() - toThirdPt.scale(r,r,r) - arcPt1 = newCent.add(CentToNewCent[0]) - arcPt2 = newCent.add(toThirdPt) - arcPt3 = newCent.add(CentToNewCent[1]) - arcPt = [arcPt1,arcPt2,arcPt3] - - arcAsEdge = [] - for i in range(2) : - toCenter[i].scale(-1,-1,-1) - delLength = arcRadius[i] * arcPt[-i].sub(arcCenter[i]).getAngle(toCenter[i]) - if delLength > arcLength[i] : - print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" - return rndEdges - V = [arcPt[-i],lVertexes[-i].Point] - arcAsEdge += [arcFrom2Pts(V[i-1],V[-i],arcCenter[i],arcAxis[i])] - - rndEdges[0] = arcAsEdge[0] - rndEdges[1] = arcAsEdge[1] - rndEdges[1:1] = [Part.Edge(Part.Arc(arcPt[0],arcPt[1],arcPt[2]))] - - return rndEdges - + print "DraftGeomUtils.fillet : Warning : arcs are coincident. Did nothing" + return rndEdges + else : + print "DraftGeomUtils.fillet : Warning : edges are already tangent. Did nothing" + return rndEdges + + if newRadius[0]+newRadius[1] < dCentToCent or \ + newRadius[0]-newRadius[1] > dCentToCent or \ + newRadius[1]-newRadius[0] > dCentToCent : + print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" + return rndEdges + + x = (dCentToCent**2+newRadius[0]**2-newRadius[1]**2)/(2*dCentToCent) + y = (newRadius[0]**2-x**2)**(0.5) + + CentToCent.normalize() ; toCenter[0].normalize() ; toCenter[1].normalize() + if abs(toCenter[0].dot(toCenter[1])) != 1 : + normVect = CentToCent.cross(CentToCent.cross(toCenter[0])) + else : + normVect = T[0] + normVect.normalize() + CentToCent.scale(x,x,x) ; normVect.scale(y,y,y) + newCent = arcCenter[0].add(CentToCent.add(normVect)) + CentToNewCent = [newCent.sub(arcCenter[0]),newCent.sub(arcCenter[1])] + for i in range(2) : + CentToNewCent[i].normalize() + if newRadius[i] == arcRadius[i]+r : + CentToNewCent[i].scale(-r,-r,-r) + else : + CentToNewCent[i].scale(r,r,r) + toThirdPt = lVertexes[1].Point.sub(newCent) ; toThirdPt.normalize() + toThirdPt.scale(r,r,r) + arcPt1 = newCent.add(CentToNewCent[0]) + arcPt2 = newCent.add(toThirdPt) + arcPt3 = newCent.add(CentToNewCent[1]) + arcPt = [arcPt1,arcPt2,arcPt3] + + arcAsEdge = [] + for i in range(2) : + toCenter[i].scale(-1,-1,-1) + delLength = arcRadius[i] * arcPt[-i].sub(arcCenter[i]).getAngle(toCenter[i]) + if delLength > arcLength[i] : + print "DraftGeomUtils.fillet : Error : radius value ", r," is too high" + return rndEdges + V = [arcPt[-i],lVertexes[-i].Point] + arcAsEdge += [arcFrom2Pts(V[i-1],V[-i],arcCenter[i],arcAxis[i])] + + rndEdges[0] = arcAsEdge[0] + rndEdges[1] = arcAsEdge[1] + rndEdges[1:1] = [Part.Edge(Part.Arc(arcPt[0],arcPt[1],arcPt[2]))] + + return rndEdges + def filletWire(aWire,r,makeClosed=True): - ''' Fillets each angle of a wire with r as radius value''' - - edges = aWire.Edges - edges = sortEdges(edges) - filEdges = [edges[0]] - for i in range(len(edges)-1): - result = fillet([filEdges[-1],edges[i+1]],r) - if len(result)>2: - filEdges[-1:] = result[0:3] - else : - filEdges[-1:] = result[0:2] - if isReallyClosed(aWire) and makeClosed : - result = fillet([filEdges[-1],filEdges[0]],r) - if len(result)>2: - filEdges[-1:] = result[0:2] - filEdges[0] = result[2] - return Part.Wire(filEdges) + ''' Fillets each angle of a wire with r as radius value''' + + edges = aWire.Edges + edges = sortEdges(edges) + filEdges = [edges[0]] + for i in range(len(edges)-1): + result = fillet([filEdges[-1],edges[i+1]],r) + if len(result)>2: + filEdges[-1:] = result[0:3] + else : + filEdges[-1:] = result[0:2] + if isReallyClosed(aWire) and makeClosed : + result = fillet([filEdges[-1],filEdges[0]],r) + if len(result)>2: + filEdges[-1:] = result[0:2] + filEdges[0] = result[2] + return Part.Wire(filEdges) # circle functions ********************************************************* @@ -1644,236 +1649,236 @@ def getBoundaryAngles(angle,alist): def circleFrom2tan1pt(tan1, tan2, point): - "circleFrom2tan1pt(edge, edge, Vector)" - if isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): - return circlefrom2Lines1Point(tan1, tan2, point) - elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): - return circlefromCircleLinePoint(tan1, tan2, point) - elif isinstance(tan2.Curve, Part.Circle) and isinstance(tan1.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): - return circlefromCircleLinePoint(tan2, tan1, point) - elif isinstance(tan2.Curve, Part.Circle) and isinstance(tan1.Curve, Part.Circle) and isinstance(point, FreeCAD.Vector): - return circlefrom2Circles1Point(tan2, tan1, point) + "circleFrom2tan1pt(edge, edge, Vector)" + if isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): + return circlefrom2Lines1Point(tan1, tan2, point) + elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): + return circlefromCircleLinePoint(tan1, tan2, point) + elif isinstance(tan2.Curve, Part.Circle) and isinstance(tan1.Curve, Part.Line) and isinstance(point, FreeCAD.Vector): + return circlefromCircleLinePoint(tan2, tan1, point) + elif isinstance(tan2.Curve, Part.Circle) and isinstance(tan1.Curve, Part.Circle) and isinstance(point, FreeCAD.Vector): + return circlefrom2Circles1Point(tan2, tan1, point) def circleFrom2tan1rad(tan1, tan2, rad): - "circleFrom2tan1rad(edge, edge, float)" - if isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Line): - return circleFrom2LinesRadius(tan1, tan2, rad) - elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Line): - return circleFromCircleLineRadius(tan1, tan2, rad) - elif isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Circle): - return circleFromCircleLineRadius(tan2, tan1, rad) - elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Circle): - return circleFrom2CirclesRadius(tan1, tan2, rad) + "circleFrom2tan1rad(edge, edge, float)" + if isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Line): + return circleFrom2LinesRadius(tan1, tan2, rad) + elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Line): + return circleFromCircleLineRadius(tan1, tan2, rad) + elif isinstance(tan1.Curve, Part.Line) and isinstance(tan2.Curve, Part.Circle): + return circleFromCircleLineRadius(tan2, tan1, rad) + elif isinstance(tan1.Curve, Part.Circle) and isinstance(tan2.Curve, Part.Circle): + return circleFrom2CirclesRadius(tan1, tan2, rad) def circleFrom1tan2pt(tan1, p1, p2): - if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector) and isinstance(p2, FreeCAD.Vector): - return circlefrom1Line2Points(tan1, p1, p2) - if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector) and isinstance(p2, FreeCAD.Vector): - return circlefrom1Circle2Points(tan1, p1, p2) + if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector) and isinstance(p2, FreeCAD.Vector): + return circlefrom1Line2Points(tan1, p1, p2) + if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector) and isinstance(p2, FreeCAD.Vector): + return circlefrom1Circle2Points(tan1, p1, p2) def circleFrom1tan1pt1rad(tan1, p1, rad): - if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector): - return circleFromPointLineRadius(p1, tan1, rad) - if isinstance(tan1.Curve, Part.Circle) and isinstance(p1, FreeCAD.Vector): - return circleFromPointCircleRadius(p1, tan1, rad) + if isinstance(tan1.Curve, Part.Line) and isinstance(p1, FreeCAD.Vector): + return circleFromPointLineRadius(p1, tan1, rad) + if isinstance(tan1.Curve, Part.Circle) and isinstance(p1, FreeCAD.Vector): + return circleFromPointCircleRadius(p1, tan1, rad) def circleFrom3tan(tan1, tan2, tan3): - tan1IsLine = isinstance(tan1.Curve, Part.Line) - tan2IsLine = isinstance(tan2.Curve, Part.Line) - tan3IsLine = isinstance(tan3.Curve, Part.Line) - tan1IsCircle = isinstance(tan1.Curve, Part.Circle) - tan2IsCircle = isinstance(tan2.Curve, Part.Circle) - tan3IsCircle = isinstance(tan3.Curve, Part.Circle) - if tan1IsLine and tan2IsLine and tan3IsLine: - return circleFrom3LineTangents(tan1, tan2, tan3) - elif tan1IsCircle and tan2IsCircle and tan3IsCircle: - return circleFrom3CircleTangents(tan1, tan2, tan3) - elif (tan1IsCircle and tan2IsLine and tan3IsLine): - return circleFrom1Circle2Lines(tan1, tan2, tan3) - elif (tan1IsLine and tan2IsCircle and tan3IsLine): - return circleFrom1Circle2Lines(tan2, tan1, tan3) - elif (tan1IsLine and tan2IsLine and tan3IsCircle): - return circleFrom1Circle2Lines(tan3, tan1, tan2) - elif (tan1IsLine and tan2IsCircle and tan3IsCircle): - return circleFrom2Circle1Lines(tan2, tan3, tan1) - elif (tan1IsCircle and tan2IsLine and tan3IsCircle): - return circleFrom2Circle1Lines(tan1, tan3, tan2) - elif (tan1IsCircle and tan2IsCircle and tan3IsLine): - return circleFrom2Circle1Lines(tan1, tan2, tan3) + tan1IsLine = isinstance(tan1.Curve, Part.Line) + tan2IsLine = isinstance(tan2.Curve, Part.Line) + tan3IsLine = isinstance(tan3.Curve, Part.Line) + tan1IsCircle = isinstance(tan1.Curve, Part.Circle) + tan2IsCircle = isinstance(tan2.Curve, Part.Circle) + tan3IsCircle = isinstance(tan3.Curve, Part.Circle) + if tan1IsLine and tan2IsLine and tan3IsLine: + return circleFrom3LineTangents(tan1, tan2, tan3) + elif tan1IsCircle and tan2IsCircle and tan3IsCircle: + return circleFrom3CircleTangents(tan1, tan2, tan3) + elif (tan1IsCircle and tan2IsLine and tan3IsLine): + return circleFrom1Circle2Lines(tan1, tan2, tan3) + elif (tan1IsLine and tan2IsCircle and tan3IsLine): + return circleFrom1Circle2Lines(tan2, tan1, tan3) + elif (tan1IsLine and tan2IsLine and tan3IsCircle): + return circleFrom1Circle2Lines(tan3, tan1, tan2) + elif (tan1IsLine and tan2IsCircle and tan3IsCircle): + return circleFrom2Circle1Lines(tan2, tan3, tan1) + elif (tan1IsCircle and tan2IsLine and tan3IsCircle): + return circleFrom2Circle1Lines(tan1, tan3, tan2) + elif (tan1IsCircle and tan2IsCircle and tan3IsLine): + return circleFrom2Circle1Lines(tan1, tan2, tan3) def circlefrom2Lines1Point(edge1, edge2, point): - "circlefrom2Lines1Point(edge, edge, Vector)" - bis = angleBisection(edge1, edge2) - if not bis: return None - mirrPoint = mirror(point, bis) - return circlefrom1Line2Points(edge1, point, mirrPoint) + "circlefrom2Lines1Point(edge, edge, Vector)" + bis = angleBisection(edge1, edge2) + if not bis: return None + mirrPoint = mirror(point, bis) + return circlefrom1Line2Points(edge1, point, mirrPoint) def circlefrom1Line2Points(edge, p1, p2): - "circlefrom1Line2Points(edge, Vector, Vector)" - p1_p2 = edg(p1, p2) - s = findIntersection(edge, p1_p2, True, True) - if not s: return None - s = s[0] - v1 = p1.sub(s) - v2 = p2.sub(s) - projectedDist = math.sqrt(abs(v1.dot(v2))) - edgeDir = vec(edge); edgeDir.normalize() - projectedCen1 = Vector.add(s, DraftVecUtils.scale(edgeDir, projectedDist)) - projectedCen2 = Vector.add(s, DraftVecUtils.scale(edgeDir, -projectedDist)) - perpEdgeDir = edgeDir.cross(Vector(0,0,1)) - perpCen1 = Vector.add(projectedCen1, perpEdgeDir) - perpCen2 = Vector.add(projectedCen2, perpEdgeDir) - mid = findMidpoint(p1_p2) - x = DraftVecUtils.crossproduct(vec(p1_p2)); x.normalize() - perp_mid = Vector.add(mid, x) - cen1 = findIntersection(edg(projectedCen1, perpCen1), edg(mid, perp_mid), True, True) - cen2 = findIntersection(edg(projectedCen2, perpCen2), edg(mid, perp_mid), True, True) - circles = [] - if cen1: - radius = DraftVecUtils.dist(projectedCen1, cen1[0]) - circles.append(Part.Circle(cen1[0], NORM, radius)) - if cen2: - radius = DraftVecUtils.dist(projectedCen2, cen2[0]) - circles.append(Part.Circle(cen2[0], NORM, radius)) + "circlefrom1Line2Points(edge, Vector, Vector)" + p1_p2 = edg(p1, p2) + s = findIntersection(edge, p1_p2, True, True) + if not s: return None + s = s[0] + v1 = p1.sub(s) + v2 = p2.sub(s) + projectedDist = math.sqrt(abs(v1.dot(v2))) + edgeDir = vec(edge); edgeDir.normalize() + projectedCen1 = Vector.add(s, DraftVecUtils.scale(edgeDir, projectedDist)) + projectedCen2 = Vector.add(s, DraftVecUtils.scale(edgeDir, -projectedDist)) + perpEdgeDir = edgeDir.cross(Vector(0,0,1)) + perpCen1 = Vector.add(projectedCen1, perpEdgeDir) + perpCen2 = Vector.add(projectedCen2, perpEdgeDir) + mid = findMidpoint(p1_p2) + x = DraftVecUtils.crossproduct(vec(p1_p2)); x.normalize() + perp_mid = Vector.add(mid, x) + cen1 = findIntersection(edg(projectedCen1, perpCen1), edg(mid, perp_mid), True, True) + cen2 = findIntersection(edg(projectedCen2, perpCen2), edg(mid, perp_mid), True, True) + circles = [] + if cen1: + radius = DraftVecUtils.dist(projectedCen1, cen1[0]) + circles.append(Part.Circle(cen1[0], NORM, radius)) + if cen2: + radius = DraftVecUtils.dist(projectedCen2, cen2[0]) + circles.append(Part.Circle(cen2[0], NORM, radius)) - if circles: return circles - else: return None + if circles: return circles + else: return None def circleFrom2LinesRadius (edge1, edge2, radius): - "circleFrom2LinesRadius(edge,edge,radius)" - int = findIntersection(edge1, edge2, True, True) - if not int: return None - int = int[0] - bis12 = angleBisection(edge1,edge2) - bis21 = Part.Line(bis12.Vertexes[0].Point,DraftVecUtils.rotate(vec(bis12), math.pi/2.0)) - ang12 = abs(DraftVecUtils.angle(vec(edge1),vec(edge2))) - ang21 = math.pi - ang12 - dist12 = radius / math.sin(ang12 * 0.5) - dist21 = radius / math.sin(ang21 * 0.5) - circles = [] - cen = Vector.add(int, DraftVecUtils.scale(vec(bis12), dist12)) - circles.append(Part.Circle(cen, NORM, radius)) - cen = Vector.add(int, DraftVecUtils.scale(vec(bis12), -dist12)) - circles.append(Part.Circle(cen, NORM, radius)) - cen = Vector.add(int, DraftVecUtils.scale(vec(bis21), dist21)) - circles.append(Part.Circle(cen, NORM, radius)) - cen = Vector.add(int, DraftVecUtils.scale(vec(bis21), -dist21)) - circles.append(Part.Circle(cen, NORM, radius)) - return circles + "circleFrom2LinesRadius(edge,edge,radius)" + int = findIntersection(edge1, edge2, True, True) + if not int: return None + int = int[0] + bis12 = angleBisection(edge1,edge2) + bis21 = Part.Line(bis12.Vertexes[0].Point,DraftVecUtils.rotate(vec(bis12), math.pi/2.0)) + ang12 = abs(DraftVecUtils.angle(vec(edge1),vec(edge2))) + ang21 = math.pi - ang12 + dist12 = radius / math.sin(ang12 * 0.5) + dist21 = radius / math.sin(ang21 * 0.5) + circles = [] + cen = Vector.add(int, DraftVecUtils.scale(vec(bis12), dist12)) + circles.append(Part.Circle(cen, NORM, radius)) + cen = Vector.add(int, DraftVecUtils.scale(vec(bis12), -dist12)) + circles.append(Part.Circle(cen, NORM, radius)) + cen = Vector.add(int, DraftVecUtils.scale(vec(bis21), dist21)) + circles.append(Part.Circle(cen, NORM, radius)) + cen = Vector.add(int, DraftVecUtils.scale(vec(bis21), -dist21)) + circles.append(Part.Circle(cen, NORM, radius)) + return circles def circleFrom3LineTangents (edge1, edge2, edge3): - "circleFrom3LineTangents(edge,edge,edge)" - def rot(ed): - return Part.Line(v1(ed),v1(ed).add(DraftVecUtils.rotate(vec(ed),math.pi/2))).toShape() - bis12 = angleBisection(edge1,edge2) - bis23 = angleBisection(edge2,edge3) - bis31 = angleBisection(edge3,edge1) - intersections = [] - int = findIntersection(bis12, bis23, True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - int = findIntersection(bis23, bis31, True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - int = findIntersection(bis31, bis12, True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - int = findIntersection(rot(bis12), rot(bis23), True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - int = findIntersection(rot(bis23), rot(bis31), True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - int = findIntersection(rot(bis31), rot(bis12), True, True) - if int: - radius = findDistance(int[0],edge1).Length - intersections.append(Part.Circle(int[0],NORM,radius)) - circles = [] - for int in intersections: - exists = False - for cir in circles: - if DraftVecUtils.equals(cir.Center, int.Center): - exists = True - break - if not exists: - circles.append(int) - if circles: - return circles - else: - return None + "circleFrom3LineTangents(edge,edge,edge)" + def rot(ed): + return Part.Line(v1(ed),v1(ed).add(DraftVecUtils.rotate(vec(ed),math.pi/2))).toShape() + bis12 = angleBisection(edge1,edge2) + bis23 = angleBisection(edge2,edge3) + bis31 = angleBisection(edge3,edge1) + intersections = [] + int = findIntersection(bis12, bis23, True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + int = findIntersection(bis23, bis31, True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + int = findIntersection(bis31, bis12, True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + int = findIntersection(rot(bis12), rot(bis23), True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + int = findIntersection(rot(bis23), rot(bis31), True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + int = findIntersection(rot(bis31), rot(bis12), True, True) + if int: + radius = findDistance(int[0],edge1).Length + intersections.append(Part.Circle(int[0],NORM,radius)) + circles = [] + for int in intersections: + exists = False + for cir in circles: + if DraftVecUtils.equals(cir.Center, int.Center): + exists = True + break + if not exists: + circles.append(int) + if circles: + return circles + else: + return None def circleFromPointLineRadius (point, edge, radius): - "circleFromPointLineRadius (point, edge, radius)" - dist = findDistance(point, edge, False) - center1 = None - center2 = None - if dist.Length == 0: - segment = vec(edge) - perpVec = DraftVecUtils.crossproduct(segment); perpVec.normalize() - normPoint_c1 = DraftVecUtils.scale(perpVec, radius) - normPoint_c2 = DraftVecUtils.scale(perpVec, -radius) - center1 = point.add(normPoint_c1) - center2 = point.add(normPoint_c2) - elif dist.Length > 2 * radius: - return None - elif dist.Length == 2 * radius: - normPoint = point.add(findDistance(point, edge, False)) - dummy = DraftVecUtils.scale(normPoint.sub(point), 0.5) - cen = point.add(dummy) - circ = Part.Circle(cen, NORM, radius) - if circ: - return [circ] - else: - return None - else: - normPoint = point.add(findDistance(point, edge, False)) - normDist = DraftVecUtils.dist(normPoint, point) - dist = math.sqrt(radius**2 - (radius - normDist)**2) - centerNormVec = DraftVecUtils.scaleTo(point.sub(normPoint), radius) - edgeDir = edge.Vertexes[0].Point.sub(normPoint); edgeDir.normalize() - center1 = centerNormVec.add(normPoint.add(DraftVecUtils.scale(edgeDir, dist))) - center2 = centerNormVec.add(normPoint.add(DraftVecUtils.scale(edgeDir, -dist))) - circles = [] - if center1: - circ = Part.Circle(center1, NORM, radius) - if circ: - circles.append(circ) - if center2: - circ = Part.Circle(center2, NORM, radius) - if circ: - circles.append(circ) + "circleFromPointLineRadius (point, edge, radius)" + dist = findDistance(point, edge, False) + center1 = None + center2 = None + if dist.Length == 0: + segment = vec(edge) + perpVec = DraftVecUtils.crossproduct(segment); perpVec.normalize() + normPoint_c1 = DraftVecUtils.scale(perpVec, radius) + normPoint_c2 = DraftVecUtils.scale(perpVec, -radius) + center1 = point.add(normPoint_c1) + center2 = point.add(normPoint_c2) + elif dist.Length > 2 * radius: + return None + elif dist.Length == 2 * radius: + normPoint = point.add(findDistance(point, edge, False)) + dummy = DraftVecUtils.scale(normPoint.sub(point), 0.5) + cen = point.add(dummy) + circ = Part.Circle(cen, NORM, radius) + if circ: + return [circ] + else: + return None + else: + normPoint = point.add(findDistance(point, edge, False)) + normDist = DraftVecUtils.dist(normPoint, point) + dist = math.sqrt(radius**2 - (radius - normDist)**2) + centerNormVec = DraftVecUtils.scaleTo(point.sub(normPoint), radius) + edgeDir = edge.Vertexes[0].Point.sub(normPoint); edgeDir.normalize() + center1 = centerNormVec.add(normPoint.add(DraftVecUtils.scale(edgeDir, dist))) + center2 = centerNormVec.add(normPoint.add(DraftVecUtils.scale(edgeDir, -dist))) + circles = [] + if center1: + circ = Part.Circle(center1, NORM, radius) + if circ: + circles.append(circ) + if center2: + circ = Part.Circle(center2, NORM, radius) + if circ: + circles.append(circ) - if len(circles): - return circles - else: - return None + if len(circles): + return circles + else: + return None def circleFrom2PointsRadius(p1, p2, radius): - "circleFrom2PointsRadiust(Vector, Vector, radius)" - if DraftVecUtils.equals(p1, p2): return None + "circleFrom2PointsRadiust(Vector, Vector, radius)" + if DraftVecUtils.equals(p1, p2): return None - p1_p2 = Part.Line(p1, p2).toShape() - dist_p1p2 = DraftVecUtils.dist(p1, p1) - mid = findMidpoint(p1_p2) - if dist_p1p2 == 2*radius: - circle = Part.Circle(mid, norm, radius) - if circle: return [circle] - else: return None - dir = vec(p1_p2); dir.normalize() - perpDir = dir.cross(Vector(0,0,1)); perpDir.normailze() - dist = math.sqrt(radius**2 - (dist_p1p2 / 2.0)**2) - cen1 = Vector.add(mid, DraftVecUtils.scale(perpDir, dist)) - cen2 = Vector.add(mid, DraftVecUtils.scale(perpDir, -dist)) - circles = [] - if cen1: circles.append(Part.Circle(cen1, norm, radius)) - if cen2: circles.append(Part.Circle(cen2, norm, radius)) - if circles: return circles - else: return None + p1_p2 = Part.Line(p1, p2).toShape() + dist_p1p2 = DraftVecUtils.dist(p1, p1) + mid = findMidpoint(p1_p2) + if dist_p1p2 == 2*radius: + circle = Part.Circle(mid, norm, radius) + if circle: return [circle] + else: return None + dir = vec(p1_p2); dir.normalize() + perpDir = dir.cross(Vector(0,0,1)); perpDir.normailze() + dist = math.sqrt(radius**2 - (dist_p1p2 / 2.0)**2) + cen1 = Vector.add(mid, DraftVecUtils.scale(perpDir, dist)) + cen2 = Vector.add(mid, DraftVecUtils.scale(perpDir, -dist)) + circles = [] + if cen1: circles.append(Part.Circle(cen1, norm, radius)) + if cen2: circles.append(Part.Circle(cen2, norm, radius)) + if circles: return circles + else: return None @@ -1886,423 +1891,423 @@ def circleFrom2PointsRadius(p1, p2, radius): def outerSoddyCircle(circle1, circle2, circle3): - ''' - Computes the outer soddy circle for three tightly packed circles. - ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): - # Original Java code Copyright (rc) 2008 Werner Randelshofer - # Converted to python by Martin Buerbaum 2009 - # http://www.randelshofer.ch/treeviz/ - # Either Creative Commons Attribution 3.0, the MIT license, or the GNU Lesser General License LGPL. + ''' + Computes the outer soddy circle for three tightly packed circles. + ''' + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): + # Original Java code Copyright (rc) 2008 Werner Randelshofer + # Converted to python by Martin Buerbaum 2009 + # http://www.randelshofer.ch/treeviz/ + # Either Creative Commons Attribution 3.0, the MIT license, or the GNU Lesser General License LGPL. - A = circle1.Curve.Center - B = circle2.Curve.Center - C = circle3.Curve.Center + A = circle1.Curve.Center + B = circle2.Curve.Center + C = circle3.Curve.Center - ra = circle1.Curve.Radius - rb = circle2.Curve.Radius - rc = circle3.Curve.Radius + ra = circle1.Curve.Radius + rb = circle2.Curve.Radius + rc = circle3.Curve.Radius - # Solution using Descartes' theorem, as described here: - # http://en.wikipedia.org/wiki/Descartes%27_theorem - k1 = 1 / ra - k2 = 1 / rb - k3 = 1 / rc - k4 = abs(k1 + k2 + k3 - 2 * math.sqrt(k1 * k2 + k2 * k3 + k3 * k1)) + # Solution using Descartes' theorem, as described here: + # http://en.wikipedia.org/wiki/Descartes%27_theorem + k1 = 1 / ra + k2 = 1 / rb + k3 = 1 / rc + k4 = abs(k1 + k2 + k3 - 2 * math.sqrt(k1 * k2 + k2 * k3 + k3 * k1)) - q1 = (k1 + 0j) * (A.x + A.y * 1j) - q2 = (k2 + 0j) * (B.x + B.y * 1j) - q3 = (k3 + 0j) * (C.x + C.y * 1j) + q1 = (k1 + 0j) * (A.x + A.y * 1j) + q2 = (k2 + 0j) * (B.x + B.y * 1j) + q3 = (k3 + 0j) * (C.x + C.y * 1j) - temp = ((q1 * q2) + (q2 * q3) + (q3 * q1)) - q4 = q1 + q2 + q3 - ((2 + 0j) * cmath.sqrt(temp) ) + temp = ((q1 * q2) + (q2 * q3) + (q3 * q1)) + q4 = q1 + q2 + q3 - ((2 + 0j) * cmath.sqrt(temp) ) - z = q4 / (k4 + 0j) + z = q4 / (k4 + 0j) - # If the formula is not solveable, we return no circle. - if (not z or not (1 / k4)): - return None + # If the formula is not solveable, we return no circle. + if (not z or not (1 / k4)): + return None - X = -z.real - Y = -z.imag - print "Outer Soddy circle: " + str(X) + " " + str(Y) + "\n" # Debug + X = -z.real + Y = -z.imag + print "Outer Soddy circle: " + str(X) + " " + str(Y) + "\n" # Debug - # The Radius of the outer soddy circle can also be calculated with the following formula: - # radiusOuter = abs(r1*r2*r3 / (r1*r2 + r1*r3 + r2*r3 - 2 * math.sqrt(r1*r2*r3 * (r1+r2+r3)))) - circ = Part.Circle(Vector(X, Y, A.z), norm, 1 / k4) - return circ + # The Radius of the outer soddy circle can also be calculated with the following formula: + # radiusOuter = abs(r1*r2*r3 / (r1*r2 + r1*r3 + r2*r3 - 2 * math.sqrt(r1*r2*r3 * (r1+r2+r3)))) + circ = Part.Circle(Vector(X, Y, A.z), norm, 1 / k4) + return circ - else: - print "debug: outerSoddyCircle bad parameters!\n" - # FreeCAD.Console.PrintMessage("debug: outerSoddyCircle bad parameters!\n") - return None + else: + print "debug: outerSoddyCircle bad parameters!\n" + # FreeCAD.Console.PrintMessage("debug: outerSoddyCircle bad parameters!\n") + return None def innerSoddyCircle(circle1, circle2, circle3): - ''' - Computes the inner soddy circle for three tightly packed circles. - ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): - # Original Java code Copyright (rc) 2008 Werner Randelshofer - # Converted to python by Martin Buerbaum 2009 - # http://www.randelshofer.ch/treeviz/ + ''' + Computes the inner soddy circle for three tightly packed circles. + ''' + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): + # Original Java code Copyright (rc) 2008 Werner Randelshofer + # Converted to python by Martin Buerbaum 2009 + # http://www.randelshofer.ch/treeviz/ - A = circle1.Curve.Center - B = circle2.Curve.Center - C = circle3.Curve.Center + A = circle1.Curve.Center + B = circle2.Curve.Center + C = circle3.Curve.Center - ra = circle1.Curve.Radius - rb = circle2.Curve.Radius - rc = circle3.Curve.Radius + ra = circle1.Curve.Radius + rb = circle2.Curve.Radius + rc = circle3.Curve.Radius - # Solution using Descartes' theorem, as described here: - # http://en.wikipedia.org/wiki/Descartes%27_theorem - k1 = 1 / ra - k2 = 1 / rb - k3 = 1 / rc - k4 = abs(k1 + k2 + k3 + 2 * math.sqrt(k1 * k2 + k2 * k3 + k3 * k1)) + # Solution using Descartes' theorem, as described here: + # http://en.wikipedia.org/wiki/Descartes%27_theorem + k1 = 1 / ra + k2 = 1 / rb + k3 = 1 / rc + k4 = abs(k1 + k2 + k3 + 2 * math.sqrt(k1 * k2 + k2 * k3 + k3 * k1)) - q1 = (k1 + 0j) * (A.x + A.y * 1j) - q2 = (k2 + 0j) * (B.x + B.y * 1j) - q3 = (k3 + 0j) * (C.x + C.y * 1j) + q1 = (k1 + 0j) * (A.x + A.y * 1j) + q2 = (k2 + 0j) * (B.x + B.y * 1j) + q3 = (k3 + 0j) * (C.x + C.y * 1j) - temp = ((q1 * q2) + (q2 * q3) + (q3 * q1)) - q4 = q1 + q2 + q3 + ((2 + 0j) * cmath.sqrt(temp) ) + temp = ((q1 * q2) + (q2 * q3) + (q3 * q1)) + q4 = q1 + q2 + q3 + ((2 + 0j) * cmath.sqrt(temp) ) - z = q4 / (k4 + 0j) + z = q4 / (k4 + 0j) - # If the formula is not solveable, we return no circle. - if (not z or not (1 / k4)): - return None + # If the formula is not solveable, we return no circle. + if (not z or not (1 / k4)): + return None - X = z.real - Y = z.imag - print "Outer Soddy circle: " + str(X) + " " + str(Y) + "\n" # Debug + X = z.real + Y = z.imag + print "Outer Soddy circle: " + str(X) + " " + str(Y) + "\n" # Debug - # The Radius of the inner soddy circle can also be calculated with the following formula: - # radiusInner = abs(r1*r2*r3 / (r1*r2 + r1*r3 + r2*r3 + 2 * math.sqrt(r1*r2*r3 * (r1+r2+r3)))) - circ = Part.Circle(Vector(X, Y, A.z), norm, 1 / k4) - return circ + # The Radius of the inner soddy circle can also be calculated with the following formula: + # radiusInner = abs(r1*r2*r3 / (r1*r2 + r1*r3 + r2*r3 + 2 * math.sqrt(r1*r2*r3 * (r1+r2+r3)))) + circ = Part.Circle(Vector(X, Y, A.z), norm, 1 / k4) + return circ - else: - print "debug: innerSoddyCircle bad parameters!\n" - # FreeCAD.Console.PrintMessage("debug: innerSoddyCircle bad parameters!\n") - return None + else: + print "debug: innerSoddyCircle bad parameters!\n" + # FreeCAD.Console.PrintMessage("debug: innerSoddyCircle bad parameters!\n") + return None def circleFrom3CircleTangents(circle1, circle2, circle3): - ''' - http://en.wikipedia.org/wiki/Problem_of_Apollonius#Inversive_methods - http://mathworld.wolfram.com/ApolloniusCircle.html - http://mathworld.wolfram.com/ApolloniusProblem.html - ''' + ''' + http://en.wikipedia.org/wiki/Problem_of_Apollonius#Inversive_methods + http://mathworld.wolfram.com/ApolloniusCircle.html + http://mathworld.wolfram.com/ApolloniusProblem.html + ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): - int12 = findIntersection(circle1, circle2, True, True) - int23 = findIntersection(circle2, circle3, True, True) - int31 = findIntersection(circle3, circle1, True, True) + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle3.Curve, Part.Circle): + int12 = findIntersection(circle1, circle2, True, True) + int23 = findIntersection(circle2, circle3, True, True) + int31 = findIntersection(circle3, circle1, True, True) - if int12 and int23 and int31: - if len(int12) == 1 and len(int23) == 1 and len(int31) == 1: - # Only one intersection with each circle. - # => "Soddy Circle" - 2 solutions. - # http://en.wikipedia.org/wiki/Problem_of_Apollonius#Mutually_tangent_given_circles:_Soddy.27s_circles_and_Descartes.27_theorem - # http://mathworld.wolfram.com/SoddyCircles.html - # http://mathworld.wolfram.com/InnerSoddyCenter.html - # http://mathworld.wolfram.com/OuterSoddyCenter.html + if int12 and int23 and int31: + if len(int12) == 1 and len(int23) == 1 and len(int31) == 1: + # Only one intersection with each circle. + # => "Soddy Circle" - 2 solutions. + # http://en.wikipedia.org/wiki/Problem_of_Apollonius#Mutually_tangent_given_circles:_Soddy.27s_circles_and_Descartes.27_theorem + # http://mathworld.wolfram.com/SoddyCircles.html + # http://mathworld.wolfram.com/InnerSoddyCenter.html + # http://mathworld.wolfram.com/OuterSoddyCenter.html - r1 = circle1.Curve.Radius - r2 = circle2.Curve.Radius - r3 = circle3.Curve.Radius - outerSoddy = outerSoddyCircle(circle1, circle2, circle3) -# print str(outerSoddy) + "\n" # Debug + r1 = circle1.Curve.Radius + r2 = circle2.Curve.Radius + r3 = circle3.Curve.Radius + outerSoddy = outerSoddyCircle(circle1, circle2, circle3) +# print str(outerSoddy) + "\n" # Debug - innerSoddy = innerSoddyCircle(circle1, circle2, circle3) -# print str(innerSoddy) + "\n" # Debug + innerSoddy = innerSoddyCircle(circle1, circle2, circle3) +# print str(innerSoddy) + "\n" # Debug - circles = [] - if outerSoddy: - circles.append(outerSoddy) - if innerSoddy: - circles.append(innerSoddy) - return circles + circles = [] + if outerSoddy: + circles.append(outerSoddy) + if innerSoddy: + circles.append(innerSoddy) + return circles - # @todo Calc all 6 homothetic centers. - # @todo Create 3 lines from the inner and 4 from the outer h. center. - # @todo Calc. the 4 inversion poles of these lines for each circle. - # @todo Calc. the radical center of the 3 circles. - # @todo Calc. the intersection points (max. 8) of 4 lines (trough each inversion pole and the radical center) with the circle. - # This gives us all the tangent points. - else: - # Some circles are inside each other or an error has occured. - return None + # @todo Calc all 6 homothetic centers. + # @todo Create 3 lines from the inner and 4 from the outer h. center. + # @todo Calc. the 4 inversion poles of these lines for each circle. + # @todo Calc. the radical center of the 3 circles. + # @todo Calc. the intersection points (max. 8) of 4 lines (trough each inversion pole and the radical center) with the circle. + # This gives us all the tangent points. + else: + # Some circles are inside each other or an error has occured. + return None - else: - print "debug: circleFrom3CircleTangents bad parameters!\n" - # FreeCAD.Console.PrintMessage("debug: circleFrom3CircleTangents bad parameters!\n") - return None + else: + print "debug: circleFrom3CircleTangents bad parameters!\n" + # FreeCAD.Console.PrintMessage("debug: circleFrom3CircleTangents bad parameters!\n") + return None def linearFromPoints (p1, p2): - ''' - Calculate linear equation from points. - Calculate the slope and offset parameters of the linear equation of a line defined by two points. + ''' + Calculate linear equation from points. + Calculate the slope and offset parameters of the linear equation of a line defined by two points. - Linear equation: - y = m * x + b - m = dy / dx - m ... Slope - b ... Offset (point where the line intersects the y axis) - dx/dy ... Delta x and y. Using both as a vector results in a non-offset direction vector. - ''' - if isinstance(p1, Vector) and isinstance(p2, Vector): - line = {} - line['dx'] = (p2.x - p1.x) - line['dy'] = (p2.y - p1.y) - line['slope'] = line['dy'] / line['dx'] - line['offset'] = p1.y - slope * p1.x - return line - else: - return None + Linear equation: + y = m * x + b + m = dy / dx + m ... Slope + b ... Offset (point where the line intersects the y axis) + dx/dy ... Delta x and y. Using both as a vector results in a non-offset direction vector. + ''' + if isinstance(p1, Vector) and isinstance(p2, Vector): + line = {} + line['dx'] = (p2.x - p1.x) + line['dy'] = (p2.y - p1.y) + line['slope'] = line['dy'] / line['dx'] + line['offset'] = p1.y - slope * p1.x + return line + else: + return None def determinant (mat,n): - ''' - determinant(matrix,int) - Determinat function. Returns the determinant - of a n-matrix. It recursively expands the minors. - ''' - matTemp = [[0.0,0.0,0.0],[0.0,0.0,0.0],[0.0,0.0,0.0]] - if (n > 1): - if n == 2: - d = mat[0][0] * mat[1][1] - mat[1][0] * mat[0][1] - else: - d = 0.0 - for j1 in range(n): - # Create minor - for i in range(1, n): - j2 = 0 - for j in range(n): - if j == j1: - continue - matTemp[i-1][j2] = mat[i][j] - j2 += 1 - d += (-1.0)**(1.0 + j1 + 1.0) * mat[0][j1] * determinant(matTemp, n-1) - return d - else: - return 0 + ''' + determinant(matrix,int) - Determinat function. Returns the determinant + of a n-matrix. It recursively expands the minors. + ''' + matTemp = [[0.0,0.0,0.0],[0.0,0.0,0.0],[0.0,0.0,0.0]] + if (n > 1): + if n == 2: + d = mat[0][0] * mat[1][1] - mat[1][0] * mat[0][1] + else: + d = 0.0 + for j1 in range(n): + # Create minor + for i in range(1, n): + j2 = 0 + for j in range(n): + if j == j1: + continue + matTemp[i-1][j2] = mat[i][j] + j2 += 1 + d += (-1.0)**(1.0 + j1 + 1.0) * mat[0][j1] * determinant(matTemp, n-1) + return d + else: + return 0 - + def findHomotheticCenterOfCircles(circle1, circle2): - ''' - findHomotheticCenterOfCircles(circle1, circle2) - Calculates the homothetic center(s) of two circles. + ''' + findHomotheticCenterOfCircles(circle1, circle2) + Calculates the homothetic center(s) of two circles. - http://en.wikipedia.org/wiki/Homothetic_center - http://mathworld.wolfram.com/HomotheticCenter.html - ''' + http://en.wikipedia.org/wiki/Homothetic_center + http://mathworld.wolfram.com/HomotheticCenter.html + ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): - if DraftVecUtils.equals(circle1.Curve.Center, circle2.Curve.Center): - return None + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): + if DraftVecUtils.equals(circle1.Curve.Center, circle2.Curve.Center): + return None - cen1_cen2 = Part.Line(circle1.Curve.Center, circle2.Curve.Center).toShape() - cenDir = vec(cen1_cen2); cenDir.normalize() + cen1_cen2 = Part.Line(circle1.Curve.Center, circle2.Curve.Center).toShape() + cenDir = vec(cen1_cen2); cenDir.normalize() - # Get the perpedicular vector. - perpCenDir = cenDir.cross(Vector(0,0,1)); perpCenDir.normalize() + # Get the perpedicular vector. + perpCenDir = cenDir.cross(Vector(0,0,1)); perpCenDir.normalize() - # Get point on first circle - p1 = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, circle1.Curve.Radius)) + # Get point on first circle + p1 = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, circle1.Curve.Radius)) - centers = [] - # Calculate inner homothetic center - # Get point on second circle - p2_inner = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, -circle1.Curve.Radius)) - hCenterInner = DraftVecUtils.intersect(circle1.Curve.Center, circle2.Curve.Center, p1, p2_inner, True, True) - if hCenterInner: - centers.append(hCenterInner) + centers = [] + # Calculate inner homothetic center + # Get point on second circle + p2_inner = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, -circle1.Curve.Radius)) + hCenterInner = DraftVecUtils.intersect(circle1.Curve.Center, circle2.Curve.Center, p1, p2_inner, True, True) + if hCenterInner: + centers.append(hCenterInner) - # Calculate outer homothetic center (only exists of the circles have different radii) - if circle1.Curve.Radius != circle2.Curve.Radius: - # Get point on second circle - p2_outer = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, circle1.Curve.Radius)) - hCenterOuter = DraftVecUtils.intersect(circle1.Curve.Center, circle2.Curve.Center, p1, p2_outer, True, True) - if hCenterOuter: - centers.append(hCenterOuter) + # Calculate outer homothetic center (only exists of the circles have different radii) + if circle1.Curve.Radius != circle2.Curve.Radius: + # Get point on second circle + p2_outer = Vector.add(circle1.Curve.Center, DraftVecUtils.scale(perpCenDir, circle1.Curve.Radius)) + hCenterOuter = DraftVecUtils.intersect(circle1.Curve.Center, circle2.Curve.Center, p1, p2_outer, True, True) + if hCenterOuter: + centers.append(hCenterOuter) - if len(centers): - return centers - else: - return None + if len(centers): + return centers + else: + return None - else: - print "debug: findHomotheticCenterOfCircles bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: findHomotheticCenterOfCirclescleFrom3tan bad parameters!\n") - return None + else: + print "debug: findHomotheticCenterOfCircles bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: findHomotheticCenterOfCirclescleFrom3tan bad parameters!\n") + return None def findRadicalAxis(circle1, circle2): - ''' - Calculates the radical axis of two circles. - On the radical axis (also called power line) of two circles any - tangents drawn from a point on the axis to both circles have the same length. + ''' + Calculates the radical axis of two circles. + On the radical axis (also called power line) of two circles any + tangents drawn from a point on the axis to both circles have the same length. - http://en.wikipedia.org/wiki/Radical_axis - http://mathworld.wolfram.com/RadicalLine.html + http://en.wikipedia.org/wiki/Radical_axis + http://mathworld.wolfram.com/RadicalLine.html - @sa findRadicalCenter - ''' + @sa findRadicalCenter + ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): - if DraftVecUtils.equals(circle1.Curve.Center, circle2.Curve.Center): - return None - r1 = circle1.Curve.Radius - r2 = circle1.Curve.Radius - cen1 = circle1.Curve.Center - # dist .. the distance from cen1 to cen2. - dist = DraftVecUtils.dist(cen1, circle2.Curve.Center) - cenDir = cen1.sub(circle2.Curve.Center); cenDir.normalize() + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): + if DraftVecUtils.equals(circle1.Curve.Center, circle2.Curve.Center): + return None + r1 = circle1.Curve.Radius + r2 = circle1.Curve.Radius + cen1 = circle1.Curve.Center + # dist .. the distance from cen1 to cen2. + dist = DraftVecUtils.dist(cen1, circle2.Curve.Center) + cenDir = cen1.sub(circle2.Curve.Center); cenDir.normalize() - # Get the perpedicular vector. - perpCenDir = cenDir.cross(Vector(0,0,1)); perpCenDir.normalize() + # Get the perpedicular vector. + perpCenDir = cenDir.cross(Vector(0,0,1)); perpCenDir.normalize() - # J ... The radical center. - # K ... The point where the cadical axis crosses the line of cen1->cen2. - # k1 ... Distance from cen1 to K. - # k2 ... Distance from cen2 to K. - # dist = k1 + k2 + # J ... The radical center. + # K ... The point where the cadical axis crosses the line of cen1->cen2. + # k1 ... Distance from cen1 to K. + # k2 ... Distance from cen2 to K. + # dist = k1 + k2 - k1 = (dist + (r1^2 - r2^2) / dist) / 2.0 - #k2 = dist - k1 + k1 = (dist + (r1^2 - r2^2) / dist) / 2.0 + #k2 = dist - k1 - K = Vector.add(cen1, DraftVecUtils.scale(cenDir, k1)) + K = Vector.add(cen1, DraftVecUtils.scale(cenDir, k1)) - # K_ .. A point somewhere between K and J (actually with a distance of 1 unit from K). - K_ = Vector,add(K, perpCenDir) + # K_ .. A point somewhere between K and J (actually with a distance of 1 unit from K). + K_ = Vector,add(K, perpCenDir) - radicalAxis = Part.Line(K, Vector.add(origin, dir)) + radicalAxis = Part.Line(K, Vector.add(origin, dir)) - if radicalAxis: - return radicalAxis - else: - return None - else: - print "debug: findRadicalAxis bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: findRadicalAxis bad parameters!\n") - return None + if radicalAxis: + return radicalAxis + else: + return None + else: + print "debug: findRadicalAxis bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: findRadicalAxis bad parameters!\n") + return None def findRadicalCenter(circle1, circle2, circle3): - ''' - findRadicalCenter(circle1, circle2, circle3): - Calculates the radical center (also called the power center) of three circles. - It is the intersection point of the three radical axes of the pairs of circles. + ''' + findRadicalCenter(circle1, circle2, circle3): + Calculates the radical center (also called the power center) of three circles. + It is the intersection point of the three radical axes of the pairs of circles. - http://en.wikipedia.org/wiki/Power_center_(geometry) - http://mathworld.wolfram.com/RadicalCenter.html + http://en.wikipedia.org/wiki/Power_center_(geometry) + http://mathworld.wolfram.com/RadicalCenter.html - @sa findRadicalAxis - ''' + @sa findRadicalAxis + ''' - if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): - radicalAxis12 = findRadicalAxis(circle1, circle2) - radicalAxis23 = findRadicalAxis(circle1, circle2) + if isinstance(circle1.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): + radicalAxis12 = findRadicalAxis(circle1, circle2) + radicalAxis23 = findRadicalAxis(circle1, circle2) - if not radicalAxis12 or not radicalAxis23: - # No radical center could be calculated. - return None + if not radicalAxis12 or not radicalAxis23: + # No radical center could be calculated. + return None - int = findIntersection(radicalAxis12, radicalAxis23, True, True) + int = findIntersection(radicalAxis12, radicalAxis23, True, True) - if int: - return int - else: - # No radical center could be calculated. - return None - else: - print "debug: findRadicalCenter bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: findRadicalCenter bad parameters!\n") - return None + if int: + return int + else: + # No radical center could be calculated. + return None + else: + print "debug: findRadicalCenter bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: findRadicalCenter bad parameters!\n") + return None def pointInversion(circle, point): - ''' - pointInversion(Circle, Vector) + ''' + pointInversion(Circle, Vector) - Circle inversion of a point. - Will calculate the inversed point an return it. - If the given point is equal to the center of the circle "None" will be returned. + Circle inversion of a point. + Will calculate the inversed point an return it. + If the given point is equal to the center of the circle "None" will be returned. - See also: - http://en.wikipedia.org/wiki/Inversive_geometry - ''' + See also: + http://en.wikipedia.org/wiki/Inversive_geometry + ''' - if isinstance(circle.Curve, Part.Circle) and isinstance(point, FreeCAD.Vector): - cen = circle.Curve.Center - rad = circle.Curve.Radius + if isinstance(circle.Curve, Part.Circle) and isinstance(point, FreeCAD.Vector): + cen = circle.Curve.Center + rad = circle.Curve.Radius - if DraftVecUtils.equals(cen, point): - return None + if DraftVecUtils.equals(cen, point): + return None - # Inverse the distance of the point - # dist(cen -> P) = r^2 / dist(cen -> invP) + # Inverse the distance of the point + # dist(cen -> P) = r^2 / dist(cen -> invP) - dist = DraftVecUtils.dist(point, cen) - invDist = rad**2 / d + dist = DraftVecUtils.dist(point, cen) + invDist = rad**2 / d - invPoint = Vector(0, 0, point.z) - invPoint.x = cen.x + (point.x - cen.x) * invDist / dist; - invPoint.y = cen.y + (point.y - cen.y) * invDist / dist; + invPoint = Vector(0, 0, point.z) + invPoint.x = cen.x + (point.x - cen.x) * invDist / dist; + invPoint.y = cen.y + (point.y - cen.y) * invDist / dist; - return invPoint + return invPoint - else: - print "debug: pointInversion bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: pointInversion bad parameters!\n") - return None + else: + print "debug: pointInversion bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: pointInversion bad parameters!\n") + return None def polarInversion(circle, edge): - ''' - polarInversion(circle, edge): - Returns the inversion pole of a line. - edge ... The polar. - i.e. The nearest point on the line is inversed. + ''' + polarInversion(circle, edge): + Returns the inversion pole of a line. + edge ... The polar. + i.e. The nearest point on the line is inversed. - http://mathworld.wolfram.com/InversionPole.html - ''' + http://mathworld.wolfram.com/InversionPole.html + ''' - if isinstance(circle.Curve, Part.Circle) and isinstance(edge.Curve, Part.Line): - nearest = circle.Curve.Center.add(findDistance(circle.Curve.Center, edge, False)) - if nearest: - inversionPole = pointInversion(circle, nearest) - if inversionPole: - return inversionPole + if isinstance(circle.Curve, Part.Circle) and isinstance(edge.Curve, Part.Line): + nearest = circle.Curve.Center.add(findDistance(circle.Curve.Center, edge, False)) + if nearest: + inversionPole = pointInversion(circle, nearest) + if inversionPole: + return inversionPole - else: - print "debug: circleInversionPole bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: circleInversionPole bad parameters!\n") - return None + else: + print "debug: circleInversionPole bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: circleInversionPole bad parameters!\n") + return None def circleInversion(circle, circle2): - ''' - pointInversion(Circle, Circle) + ''' + pointInversion(Circle, Circle) - Circle inversion of a circle. - ''' - if isinstance(circle.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): - cen1 = circle.Curve.Center - rad1 = circle.Curve.Radius + Circle inversion of a circle. + ''' + if isinstance(circle.Curve, Part.Circle) and isinstance(circle2.Curve, Part.Circle): + cen1 = circle.Curve.Center + rad1 = circle.Curve.Radius - if DraftVecUtils.equals(cen1, point): - return None + if DraftVecUtils.equals(cen1, point): + return None - invCen2 = Inversion(circle, circle2.Curve.Center) + invCen2 = Inversion(circle, circle2.Curve.Center) - pointOnCircle2 = Vector.add(circle2.Curve.Center, Vector(circle2.Curve.Radius, 0, 0)) - invPointOnCircle2 = Inversion(circle, pointOnCircle2) + pointOnCircle2 = Vector.add(circle2.Curve.Center, Vector(circle2.Curve.Radius, 0, 0)) + invPointOnCircle2 = Inversion(circle, pointOnCircle2) - return Part.Circle(invCen2, norm, DraftVecUtils.dist(invCen2, invPointOnCircle2)) + return Part.Circle(invCen2, norm, DraftVecUtils.dist(invCen2, invPointOnCircle2)) - else: - print "debug: circleInversion bad parameters!\n" - FreeCAD.Console.PrintMessage("debug: circleInversion bad parameters!\n") - return None + else: + print "debug: circleInversion bad parameters!\n" + FreeCAD.Console.PrintMessage("debug: circleInversion bad parameters!\n") + return None