Skip to content

P
pyMKcam
Commit 73cb1bb0 authored by Whitham D. Reeve II's avatar Whitham D. Reeve II
Browse files
Options

Removed old commented out code and converted PolygonExtractor.py to new Point tuple style.

parent d45b6dcb
Show whitespace changes
Inline Side-by-side
Showing with 90 additions and 270 deletions
pycam/Geometry/Letters.py
View file @ 73cb1bb0
......@@ -119,7 +119,6 @@ class Charset(object):
for character in line:
if character == " ":
base = padd(base, (word_spacing, 0, 0))
#base = base.add(Point(word_spacing, 0, 0))
elif character in self.letters.keys():
charset_letter = self.letters[character]
new_model = ContourModel()
......@@ -133,14 +132,11 @@ class Charset(object):
line_height = max(line_height, charset_letter.maxy())
# shift the base position
base = padd(base, (charset_letter.maxx() + letter_spacing, 0, 0))
#base = base.add((charset_letter.maxx() + letter_spacing, 0, 0))
else:
# unknown character - add a small whitespace
base = padd(base, (letter_spacing, 0, 0))
#base = base.add(Point(letter_spacing, 0, 0))
# go to the next line
base = (origin[0], base[1] - line_height * line_factor, origin[2])
#base = Point(origin.x, base.y - line_height * line_factor, origin.z)
if not current_line.maxx is None:
if align == TEXT_ALIGN_CENTER:
current_line.shift(-current_line.maxx / 2, 0, 0)
......
pycam/Geometry/Line.py
View file @ 73cb1bb0
......@@ -54,18 +54,15 @@ class Line(IDGenerator, TransformableContainer):
def vector(self):
if self._vector is None:
self._vector = psub(self.p2, self.p1)
#self._vector = self.p2.sub(self.p1)
return self._vector
@property
def dir(self):
return pnormalized(self.vector)
#return self.vector.normalized()
@property
def len(self):
return pnorm(self.vector)
#return self.vector.norm
@property
def minx(self):
......@@ -125,8 +122,6 @@ class Line(IDGenerator, TransformableContainer):
def next(self):
yield "p1"
#yield self.p1
#yield lambda x: self.p2 = x
yield "p2"
def get_children_count(self):
......@@ -147,13 +142,11 @@ class Line(IDGenerator, TransformableContainer):
def point_with_length_multiply(self, l):
return padd(self.p1, pmul(self.dir, l*self.len))
#return self.p1.add(self.dir.mul(l*self.len))
def get_length_line(self, length):
""" return a line with the same direction and the specified length
"""
return Line(self.p1, padd(self.p1, pmul(self.dir, length)))
#return Line(self.p1, self.p1.add(self.dir.mul(length)))
def closest_point(self, p):
v = self.dir
......@@ -161,13 +154,10 @@ class Line(IDGenerator, TransformableContainer):
# for zero-length lines
return self.p1
l = pdot(self.p1, v) - pdot(p, v)
#l = self.p1.dot(v) - p.dot(v)
return psub(self.p1, pmul(v, l))
#return self.p1.sub(v.mul(l))
def dist_to_point_sq(self, p):
return pnormsq(psub(p, self.closes_point(p)))
#return p.sub(self.closest_point(p)).normsq
def dist_to_point(self, p):
return sqrt(self.dist_to_point_sq(p))
......@@ -178,9 +168,7 @@ class Line(IDGenerator, TransformableContainer):
return True
dir1 = pnormalized(psub(p, self.p1))
#dir1 = p.sub(self.p1).normalized()
dir2 = pnormalized(psub(self.p2, p))
#dir2 = self.p2.sub(p).normalized()
# True if the two parts of the line have the same direction or if the
# point is self.p1 or self.p2.
return (dir1 == dir2 == self.dir) or (dir1 is None) or (dir2 is None)
......@@ -210,19 +198,14 @@ class Line(IDGenerator, TransformableContainer):
"""
x1, x2, x3, x4 = self.p1, self.p2, line.p1, line.p2
a = psub(x2, x1)
#a = x2.sub(x1)
b = psub(x4, x3)
#b = x4.sub(x3)
c = psub(x3, x1)
#c = x3.sub(x1)
# see http://mathworld.wolfram.com/Line-LineIntersection.html (24)
try:
factor = pdot(pcross(c, b), pcross(a, b)) / pnormsq(pcross(a, b))
#factor = c.cross(b).dot(a.cross(b)) / a.cross(b).normsq
except ZeroDivisionError:
# lines are parallel
# check if they are _one_ line
#if a.cross(c).norm != 0:
if pnorm(pcross(a,c)) != 0:
# the lines are parallel with a distance
return None, None
......@@ -232,7 +215,6 @@ class Line(IDGenerator, TransformableContainer):
candidates.append((x3, pnorm(c) / pnorm(a)))
elif self.is_point_inside(x4):
candidates.append((x4, pnorm(psub(line.p2, self.p1)) / pnorm(a)))
#candidates.append((x4, line.p2.sub(self.p1).norm / a.norm))
elif line.is_point_inside(x1):
candidates.append((x1, 0))
elif line.is_point_inside(x2):
......@@ -244,7 +226,6 @@ class Line(IDGenerator, TransformableContainer):
return candidates[0]
if infinite_lines or (-epsilon <= factor <= 1 + epsilon):
intersection = padd(x1, pmul(a, factor))
#intersection = x1.add(a.mul(factor))
# check if the intersection is between x3 and x4
if infinite_lines:
return intersection, factor
......
pycam/Geometry/Model.py
View file @ 73cb1bb0
......@@ -980,7 +980,6 @@ class PolygonGroup(object):
line_distances = []
for line in self.lines:
cross_product = pcross(line.dir, psub(point, line.p1))
#cross_product = line.dir.cross(point.sub(line.p1))
if cross_product[2] > 0:
close_points = []
close_point = line.closest_point(point)
......@@ -991,9 +990,7 @@ class PolygonGroup(object):
close_points.append(close_point)
for p in close_points:
direction = psub(point, p)
#direction = point.sub(p)
dist = pnorm(direction)
#dist = direction.norm
line_distances.append(dist)
elif cross_product.z == 0:
# the point is on the line
......@@ -1071,7 +1068,6 @@ class Rectangle(IDGenerator, TransformableContainer):
orders = ((p1, p2, p3, p4), (p1, p2, p4, p3), (p1, p3, p2, p4),
(p1, p3, p4, p2), (p1, p4, p2, p3), (p1, p4, p3, p2))
for order in orders:
#if abs(order[0].sub(order[2]).norm - order[1].sub(order[3]).norm) < epsilon:
if abs(pnorm(psub(order[0], order[2])) - pnorm(psub(order[1], order[3]))) < epsilon:
t1 = Triangle(order[0], order[1], order[2])
t2 = Triangle(order[2], order[3], order[0])
......@@ -1101,10 +1097,6 @@ class Rectangle(IDGenerator, TransformableContainer):
return (self.p1, self.p2, self.p3, self.p4)
def next(self):
#yield self.p1
#yield self.p2
#yield self.p3
#yield self.p4
yield "p1"
yield "p2"
yield "p3"
......@@ -1141,7 +1133,6 @@ class Rectangle(IDGenerator, TransformableContainer):
log.error("Invalid number of vertices: %s" % unique_vertices)
return None
if abs(pnorm(psub(unique_verticies[0], unique_verticies[1])) - pnorm(psub(shared_vertices[0], shared_vertices[1]))) < epsilon:
#if abs(unique_vertices[0].sub(unique_vertices[1]).norm - shared_vertices[0].sub(shared_vertices[1]).norm) < epsilon:
try:
return Rectangle(unique_vertices[0], unique_vertices[1],
shared_vertices[0], shared_vertices[1],
......
pycam/Geometry/Plane.py
View file @ 73cb1bb0
......@@ -77,13 +77,10 @@ class Plane(IDGenerator, TransformableContainer):
if direction is None:
return (None, INFINITE)
denom = pdot(self.n, direction)
#denom = self.n.dot(direction)
if denom == 0:
return (None, INFINITE)
l = -(pdot(self.n, point) - pdot(self.n, self.p)) / denom
#l = -(self.n.dot(point) - self.n.dot(self.p)) / denom
cp = padd(point, pmul(direction, l))
#cp = point.add(direction.mul(l))
return (cp, l)
def intersect_triangle(self, triangle, counter_clockwise=False):
......@@ -121,8 +118,6 @@ class Plane(IDGenerator, TransformableContainer):
if collision_line.len == 0:
return collision_line
cross = pcross(self.n, collision_line.dir)
#cross = self.n.cross(collision_line.dir)
#if (cross.dot(triangle.normal) < 0) == bool(not counter_clockwise):
if (pdot(cross, triangle.normal) < 0) == bool(not counter_clockwise):
# anti-clockwise direction -> revert the direction of the line
collision_line = Line(collision_line.p2, collision_line.p1)
......
pycam/Geometry/PointKdtree.py
View file @ 73cb1bb0
......@@ -47,7 +47,6 @@ class PointKdtree(kdtree):
return dx*dx+dy*dy+dz*dz
def Point(self, x, y, z):
#return Point(x,y,z)
if self._n:
n = self._n
n.bound = (x, y, z)
......
pycam/Geometry/Polygon.py
View file @ 73cb1bb0
......@@ -72,7 +72,6 @@ class PolygonInTree(IDGenerator):
def get_cost(self, other):
return pnorm(psub(other.start, self.end))
#return other.start.sub(self.end).norm
class PolygonPositionSorter(object):
......@@ -251,7 +250,6 @@ class Polygon(TransformableContainer):
self._update_limits(line.p2)
elif self._points[-1] == line.p1:
# the new Line can be added to the end of the polygon
#if line.dir == self._points[-1].sub(self._points[-2]).normalized():
if line.dir == pnormalized(psub(self._points[-1], self._points[-2])):
# Remove the last point, if the previous point combination
# is in line with the new Line. This avoids unnecessary
......@@ -266,7 +264,6 @@ class Polygon(TransformableContainer):
self.reset_cache()
else:
# the new Line can be added to the beginning of the polygon
#if (len(self._points) > 1) and (line.dir == self._points[1].sub(self._points[0]).normalized()):
if (len(self._points) > 1) and (line.dir == pnormalized(psub(self._points[1], self._points[0]))):
# Avoid points on straight lines - see above.
self._points.pop(0)
......@@ -323,8 +320,6 @@ class Polygon(TransformableContainer):
return False
def next(self):
#for idx,point in enumerate(self._points):
# yield (idx, "_points")
yield "_points"
yield self.plane
......@@ -418,16 +413,13 @@ class Polygon(TransformableContainer):
return None
else:
return pdiv(padd(self._points[index], self._points[(index + 1) % len(self._points)]), 2)
#return self._points[index].add(self._points[(index + 1) % len(self._points)]).div(2)
def get_lengths(self):
result = []
for index in range(len(self._points) - 1):
result.append(pnorm(psub(self._points[index + 1], self._points[index])))
#result.append(self._points[index + 1].sub(self._points[index]).norm)
if self.is_closed:
result.append(pnorm(psub(self._points[0], self._points[-1])))
#result.append(self._points[0].sub(self._points[-1]).norm)
return result
def get_max_inside_distance(self):
......@@ -436,13 +428,11 @@ class Polygon(TransformableContainer):
if len(self._points) < 2:
return None
distance = pnorm(psub(self._points[1], self._points[0]))
#distance = self._points[1].sub(self._points[0]).norm
for p1 in self._points:
for p2 in self._points:
if p1 is p2:
continue
distance = max(distance, pnorm(psub(p2, p1)))
#distance = max(distance, p2.sub(p1).norm)
return distance
def is_outer(self):
......@@ -595,15 +585,11 @@ class Polygon(TransformableContainer):
p1 = self._points[index]
p2 = self._points[(index + 1) % len(self._points)]
cross_offset = pnormalized(pcross(psub(p2, p1), self.plane.n))
#cross_offset = p2.sub(p1).cross(self.plane.n).normalized()
bisector_normalized = self.get_bisector(index)
factor = pdot(cross_offset, bisector_normalized)
#factor = cross_offset.dot(bisector_normalized)
if factor != 0:
bisector_sized = pmul(bisector_normalized, offset / factor)
#bisector_sized = bisector_normalized.mul(offset / factor)
return padd(p1, bisector_sized)
#return p1.add(bisector_sized)
else:
return p2
if offset * 2 >= self.get_max_inside_distance():
......@@ -615,7 +601,6 @@ class Polygon(TransformableContainer):
max_dist = 1000 * epsilon
def test_point_near(p, others):
for o in others:
#if p.sub(o).norm < max_dist:
if pnorm(psub(p, o)) < max_dist:
return True
return False
......@@ -626,10 +611,7 @@ class Polygon(TransformableContainer):
p1 = points[index]
p2 = points[next_index]
diff = psub(p2, p1)
#diff = p2.sub(p1)
old_dir = pnormalized(psub(self._points[next_index], self._points[index]))
#old_dir = self._points[next_index].sub(self._points[index]).normalized()
#if diff.normalized() != old_dir:
if pnormalized(diff) != old_dir:
# the direction turned around
if pnorm(diff) > max_dist:
......@@ -663,7 +645,6 @@ class Polygon(TransformableContainer):
# no lines are left
print "out 2"
return []
#if prev_line.p2.sub(next_line.p1).norm > max_dist:
if pnorm(psub(prev_line.p2, next_line.p1)) > max_dist:
cp, dist = prev_line.get_intersection(next_line)
else:
......@@ -711,9 +692,7 @@ class Polygon(TransformableContainer):
# maybe we have been here before
if not cp in split_points:
split_points.append(cp)
#elif (cp.sub(line.p1).norm < max_dist) or (cp.sub(line.p2).norm < max_dist):
elif (pnorm(psub(cp, line.p1)) < max_dist) or (pnorm(psub(cp, line.p2)) < max_dist):
#if cp.sub(line.p1).norm < cp.sub(line.p2).norm:
if pnorm(psub(cp, lines.p1)) < pnorm(psub(cp, line.p2)):
non_reversed[index] = Line(cp, line.p2)
else:
......@@ -790,7 +769,6 @@ class Polygon(TransformableContainer):
if len(group) <= 2:
continue
poly = Polygon(self.plane)
#print "**************************************"
for line in group:
try:
poly.append(line)
......@@ -871,15 +849,11 @@ class Polygon(TransformableContainer):
p1 = self._points[index]
p2 = self._points[(index + 1) % len(self._points)]
cross_offset = pnormalized(pcross(psub(p2, p1), self.plane.n))
#cross_offset = p2.sub(p1).cross(self.plane.n).normalized()
bisector_normalized = self.get_bisector(index)
factor = pdot(cross_offset, bisector_normalized)
#factor = cross_offset.dot(bisector_normalized)
if factor != 0:
bisector_sized = pmul(bisector_normalized, offset / factor)
#bisector_sized = bisector_normalized.mul(offset / factor)
return padd(p1, bisector_sized)
#return p1.add(bisector_sized)
else:
return p2
def simplify_polygon_intersections(lines):
......@@ -1056,16 +1030,13 @@ class Polygon(TransformableContainer):
return Line(line.p1, line.p2)
else:
cross_offset = pmul(pnormalized(pcross(line.dir, self.plane.n)), offset)
#cross_offset = line.dir.cross(self.plane.n).normalized().mul(offset)
# Prolong the line at the beginning and at the end - to allow
# overlaps. Use factor "2" to take care for star-like structure
# where a complete convex triangle would get cropped (two lines
# get lost instead of just one). Use the "abs" value to
# compensate negative offsets.
in_line = pmul(line.dir, 2 * abs(offset))
#in_line = line.dir.mul(2 * abs(offset))
return Line(psub(padd(line.p1, cross_offset), in_line), padd(padd(line.p2, cross_offset), in_line))
#return Line(line.p1.add(cross_offset).sub(in_line), line.p2.add(cross_offset).add(in_line))
def do_lines_intersection(l1, l2):
""" calculate the new intersection between two neighbouring lines
"""
......@@ -1078,15 +1049,11 @@ class Polygon(TransformableContainer):
return
x1, x2, x3, x4 = l2.p1, l2.p2, l1.p1, l1.p2
a = psub(x2, x1)
#a = x2.sub(x1)
b = psub(x4, x3)
#b = x4.sub(x3)
c = psub(x3, x1)
#c = x3.sub(x1)
# see http://mathworld.wolfram.com/Line-LineIntersection.html (24)
try:
factor = pdot(pcross(c, b), pcross(a, b)) / pnormsq(pcross(a, b))
#factor = c.cross(b).dot(a.cross(b)) / a.cross(b).normsq
except ZeroDivisionError:
l2.p1 = None
return
......@@ -1095,7 +1062,6 @@ class Polygon(TransformableContainer):
l2.p1 = None
else:
intersection = padd(x1, pmul(a, factor))
#intersection = x1.add(a.mul(factor))
if Line(l1.p1, intersection).dir != l1.dir:
# Remove lines that would change their direction due to the
# new intersection. These are usually lines that become
......@@ -1338,14 +1304,12 @@ class Polygon(TransformableContainer):
for index in range(len(collisions) - 1):
p1 = collisions[index][0]
p2 = collisions[index + 1][0]
#if p1.sub(p2).norm < epsilon:
if pnorm(psub(p1, p2)) < epsilon:
# ignore zero-length lines
continue
# Use the middle between p1 and p2 to check the
# inner/outer state.
p_middle = pdiv(padd(p1, p2), 2)
#p_middle = p1.add(p2).div(2)
p_inside = poly2.is_point_inside(p_middle) \
and not poly2.is_point_on_outline(p_middle)
if not p_inside:
......@@ -1391,7 +1355,6 @@ class Polygon(TransformableContainer):
p2, d2 = intersections[index + 1]
if p1 != p2:
middle = pdiv(padd(p1, p2), 2)
#middle = p1.add(p2).div(2)
new_line = Line(get_original_point(d1), get_original_point(d2))
if self.is_point_inside(middle):
inner.append(new_line)
......
pycam/Geometry/PolygonExtractor.py
View file @ 73cb1bb0
......@@ -80,7 +80,7 @@ class PolygonExtractor(object):
print "points=", path.points
i = 0
while i < len(path.points)-1:
if path.points[i].x > path.points[i+1].x:
if path.points[i][0] > path.points[i+1][0]:
if DEBUG_POLYGONEXTRACTOR2:
print "drop point %d:" % path.points[i].id
path.points = path.points[:i] + path.points[i+1:]
......@@ -95,7 +95,7 @@ class PolygonExtractor(object):
print "%d:" % path.id,
print "%d ->" % path.top_join.id
for point in path.points:
print "%d(%g,%g)" % (point.id, point.x, point.y),
print "%d(%g,%g)" % (point.id, point[0], point[1]),
print "->%d" % path.bot_join.id
path_list = []
......@@ -133,7 +133,7 @@ class PolygonExtractor(object):
for path in path_list:
print "path %d(w=%d): " % (path.id, path.winding),
for point in path.points:
print "%d(%g,%g)" % (point.id, point.x, point.y),
print "%d(%g,%g)" % (point.id, point[0], point[1]),
print
if self.current_dir == 0:
......@@ -156,13 +156,13 @@ class PolygonExtractor(object):
self.svg.fill("red")
else:
self.svg.fill("blue")
self.svg.AddDot(p.x, p.y)
self.svg.AddText(p.x, p.y, str(p.id))
self.svg.AddDot(p[0], p[1])
self.svg.AddText(p[0], p[1], str(p.id))
if prev:
self.svg.AddLine(p.x, p.y, prev.x, prev.y)
self.svg.AddLine(p[0], p[1], prev[0], prev[1])
prev = p
p = path.points[0]
self.svg.AddLine(p.x, p.y, prev.x, prev.y)
self.svg.AddLine(p[0], p[1], prev[0], prev[1])
self.svg.close()
self.cont.close()
......@@ -176,8 +176,8 @@ class PolygonExtractor(object):
def append(self, p):
if DEBUG_POLYGONEXTRACTOR3:
p.dir = self.current_dir
self.svg.AddDot(p.x, p.y)
self.svg.AddText(p.x, p.y, str(p.id))
self.svg.AddDot(p[0], p[1])
self.svg.AddText(p[0], p[1], str(p.id))
self.curr_line.append(p)
def end_scanline(self):
......@@ -187,7 +187,7 @@ class PolygonExtractor(object):
if self.policy == PolygonExtractor.CONTOUR and self.hor_path_list:
next_x = -INFINITE
if len(self.curr_line) > 0:
next_x = self.curr_line[0].x
next_x = self.curr_line[0][0]
self.delta_x = next_x - self.last_x
self.last_x = next_x
else:
......@@ -204,7 +204,7 @@ class PolygonExtractor(object):
inside = False
s = ""
for point in scanline:
next_x = point.x
next_x = point[0]
if inside:
s += "*" * int(next_x - last)
else:
......@@ -217,17 +217,17 @@ class PolygonExtractor(object):
print "active paths: ",
for path in self.curr_path_list:
print "%d(%g,%g)" \
% (path.id, path.points[-1].x, path.points[-1].y),
% (path.id, path.points[-1][0], path.points[-1][1]),
print
print "prev points: ",
for point in self.prev_line:
print "(%g,%g)" % (point.x, point.y),
print "(%g,%g)" % (point[0], point[1]),
print
print "active points: ",
for point in scanline:
print "%d(%g,%g)" % (point.id, point.x, point.y),
print "%d(%g,%g)" % (point.id, point[0], point.[1]),
print
prev_point = Iterator(self.prev_line)
......@@ -246,13 +246,13 @@ class PolygonExtractor(object):
p0 = Path()
p0.winding = winding + 1
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g)" % (p0.id, c0.x, c0.y)
print "new path %d(%g,%g)" % (p0.id, c0[0], c0[1])
p0.append(c0)
self.curr_path_list.append(p0)
p1 = Path()
p1.winding = winding
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g)" % (p1.id, c1.x, c1.y)
print "new path %d(%g,%g)" % (p1.id, c1[0], c1[1])
p1.append(c1)
self.curr_path_list.append(p1)
p0.top_join = p1
......@@ -282,16 +282,16 @@ class PolygonExtractor(object):
c1 = curr_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "overlap test: p0=%g p1=%g" % (p0.x, p1.x)
print "overlap test: c0=%g c1=%g" % (c0.x, c1.x)
print "overlap test: p0=%g p1=%g" % (p0[0], p1[0])
print "overlap test: c0=%g c1=%g" % (c0[0], c1[0])
if c1.x < p0.x:
if c1[0] < p0[0]:
# new segment is completely to the left
# new path starts
s0 = Path()
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g) w=%d" \
% (s0.id, c0.x, c0.y, winding + 1)
% (s0.id, c0[0], c0[0], winding + 1)
s0.append(c0)
curr_path.insert(s0)
s1 = Path()
......@@ -299,14 +299,14 @@ class PolygonExtractor(object):
s1.winding = winding
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g) w=%d" \
% (s1.id, c1.x, c1.y, winding)
% (s1.id, c1[0], c1[1], winding)
s1.append(c1)
curr_path.insert(s1)
curr_point.next()
curr_point.next()
s0.top_join = s1
s1.top_join = s0
elif c0.x > p1.x:
elif c0[0] > p1[0]:
# new segment is completely to the right
# old path ends
s0 = curr_path.takeNext()
......@@ -342,9 +342,9 @@ class PolygonExtractor(object):
p2 = prev_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "join test: p0=%g p1=%g p2=%g" \
% (p0.x, p1.x, p2.x)
print "join test: c0=%g c1=%g" % (c0.x, c1.x)
if p2.x <= c1.x:
% (p0[0], p1[0], p2[0])
print "join test: c0=%g c1=%g" % (c0[0], c1[0])
if p2[0] <= c1[0]:
overlap_p = True
if self.policy == PolygonExtractor.CONTOUR:
s0 = curr_path.takeNext()
......@@ -384,10 +384,10 @@ class PolygonExtractor(object):
if curr_point.remains()>=2:
c2 = curr_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "split test: p0=%g p1=%g" % (p0.x, p1.x)
print "split test: p0=%g p1=%g" % (p0[0], p1[0])
print "split test: c0=%g c1=%g c2=%g" \
% (c0.x, c1.x, c2.x)
if c2.x <= p1.x:
% (c0[0], c1[0], c2[0])
if c2[0] <= p1[0]:
overlap_c = True
s0 = Path()
s1 = Path()
......@@ -419,14 +419,14 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR:
print "add to path %d(%g,%g)" \
% (left_path.id, left_point.x, left_point.y)
% (left_path.id, left_point[0], left_point[1])
left_path.append(left_point)
right_path.append(right_point)
if right_path == curr_path.peek():
curr_path.next()
if DEBUG_POLYGONEXTRACTOR:
print "add to path %d(%g,%g)" \
% (right_path.id, right_point.x, right_point.y)
% (right_path.id, right_point[0], right_point[1])
winding = right_path.winding
prev_point.next()
curr_point.next()
......@@ -434,8 +434,8 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR:
print "active paths: ",
for path in self.curr_path_list:
print "%d(%g,%g,w=%d)" % (path.id, path.points[-1].x,
path.points[-1].y, path.winding),
print "%d(%g,%g,w=%d)" % (path.id, path.points[-1][0],
path.points[-1][1], path.winding),
print
self.prev_line = scanline
......@@ -448,11 +448,11 @@ class PolygonExtractor(object):
self.cont.fill("red")
else:
self.cont.fill("blue")
self.cont.AddDot(p.x, p.y)
self.cont.AddDot(p[0], p[1])
self.cont.fill("black")
self.cont.AddText(p.x, p.y, str(p.id))
self.cont.AddText(p[0], p[1], str(p.id))
if prev:
self.cont.AddLine(prev.x, prev.y, p.x, p.y)
self.cont.AddLine(prev[0], prev[1], p[0], p[1])
prev = p
if DEBUG_POLYGONEXTRACTOR:
......@@ -460,7 +460,7 @@ class PolygonExtractor(object):
inside = False
s = ""
for point in scanline:
next_y = point.y
next_y = point[1]
if inside:
s += "*" * int(next_y - last)
else:
......@@ -473,17 +473,17 @@ class PolygonExtractor(object):
print "active paths: ",
for path in self.curr_path_list:
print "%d(%g,%g)" \
% (path.id, path.points[-1].x, path.points[-1].y),
% (path.id, path.points[-1][0], path.points[-1][1]),
print
print "prev points: ",
for point in self.prev_line:
print "(%g,%g)" % (point.x, point.y),
print "(%g,%g)" % (point[0], point[1]),
print
print "active points: ",
for point in scanline:
print "%d(%g,%g)" % (point.id, point.x, point.y),
print "%d(%g,%g)" % (point.id, point[0], point[1]),
print
prev_point = Iterator(self.prev_line)
......@@ -502,13 +502,13 @@ class PolygonExtractor(object):
p0 = Path()
p0.winding = winding + 1
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g)" % (p0.id, c0.x, c0.y)
print "new path %d(%g,%g)" % (p0.id, c0[0], c0[1])
p0.append(c0)
self.curr_path_list.append(p0)
p1 = Path()
p1.winding = winding
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g)" % (p1.id, c1.x, c1.y)
print "new path %d(%g,%g)" % (p1.id, c1[0], c1[1])
p1.append(c1)
self.curr_path_list.append(p1)
p0.top_join = p1
......@@ -538,16 +538,16 @@ class PolygonExtractor(object):
c1 = curr_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "overlap test: p0=%g p1=%g" % (p0.x, p1.x)
print "overlap test: c0=%g c1=%g" % (c0.x, c1.x)
print "overlap test: p0=%g p1=%g" % (p0[0], p1[0])
print "overlap test: c0=%g c1=%g" % (c0[0], c1[0])
if c1.y < p0.y:
if c1[1] < p0[1]:
# new segment is completely to the left
# new path starts
s0 = Path()
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g) w=%d" \
% (s0.id, c0.x, c0.y, winding + 1)
% (s0.id, c0[0], c0[1], winding + 1)
s0.append(c0)
curr_path.insert(s0)
s1 = Path()
......@@ -555,14 +555,14 @@ class PolygonExtractor(object):
s1.winding = winding
if DEBUG_POLYGONEXTRACTOR:
print "new path %d(%g,%g) w=%d" \
% (s1.id, c1.x, c1.y, winding)
% (s1.id, c1[0], c1[1], winding)
s1.append(c1)
curr_path.insert(s1)
curr_point.next()
curr_point.next()
s0.top_join = s1
s1.top_join = s0
elif c0.y > p1.y:
elif c0[1] > p1[1]:
# new segment is completely to the right
# old path ends
s0 = curr_path.takeNext()
......@@ -598,9 +598,9 @@ class PolygonExtractor(object):
p2 = prev_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "join test: p0=%g p1=%g p2=%g" \
% (p0.x, p1.x, p2.x)
print "join test: c0=%g c1=%g" % (c0.x, c1.x)
if p2.y <= c1.y:
% (p0[0], p1[0], p2[0])
print "join test: c0=%g c1=%g" % (c0[0], c1[0])
if p2[1] <= c1[1]:
overlap_p = True
if self.policy == PolygonExtractor.CONTOUR:
s0 = curr_path.takeNext()
......@@ -640,10 +640,10 @@ class PolygonExtractor(object):
if curr_point.remains()>=2:
c2 = curr_point.peek(1)
if DEBUG_POLYGONEXTRACTOR:
print "split test: p0=%g p1=%g" % (p0.x, p1.x)
print "split test: p0=%g p1=%g" % (p0[0], p1[0])
print "split test: c0=%g c1=%g c2=%g" \
% (c0.x, c1.x, c2.x)
if c2.y <= p1.y:
% (c0[0], c1[0], c2[0])
if c2[1] <= p1[1]:
overlap_c = True
s0 = Path()
s1 = Path()
......@@ -675,14 +675,14 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR:
print "add to path %d(%g,%g)" \
% (left_path.id, left_point.x, left_point.y)
% (left_path.id, left_point[0], left_point[1])
left_path.append(left_point)
right_path.append(right_point)
if right_path == curr_path.peek():
curr_path.next()
if DEBUG_POLYGONEXTRACTOR:
print "add to path %d(%g,%g)" \
% (right_path.id, right_point.x, right_point.y)
% (right_path.id, right_point[0], right_point[1])
winding = right_path.winding
prev_point.next()
curr_point.next()
......@@ -690,8 +690,8 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR:
print "active paths: ",
for path in self.curr_path_list:
print "%d(%g,%g,w=%d)" % (path.id, path.points[-1].x,
path.points[-1].y, path.winding),
print "%d(%g,%g,w=%d)" % (path.id, path.points[-1][0],
path.points[-1][1], path.winding),
print
self.prev_line = scanline
......@@ -702,26 +702,26 @@ class PolygonExtractor(object):
hor_path_list = []
for s in self.hor_path_list:
allsame = True
miny = s.points[0].y
maxy = s.points[0].y
miny = s.points[0][1]
maxy = s.points[0][1]
for p in s.points:
if not p.x == s.points[0].x:
if not p[0] == s.points[0][0]:
allsame = False
if p.y < miny:
miny = p.y
if p.y > maxy:
maxy = p.y
if p[1] < miny:
miny = p[1]
if p[1] > maxy:
maxy = p[1]
if allsame:
if DEBUG_POLYGONEXTRACTOR2:
print "all same !"
s0 = Path()
for p in s.points:
if p.y == miny:
if p[1] == miny:
s0.append(p)
hor_path_list.append(s0)
s1 = Path()
for p in s.points:
if p.y == maxy:
if p[1] == maxy:
s1.append(p)
hor_path_list.append(s1)
continue
......@@ -729,28 +729,28 @@ class PolygonExtractor(object):
p_iter = CyclicIterator(s.points)
p = s.points[0]
next_p = p_iter.next()
while not ((prev.x >= p.x) and (next_p.x > p.x)):
while not ((prev[0] >= p[0]) and (next_p[0] > p[0])):
p = next_p
next_p = p_iter.next()
count = 0
while count < len(s.points):
s0 = Path()
while next_p.x >= p.x:
while next_p[0] >= p[0]:
s0.append(p)
p = next_p
next_p = p_iter.next()
count += 1
s0.append(p)
while (len(s0.points) > 1) \
and (s0.points[0].x == s0.points[1].x):
and (s0.points[0][0] == s0.points[1][0]):
s0.points = s0.points[1:]
while (len(s0.points) > 1) \
and (s0.points[-2].x == s0.points[-1].x):
and (s0.points[-2][0] == s0.points[-1][0]):
s0.points = s0.points[0:-1]
hor_path_list.append(s0)
s1 = Path()
while next_p.x <= p.x:
while next_p[0] <= p[0]:
s1.append(p)
p = next_p
next_p = p_iter.next()
......@@ -758,13 +758,13 @@ class PolygonExtractor(object):
s1.append(p)
s1.reverse()
while (len(s1.points) > 1) \
and (s1.points[0].x == s1.points[1].x):
and (s1.points[0][0] == s1.points[1][0]):
s1.points = s1.points[1:]
while (len(s1.points) > 1) \
and (s1.points[-2].x == s1.points[-1].x):
and (s1.points[-2][0] == s1.points[-1][0]):
s1.points = s1.points[:-1]
hor_path_list.append(s1)
hor_path_list.sort(cmp=lambda a, b: cmp(a.points[0].x, b.points[0].x))
hor_path_list.sort(cmp=lambda a, b: cmp(a.points[0][0], b.points[0][0]))
if DEBUG_POLYGONEXTRACTOR2:
print "ver_hor_path_list = ", hor_path_list
for s in hor_path_list:
......@@ -785,12 +785,12 @@ class PolygonExtractor(object):
next_x = INFINITE
if self.ver_hor_path_list \
and (self.ver_hor_path_list[0].points[0].x < next_x):
next_x = self.ver_hor_path_list[0].points[0].x
and (self.ver_hor_path_list[0].points[0][0] < next_x):
next_x = self.ver_hor_path_list[0].points[0][0]
if self.act_hor_path_list \
and (self.act_hor_path_list[0].points[0].x < next_x):
next_x = self.act_hor_path_list[0].points[0].x
and (self.act_hor_path_list[0].points[0][0] < next_x):
next_x = self.act_hor_path_list[0].points[0][0]
if next_x >= _next_x:
return
......@@ -801,13 +801,13 @@ class PolygonExtractor(object):
print "next_x =", next_x
if self.ver_hor_path_list \
and (self.ver_hor_path_list[0].points[0].x <= next_x):
and (self.ver_hor_path_list[0].points[0][0] <= next_x):
while self.ver_hor_path_list \
and (self.ver_hor_path_list[0].points[0].x <= next_x):
and (self.ver_hor_path_list[0].points[0][0] <= next_x):
self.act_hor_path_list.append(self.ver_hor_path_list[0])
self.ver_hor_path_list = self.ver_hor_path_list[1:]
self.act_hor_path_list.sort(cmp=lambda a, b:
cmp(a.points[0].x, b.points[0].x))
cmp(a.points[0][0], b.points[0][0]))
scanline = []
i = 0
......@@ -816,7 +816,7 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR2:
print "s =", s
scanline.append(s.points[0])
if s.points[0].x <= next_x:
if s.points[0][0] <= next_x:
if len(s.points) <= 1:
if DEBUG_POLYGONEXTRACTOR2:
print "remove list"
......@@ -830,17 +830,17 @@ class PolygonExtractor(object):
if DEBUG_POLYGONEXTRACTOR2:
print "remove point", s.points[0]
s.points = s.points[1:]
if len(s.points)> 0 and s.points[0].x == next_x:
if len(s.points)> 0 and s.points[0][0] == next_x:
# TODO: the variable "repeat" is never used.
# Any idea?
repeat = True
i += 1
self.act_hor_path_list.sort(cmp=lambda a, b:
cmp(a.points[0].x, b.points[0].x))
cmp(a.points[0][0], b.points[0][0]))
if len(scanline) == 0:
return
scanline.sort(cmp=lambda a, b: cmp(a.y, b.y))
scanline.sort(cmp=lambda a, b: cmp(a[1], b[1]))
if DEBUG_POLYGONEXTRACTOR2:
print "scanline' =", scanline
print "ver_hor_path_list =", self.ver_hor_path_list
......
pycam/Geometry/Triangle.py
View file @ 73cb1bb0
......@@ -65,26 +65,18 @@ class Triangle(IDGenerator, TransformableContainer):
# calculate normal, if p1-p2-pe are in clockwise order
if self.normal is None:
self.normal = pnormalized(pcross(psub(self.p3, self.p1), psub(self.p2, self.p1)))
#self.normal = self.p3.sub(self.p1).cross(self.p2.sub( \
# self.p1)).normalized()
if not len(self.normal) > 3:
self.normal = (self.normal[0], self.normal[1], self.normal[2], 'v')
self.center = pdiv(padd(padd(self.p1, self.p2), self.p3), 3)
# self.center = self.p1.add(self.p2).add(self.p3).div(3)
self.plane = Plane(self.center, self.normal)
# calculate circumcircle (resulting in radius and middle)
denom = pnorm(pcross(psub(self.p2, self.p1), psub(self.p3, self.p2)))
#denom = self.p2.sub(self.p1).cross(self.p3.sub(self.p2)).norm
self.radius = (pnorm(psub(self.p2, self.p1)) * pnorm(psub(self.p3, self.p2)) * pnorm(psub(self.p3, self.p1))) / (2 * denom)
#self.radius = (self.p2.sub(self.p1).norm * self.p3.sub(self.p2).norm * self.p3.sub(self.p1).norm) / (2 * denom)
self.radiussq = self.radius ** 2
denom2 = 2 * denom * denom
alpha = pnormsq(psub(self.p3, self.p2)) * pdot(psub(self.p1, self.p2), psub(self.p1, self.p3)) / denom2
#alpha = self.p3.sub(self.p2).normsq * self.p1.sub(self.p2).dot(self.p1.sub(self.p3)) / denom2
beta = pnormsq(psub(self.p1, self.p3)) * pdot(psub(self.p2, self.p1), psub(self.p2, self.p3)) / denom2
#beta = self.p1.sub(self.p3).normsq * self.p2.sub(self.p1).dot(self.p2.sub(self.p3)) / denom2
gamma = pnormsq(psub(self.p1, self.p2)) * pdot(psub(self.p3, self.p1), psub(self.p3, self.p2)) / denom2
#gamma = self.p1.sub(self.p2).normsq * self.p3.sub(self.p1).dot(self.p3.sub(self.p2)) / denom2
self.middle = (self.p1[0] * alpha + self.p2[0] * beta + self.p3[0] * gamma,
self.p1[1] * alpha + self.p2[1] * beta + self.p3[1] * gamma,
self.p1[2] * alpha + self.p2[2] * beta + self.p3[2] * gamma)
......@@ -145,17 +137,12 @@ class Triangle(IDGenerator, TransformableContainer):
c = self.center
GL.glTranslate(c[0], c[1], c[2])
p12 = pmul(padd(self.p1, self.p2), 0.5)
#p12 = self.p1.add(self.p2).mul(0.5)
p3_12 = pnormalized(psub(self.p3, p12))
#p3_12 = self.p3.sub(p12).normalized()
p2_1 = pnormalized(psub(self.p1, self.p2))
#p2_1 = self.p1.sub(self.p2).normalized()
pn = pcross(p2_1, p3_12)
#pn = p2_1.cross(p3_12)
GL.glMultMatrixf((p2_1[0], p2_1[1], p2_1[2], 0, p3_12[0], p3_12[1],
p3_12[2], 0, pn[0], pn[1], pn[2], 0, 0, 0, 0, 1))
n = pmul(self.normal, 0.01)
#n = self.normal.mul(0.01)
GL.glTranslatef(n[0], n[1], n[2])
maxdim = max((self.maxx - self.minx), (self.maxy - self.miny),
(self.maxz - self.minz))
......@@ -172,19 +159,13 @@ class Triangle(IDGenerator, TransformableContainer):
if False: # draw point id on triangle face
c = self.center
p12 = pmul(padd(self.p1, self.p2), 0.5)
#p12 = self.p1.add(self.p2).mul(0.5)
p3_12 = pnormalized(psub(self.p3, p12))
#p3_12 = self.p3.sub(p12).normalized()
p2_1 = pnormalized(psub(self.p1, self.p2))
#p2_1 = self.p1.sub(self.p2).normalized()
pn = pcross(p2_1, p3_12)
#pn = p2_1.cross(p3_12)
n = pmul(self.normal, 0.01)
#n = self.normal.mul(0.01)
for p in (self.p1, self.p2, self.p3):
GL.glPushMatrix()
pp = psub(p, pmul(psub(p, c), 0.3))
#pp = p.sub(p.sub(c).mul(0.3))
GL.glTranslate(pp[0], pp[1], pp[2])
GL.glMultMatrixf((p2_1[0], p2_1[1], p2_1[2], 0, p3_12[0], p3_12[1],
p3_12[2], 0, pn[0], pn[1], pn[2], 0, 0, 0, 0, 1))
......@@ -202,17 +183,14 @@ class Triangle(IDGenerator, TransformableContainer):
# http://www.blackpawn.com/texts/pointinpoly/default.html
# Compute vectors
v0 = psub(self.p3, self.p1)
#v0 = self.p3.sub(self.p1)
v1 = psub(self.p2, self.p1)
#v1 = self.p2.sub(self.p1)
v2 = psub(p, self.p1)
#v2 = p.sub(self.p1)
# Compute dot products
dot00 = pdot(v0, v0) # dot00 = v0.dot(v0)
dot01 = pdot(v0, v1) # dot01 = v0.dot(v1)
dot02 = pdot(v0, v2) # dot02 = v0.dot(v2)
dot11 = pdot(v1, v1) # dot11 = v1.dot(v1)
dot12 = pdot(v1, v2) # dot12 = v1.dot(v2)
dot00 = pdot(v0, v0)
dot01 = pdot(v0, v1)
dot02 = pdot(v0, v2)
dot11 = pdot(v1, v1)
dot12 = pdot(v1, v2)
# Compute barycentric coordinates
denom = dot00 * dot11 - dot01 * dot01
if denom == 0:
......@@ -232,11 +210,8 @@ class Triangle(IDGenerator, TransformableContainer):
sub.append(self)
else:
p4 = pdiv(padd(self.p1, self.p2), 2)
#p4 = self.p1.add(self.p2).div(2)
p5 = pdiv(padd(self.p2, self.p3), 2)
#p5 = self.p2.add(self.p3).div(2)
p6 = pdiv(padd(self.p3, self.p1), 2)
#p6 = self.p3.add(self.p1).div(2)
sub += Triangle(self.p1, p4, p6).subdivide(depth - 1)
sub += Triangle(p6, p5, self.p3).subdivide(depth - 1)
sub += Triangle(p6, p4, p5).subdivide(depth - 1)
......@@ -245,6 +220,5 @@ class Triangle(IDGenerator, TransformableContainer):
def get_area(self):
cross = pcross(psub(self.p2, self.p1), psub(self.p3, self.p1))
#cross = self.p2.sub(self.p1).cross(self.p3.sub(self.p1))
return pnorm(cross) / 2
pycam/Geometry/__init__.py
View file @ 73cb1bb0
......@@ -38,23 +38,16 @@ def get_bisector(p1, p2, p3, up_vector):
of the angle.
"""
d1 = pnormalized(psub(p2, p1))
#d1 = p2.sub(p1).normalized()
d2 = pnormalized(psub(p2, p3))
#d2 = p2.sub(p3).normalized()
bisector_dir = pnormalized(padd(d1, d2))
#bisector_dir = d1.add(d2).normalized()
if bisector_dir is None:
# the two vectors pointed to opposite directions
bisector_dir = pnormalized(pcross(d1, up_vector))
#bisector_dir = d1.cross(up_vector).normalized()
else:
skel_up_vector = pcross(bisector_dir, psub(p2, p1))
#skel_up_vector = bisector_dir.cross(p2.sub(p1))
#if up_vector.dot(skel_up_vector) < 0:
if pdot(up_vector, skel_up_vector) < 0:
# reverse the skeleton vector to point outwards
bisector_dir = pmul(bisector_dir, -1)
#bisector_dir = bisector_dir.mul(-1)
return bisector_dir
def get_angle_pi(p1, p2, p3, up_vector, pi_factor=False):
......@@ -69,13 +62,10 @@ def get_angle_pi(p1, p2, p3, up_vector, pi_factor=False):
The result is in a range between 0 and 2*PI.
"""
d1 = pnormalized(psub(p2, p1))
#d1 = p2.sub(p1).normalized()
d2 = pnormalized(psub(p2, p3))
#d2 = p2.sub(p3).normalized()
if (d1 is None) or (d2 is None):
return 2 * math.pi
angle = math.acos(pdot(d1, d2))
#angle = math.acos(d1.dot(d2))
# check the direction of the points (clockwise/anti)
# The code is taken from Polygon.get_area
value = [0, 0, 0]
......@@ -145,7 +135,6 @@ def get_bezier_lines(points_with_bulge, segments=32):
# straight line
return [Line.Line(p1, p2)]
straight_dir = pnormalized(psub(p2, p1))
#straight_dir = p2.sub(p1).normalized()
#bulge1 = max(-1.0, min(1.0, bulge1))
bulge1 = math.atan(bulge1)
rot_matrix = Matrix.get_rotation_matrix_axis_angle((0, 0, 1),
......@@ -170,7 +159,6 @@ def get_bezier_lines(points_with_bulge, segments=32):
# and a bulge of 1 is a semicircle.
alpha = 2 * (abs(bulge1) + abs(bulge2))
dist = pnorm(psub(p2, p1))
#dist = p2.sub(p1).norm
# calculate the radius of the circumcircle - avoiding divide-by-zero
if (abs(alpha) < epsilon) or (abs(math.pi - alpha) < epsilon):
radius = dist / 2.0
......@@ -181,20 +169,11 @@ def get_bezier_lines(points_with_bulge, segments=32):
# seems to work well.
factor = 4 * radius * math.tan(alpha / 4.0)
dir1 = pmul(dir1, factor)
#dir1 = dir1.mul(factor)
dir2 = pmul(dir2, factor)
#dir2 = dir2.mul(factor)
for index in range(segments + 1):
# t: 0..1
t = float(index) / segments
# see: http://en.wikipedia.org/wiki/Cubic_Hermite_spline
#p = p1.mul(2 * t ** 3 - 3 * t ** 2 + 1).add(
# dir1.mul(t ** 3 - 2 * t ** 2 + t).add(
# p2.mul(-2 * t ** 3 + 3 * t ** 2).add(
# dir2.mul(t ** 3 - t ** 2)
# )
# )
#)
p = padd( pmul(p1, 2 * t ** 3 - 3 * t ** 2 + 1) ,padd( pmul(dir1, t ** 3 - 2 * t ** 2 + t), padd(pmul(p2, -2 * t ** 3 + 3 * t ** 2) ,pmul(dir2, t ** 3 - t ** 2))))
result_points.append(p)
# create lines
......
pycam/Geometry/intersection.py
View file @ 73cb1bb0
......@@ -61,16 +61,13 @@ def intersect_lines(xl, zl, nxl, nzl, xm, zm, nxm, nzm):
def intersect_cylinder_point(center, axis, radius, radiussq, direction, point):
# take a plane along direction and axis
n = pnormalized(pcross(direction, axis))
#n = direction.cross(axis).normalized()
# distance of the point to this plane
d = pdot(n, point) - pdot(n, center)
#d = n.dot(point) - n.dot(center)
if abs(d) > radius - epsilon:
return (None, None, INFINITE)
# ccl is on cylinder
d2 = sqrt(radiussq-d*d)
ccl = padd( padd(center, pmul(n, d)), pmul(direction, d2))
#ccl = center.add(n.mul(d)).add(direction.mul(d2))
# take plane through ccl and axis
plane = Plane(ccl, direction)
# intersect point with plane
......@@ -81,7 +78,6 @@ def intersect_cylinder_line(center, axis, radius, radiussq, direction, edge):
d = edge.dir
# take a plane throught the line and along the cylinder axis (1)
n = pcross(d, axis)
#n = d.cross(axis)
if pnorm(n) == 0:
# no contact point, but should check here if cylinder *always*
# intersects line...
......@@ -90,22 +86,17 @@ def intersect_cylinder_line(center, axis, radius, radiussq, direction, edge):
# the contact line between the cylinder and this plane (1)
# is where the surface normal is perpendicular to the plane
# so line := ccl + \lambda * axis
#if n.dot(direction) < 0:
if pdot(n, direction) < 0:
ccl = psub(center, pmul(n, radius))
#ccl = center.sub(n.mul(radius))
else:
ccl = padd(center, pmul(n, radius))
#ccl = center.add(n.mul(radius))
# now extrude the contact line along the direction, this is a plane (2)
n2 = pcross(direction, axis)
#n2 = direction.cross(axis)
if pnorm(n2) == 0:
# no contact point, but should check here if cylinder *always*
# intersects line...
return (None, None, INFINITE)
n2 = pnormalized(n2)
#n2 = n2.normalized()
plane1 = Plane(ccl, n2)
# intersect this plane with the line, this gives us the contact point
(cp, l) = plane1.intersect_point(d, edge.p1)
......@@ -118,7 +109,6 @@ def intersect_cylinder_line(center, axis, radius, radiussq, direction, edge):
# gives us the cutter contact point
(ccp, l) = plane2.intersect_point(direction, cp)
cp = padd(ccp, pmul(direction, -l))
#cp = ccp.add(direction.mul(-l))
return (ccp, cp, -l)
def intersect_circle_plane(center, radius, direction, triangle):
......@@ -131,13 +121,10 @@ def intersect_circle_plane(center, radius, direction, triangle):
if pnorm(n2) == 0:
(cp, d) = triangle.plane.intersect_point(direction, center)
ccp = psub(cp, pmul(direction, d))
#ccp = cp.sub(direction.mul(d))
return (ccp, cp, d)
n2 = pnormalized(n2)
#n2 = n2.normalized()
# the cutter contact point is on the circle, where the surface normal is n
ccp = padd(center, pmul(n2, -radius))
#ccp = center.add(n2.mul(-radius))
# intersect the plane with a line through the contact point
(cp, d) = triangle.plane.intersect_point(direction, ccp)
return (ccp, cp, d)
......@@ -148,7 +135,6 @@ def intersect_circle_point(center, axis, radius, radiussq, direction, point):
# intersect with line gives ccp
(ccp, l) = plane.intersect_point(direction, point)
# check if inside circle
#if ccp and (center.sub(ccp).normsq < radiussq - epsilon):
if ccp and (pnormsq(psub(center, ccp)) < radiussq - epsilon):
return (ccp, point, -l)
return (None, None, INFINITE)
......@@ -164,38 +150,30 @@ def intersect_circle_line(center, axis, radius, radiussq, direction, edge):
(p2, l) = plane.intersect_point(direction, edge.p2)
pc = Line(p1, p2).closest_point(center)
d_sq = pnormsq(psub(pc, center))
#d_sq = pc.sub(center).normsq
if d_sq >= radiussq:
return (None, None, INFINITE)
a = sqrt(radiussq - d_sq)
d1 = pdot(psub(p1, pc), d)
#d1 = p1.sub(pc).dot(d)
d2 = pdot(psub(p2, pc), d)
#d2 = p2.sub(pc).dot(d)
ccp = None
cp = None
if abs(d1) < a - epsilon:
ccp = p1
cp = psub(p1, pmul(direction, l))
#cp = p1.sub(direction.mul(l))
elif abs(d2) < a - epsilon:
ccp = p2
cp = psub(p2, pmul(direction, l))
#cp = p2.sub(direction.mul(l))
elif ((d1 < -a + epsilon) and (d2 > a - epsilon)) \
or ((d2 < -a + epsilon) and (d1 > a - epsilon)):
ccp = pc
cp = psub(pc, pmul(direction, l))
#cp = pc.sub(direction.mul(l))
return (ccp, cp, -l)
n = pcross(d, direction)
#n = d.cross(direction)
if pnorm(n)== 0:
# no contact point, but should check here if circle *always* intersects
# line...
return (None, None, INFINITE)
n = pnormalized(n)
#n = n.normalized()
# take a plane through the base
plane = Plane(center, axis)
# intersect base with line
......@@ -204,21 +182,16 @@ def intersect_circle_line(center, axis, radius, radiussq, direction, edge):
return (None, None, INFINITE)
# intersection of 2 planes: lp + \lambda v
v = pcross(axis, n)
#v = axis.cross(n)
if pnorm(v) == 0:
return (None, None, INFINITE)
v = pnormalized(v)
#v = v.normalized()
# take plane through intersection line and parallel to axis
n2 = pcross(v, axis)
#n2 = v.cross(axis)
if pnorm(n2) == 0:
return (None, None, INFINITE)
n2 = pnormalized(n2)
#n2 = n2.normalized()
# distance from center to this plane
dist = pdot(n2, center) - pdot(n2, lp)
#dist = n2.dot(center) - n2.dot(lp)
distsq = dist * dist
if distsq > radiussq - epsilon:
return (None, None, INFINITE)
......@@ -227,9 +200,7 @@ def intersect_circle_line(center, axis, radius, radiussq, direction, edge):
if pdot(d, axis) < 0:
dist2 = -dist2
ccp = psub(center, psub(pmul(n2, dist), pmul(v, dist2)))
#ccp = center.sub(n2.mul(dist)).sub(v.mul(dist2))
plane = Plane(edge.p1, pcross(pcross(d, direction), d))
#plane = Plane(edge.p1, d.cross(direction).cross(d))
(cp, l) = plane.intersect_point(direction, ccp)
return (ccp, cp, l)
......@@ -241,10 +212,8 @@ def intersect_sphere_plane(center, radius, direction, triangle):
# the cutter contact point is on the sphere, where the surface normal is n
if pdot(n, direction) < 0:
ccp = psub(center, pmul(n, radius))
#ccp = center.sub(n.mul(radius))
else:
ccp = padd(center, pmul(n, radius))
#ccp = center.add(n.mul(radius))
# intersect the plane with a line through the contact point
(cp, d) = triangle.plane.intersect_point(direction, ccp)
return (ccp, cp, d)
......@@ -256,9 +225,7 @@ def intersect_sphere_point(center, radius, radiussq, direction, point):
# (2) (x-x_0)^2 = R^2
# (1) in (2) gives a quadratic in \lambda
p0_x0 = psub(center, point)
#p0_x0 = center.sub(point)
a = pnormsq(direction)
#a = direction.normsq
b = 2 * pdot(p0_x0, direction)
c = pnormsq(p0_x0) - radiussq
d = b * b - 4 * a * c
......@@ -270,24 +237,20 @@ def intersect_sphere_point(center, radius, radiussq, direction, point):
l = (-b - sqrt(d)) / (2 * a)
# cutter contact point
ccp = padd(point, pmul(direction, -l))
#ccp = point.add(direction.mul(-l))
return (ccp, point, l)
def intersect_sphere_line(center, radius, radiussq, direction, edge):
# make a plane by sliding the line along the direction (1)
d = edge.dir
n = pcross(n, direction)
#n = d.cross(direction)
if pnorm(n) == 0:
# no contact point, but should check here if sphere *always* intersects
# line...
return (None, None, INFINITE)
n = pnormalized(n)
#n = n.normalized()
# calculate the distance from the sphere center to the plane
dist = - pdot(center, n) + pdot(edge.p1, n)
#dist = - center.dot(n) + edge.p1.dot(n)
if abs(dist) > radius - epsilon:
return (None, None, INFINITE)
# this gives us the intersection circle on the sphere
......@@ -297,14 +260,12 @@ def intersect_sphere_line(center, radius, radiussq, direction, edge):
# which means the other component is perpendicular to this plane (2)
n2 = pnormalized(pcross(n, d))
#n2 = n.cross(d).normalized()
# the contact point is on a big circle through the sphere...
dist2 = sqrt(radiussq - dist * dist)
# ... and it's on the plane (1)
ccp = padd(center, padd(pmul(n, dist), pmul(n2, dist2)))
#ccp = center.add(n.mul(dist)).add(n2.mul(dist2))
# now intersect a line through this point with the plane (2)
plane = Plane(edge.p1, n2)
......@@ -321,17 +282,13 @@ def intersect_torus_plane(center, axis, majorradius, minorradius, direction,
return (None, None, INFINITE)
# find place on torus where surface normal is n
b = pmul(n, -1)
#b = n.mul(-1)
z = axis
a = psub(b, pmul(z,pdot(z, b)))
#a = b.sub(z.mul(z.dot(b)))
a_sq = pnormsq(a)
if a_sq <= 0:
return (None, None, INFINITE)
a = pdiv(a, sqrt(a_sq))
#a = a.div(sqrt(a_sq))
ccp = padd(padd(center, pmul(a, majorradius)), pmul(b, minorradius))
#ccp = center.add(a.mul(majorradius)).add(b.mul(minorradius))
# find intersection with plane
(cp, l) = triangle.plane.intersect_point(direction, ccp)
return (ccp, cp, l)
......@@ -360,38 +317,24 @@ def intersect_torus_point(center, axis, majorradius, minorradius, majorradiussq,
return (None, None, INFINITE)
l = majorradius + sqrt(minorradiussq - z * z)
n = pcross(axis, direction)
#n = axis.cross(direction)
d = pdot(n, point) - pdot(n, center)
#d = n.dot(point) - n.dot(center)
if abs(d) > l - epsilon:
return (None, None, INFINITE)
a = sqrt(l * l - d * d)
ccp = padd(padd(center, pmul(n, d)), pmul(direction, a))
#ccp = center.add(n.mul(d).add(direction.mul(a)))
ccp = (ccp[0], ccp[1], point[2])
#ccp.z = point.z
dist = pdot(psub(point, ccp), direction)
#dist = point.sub(ccp).dot(direction)
else:
# general case
x = psub(point, center)
#x = point.sub(center)
v = pmul(direction, -1)
#v = direction.mul(-1)
x_x = pdot(x, x)
#x_x = x.dot(x)
x_v = pdot(x, v)
#x_v = x.dot(v)
x1 = (x[0], x[1], 0)
#x1 = Point(x.x, x.y, 0)
v1 = (v[0], v[1], 0)
#v1 = Point(v.x, v.y, 0)
x1_x1 = pdot(x1, x1)
#x1_x1 = x1.dot(x1)
x1_v1 = pdot(x1, v1)
#x1_v1 = x1.dot(v1)
v1_v1 = pdot(v1, v1)
#v1_v1 = v1.dot(v1)
R2 = majorradiussq
r2 = minorradiussq
a = 1.0
......@@ -405,7 +348,6 @@ def intersect_torus_point(center, axis, majorradius, minorradius, majorradiussq,
else:
l = min(r)
ccp = padd(point, pmul(direction, -l))
#ccp = point.add(direction.mul(-l))
dist = l
return (ccp, point, dist)
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment