diff --git a/h2/MBC_carlsen.py b/h2/MBC_carlsen.py new file mode 100644 index 0000000..402b27e --- /dev/null +++ b/h2/MBC_carlsen.py @@ -0,0 +1,193 @@ + +import matplotlib.pyplot as plt +import collections +from enum import Enum, auto +from random import randint + +Point = collections.namedtuple("Point", "x y") + +class Side(Enum): + ON = auto() + ABOVE = auto() + BELOW = auto() + + +def sidedness(slope, intersection, p3, eps=0.0000001): + print(slope * p3[0] + intersection ) + """ finds where a point is in regards to a line """ + if p3[1] - eps <= slope * p3[0] + intersection <= p3[1] + eps or p3[0] - eps <= (p3[1] - intersection)/slope <= p3[0] + eps: + return Side.ON + elif p3[1] > slope * p3[0] + intersection: + return Side.ABOVE + return Side.BELOW + +def diplay_prune_points(points, p1, p2): + xs = [p[0] for p in points] + ys = [p[1] for p in points] + + plt.plot(xs, ys, 'ro') + plt.plot([p1[0], p2[0]], [p1[1], p2[1]]) + plt.show() + + +def solve1D(points, xm, iteration_num): + #print("iter:", iteration_num) + point = points[iteration_num] + + if iteration_num == 0: + return -float('Inf'), -float('Inf') + + # lad point[1] = point[0] * a + b <=> y = x * a + b + # isolere b og sæt ind i constraints + a = None + b = point[1] - point[0] # * a + + # minimere xm*a + b, hvor b har ny værdi + # Vi regner kun med koefficienterne + # obj_fun = (xm - point[0]) + points[1] = (xm*a - xi*a) + y + + # max eller min + #print("XM og point[0]:", xm, point[0]) + a = xm - point[0] + #print("a lige her:", a) + + a_constraint_list = [] + # looping over the i constraints + for p in points[:iteration_num]: + # we can't make a straigt vertical line + if p[0] == point[0]: + if p[1] > point[1]: + p[0] = p[0] + 0.0001 + else: + point[0] = point[0] + 0.0001 + + print(p, point) + # Spring over den i'te constraint + if p != point: + # y_j - yi + c = p[1] - point[1] + # x_j * a - xi * a + a_diff = p[0] - point[0] + + print(c, a_diff, c/a_diff) + # det her er forkert og skal fikses + if a >= 0 and a_diff < 0: + a_constraint_list.append(-float('Inf')) + else: + a_constraint_list.append(c/a_diff) + + # hvis a > 0 så min + # hvis a < 0 så max. + if a >= 0: + v1 = max(a_constraint_list) + v2 = point[1] - point[0] * v1 + v = (v1, v2) + elif a < 0: + v1 = min(a_constraint_list) + v2 = point[1] - point[0] * v1 + v = (v1, v2) + + return v + + + +def findBridge(points, xm): + if xm > 0: + v = (-float('Inf'), -float('Inf')) + elif xm < 0: + 4 # noget + else: + # TODO: xm == 0 + 4 + # looping over constraints + for point in points: + # checking for violation + if not point[1] <= point[0] * v[0] + v[1]: + v = solve1D(points, xm, points.index(point)) + #print("HER OVER: ", v) + + slope = v[0] + intercept = v[1] + + line_points = [p for p in points if sidedness(slope, intercept, p) == Side.ON] + + if len(line_points) > 3: + print("Halli HAllo") + print(line_points) + + return line_points[0], line_points[1] + + +def find_median(points): + if len(points) % 2 == 0: + first_med_idx = int(len(points) / 2 - 1) + second_med_idx = int(len(points) / 2) + return (points[first_med_idx][0] + points[second_med_idx][0]) / 2 + else: + idx = int((len(points)-1) / 2) + return points[idx][0] + + +def upperHull(points, all_points): + print("Punkter:", points) + if len(points) < 2: + return [] + xm = find_median(points) + + # end-points of bridge + (xi, yi), (xj, yj) = findBridge(points, xm) + + #print(xm) + print((xi, yi), (xj, yj)) + prune_points = [p for p in points if p[0] < xi or xj < p[0]] + [(xi, yi), (xj, yj)] + + # Neden for er mest for visualisering + prune_all_points = [p for p in all_points if p[0] < xi or xj < p[0]] + [(xi, yi), (xj, yj)] + diplay_prune_points(prune_all_points, (xi, yi), (xj, yj)) + + Pl = [p for p in prune_points if p[0] < xm] + Pr = [p for p in prune_points if p[0] >= xm] + #print("Pl:", Pl, "Pr", Pr) + + print("\n") + # recurse results and return + ret = [(xi, yi), (xj, yj)] + ret = ret + [p for p in upperHull(Pl, prune_all_points) if p not in ret] + ret = ret + [p for p in upperHull(Pr, prune_all_points) if p not in ret] + return ret + + + +p1 = (2, 1) +p2 = (3, 4) +p3 = (5, 2) +p4 = (6, 5) + +list_of_points = [] + +list_of_points.append(p1) +list_of_points.append(p2) +list_of_points.append(p3) +list_of_points.append(p4) + + + + +xs = [p[0] for p in list_of_points] +ys = [p[1] for p in list_of_points] + +plt.plot(xs, ys, 'ro') +plt.show() + +result = list(sorted(upperHull(list_of_points, list_of_points))) + +result +res_xs = [p[0] for p in result] +res_ys = [p[1] for p in result] + + +#print("result", result) +plt.plot(res_xs, res_ys) +plt.plot(xs, ys, 'ro') +plt.show() +