#!/usr/bin/env python3 from sympy import * from gcode2contour import Position, contour from sympy.plotting import plot3d from mpl_toolkits.mplot3d import Axes3D import matplotlib.colors as mcolors from copy import deepcopy, copy import matplotlib.pyplot as plt import numpy as np def set_aspect_equal_3d(ax): """Fix equal aspect bug for 3D plots.""" xlim = ax.get_xlim3d() ylim = ax.get_ylim3d() zlim = ax.get_zlim3d() from numpy import mean xmean = mean(xlim) ymean = mean(ylim) zmean = mean(zlim) plot_radius = max([abs(lim - mean_) for lims, mean_ in ((xlim, xmean), (ylim, ymean), (zlim, zmean)) for lim in lims]) ax.set_xlim3d([xmean - plot_radius, xmean + plot_radius]) ax.set_ylim3d([ymean - plot_radius, ymean + plot_radius]) ax.set_zlim3d([zmean - plot_radius, zmean + plot_radius]) def plot_contours(*args): """ Plots a list of contours Each input arguement is a list of contours All contours within each list will be the same color Contours in different lists will be different colors """ fig = plt.figure() ax = Axes3D(fig) # colors = [k for k in mcolors.cnames] colors = ['blue', 'red', 'green'] for i, contours in enumerate(args): for c in contours: xs = [pos[0] for pos in c.pos] ys = [pos[1] for pos in c.pos] zs = [pos[2] for pos in c.pos] ax.plot(xs, ys, zs, color=colors[i]) set_aspect_equal_3d(ax) plt.show() return class solver: """ Handles symbolic variables to solve for layer positions in various planes """ def __init__(self, cl, cx, cy, dz): self.x, self.y, self.z, self.cz = symbols('x y z cz') self.layer = cl * sin(cx*self.x)*sin(cy*self.y) + self.cz self.dz = dz self.def_prism() def get_z(self, x, y, layer=0): return float(self.layer.subs([(self.x, x), (self.y, y), (self.cz, layer*self.dz)])) def def_prism(self, x_min = -0.05, x_max = 0.05, y_min = -0.05, y_max = 0.05, z_min = 0., z_max = 0.1): """ save the prism size """ self.range = { self.x: (x_min, x_max), self.y: (y_min, y_max), self.z: (z_min, z_max), } def plane_intersection(self, sym, val, layer = 0): """ Takes a plane (not any plane, only x,y or z=val) and interesects it with the nth layer Returns symbolic expression for the intersection. Need to sample """ return self.layer.subs([(sym, val), (self.cz, layer*self.dz)]) def sample(self, expression, sym_in, res = 0.001): """ sample across one variable, get values of second variable in an expression """ v1 = np.arange(self.range[sym_in][0], self.range[sym_in][1]+ res, res) zs = [] for v in v1: zs.append(float(expression.subs(sym_in, v))) return v1, zs def trim_contour(self, c, zlims): """ take a contour and split it at z limits returns a list of contours. This list will be the original contour if it doesnt reach out of the limits If it doesn, regions outside the limits will be cut off, and regions inside will be split """ c_out = [] j = 0 # index of last intersection with limit state = zlims[0] < c.pos[0][2] < zlims[1] for i, pos in enumerate(c): if (zlims[0] < pos[2] < zlims[1]) is not state: # State change has happened if state: # Leaving zlim, save contour so far c_out.append(contour(c.pos[j:i],0)) else: # Entering zlim, set start point j = i state = not state #save new state if state: # the last stretch was in range c_out.append(contour(c.pos[j:i],0)) return c_out def contour_n(self, n): """ Get 4 contours for the nth layer 4 sides of the prism unlinked """ contours = [] expr1 = self.plane_intersection(self.x, self.range[self.x][0], layer = n) ys, zs = self.sample(expr1, self.y) contours += self.trim_contour( contour([Position(self.range[self.x][0],ys[i],zs[i]) for i in range(len(ys))], 0), self.range[self.z] ) expr2 = self.plane_intersection(self.x, self.range[self.x][1], layer = n) ys, zs = self.sample(expr2, self.y) contours += self.trim_contour( contour([Position(self.range[self.x][1],ys[i],zs[i]) for i in range(len(ys))], 0), self.range[self.z] ) expr3 = self.plane_intersection(self.y, self.range[self.y][0], layer = n) xs, zs = self.sample(expr3, self.x) contours += self.trim_contour( contour([Position(xs[i],self.range[self.y][0],zs[i]) for i in range(len(ys))], 0), self.range[self.z] ) expr4 = self.plane_intersection(self.y, self.range[self.y][1], layer = n) xs, zs = self.sample(expr4, self.x) contours += self.trim_contour( contour([Position(xs[i],self.range[self.y][1],zs[i]) for i in range(len(ys))], 0), self.range[self.z] ) # TODO Join adjacent contours into one return contours def infill_contour(self, n, density): """ return one contour for the infill if n is even, the infill for that layer is parallel to the x-axis """ # Get the direction to print in if n%2 == 0: main = self.x other = self.y else: main = self.y other = self.x # get the list of rows to do main_vals = np.arange(self.range[main][0],self.range[main][1], density)[1:] # do the first row expr1 = self.plane_intersection(main, main_vals[0], layer = n) os, zs = self.sample(expr1, other) poses = [Position(self.range[main][0],os[i],zs[i]) for i in range(len(os))] # loop through th remaining rows, linking them for row in main_vals[1:]: # TODO link to the next row # TODO do the next row (need to reverse directions each time) expr1 = self.plane_intersection(main, row, layer = n) os, zs = self.sample(expr1, other) poses = [Position(self.range[main][0],os[i],zs[i]) for i in range(len(os))] def show(self): """ Show the surface of the layer in the range of the cube at layer 0 """ plot3d(self.layer.subs(self.cz, 0), (self.x, self.range[self.x][0], self.range[self.x][1]), (self.y, self.range[self.y][0], self.range[self.y][1])) if __name__ == '__main__': s = solver(0.02, 100, 100, 0.003) contours = [] for i in np.arange(-15,55): contours += s.contour_n(i) contours += s.infill_contour(i, 0.02) plot_contours(contours)