2.1.1 Make CAD model

The part we have chosen to print is “3 sides of a cube”. This is a very simple 3D geometry that can easily have errors classified through an image classification neural network. Only 3 sides for speed, less material, and stackability.

OpenSCAD to create the model.

From OpenSCAD, the object was exported as an .STL file, which was then sliced into gcode by CURA Lulzbot Edition for the Lulzbot Mini.

2.1.2 Material

Lulzbot
Natural ABS 3mm Filament 1 kg reel
SKU: (3) 1775201312 (3)
Aleph Objects Incorporated

Reel 1: Prints 1-26, 167, 3b, 7b-13b, 314(Need to drop this print see comments section at end!)
Reel 2: Prints 27-112
Reel 3: Prints 113-127, 146, 161, 173, 186, 201, 202, 211, 239, 241, 259, 270, 272, 274

2.1.3 Choosing Various Printing Parameters

The parameters of interest are the Nozzle Temperature, Bed Temperature, Print Speed, Extrusion Multiplier and Fan Speed.

• The ideal bed temperature is suggested to be 110 C. Therefore, we are testing 100-120 C.
  
• Print speed was chosen to vary between 20 ad 80 mm/sec based on print speeds suggested here.
  
• Rao 2015 suggested the extrusion multiplier does best around 100%, thus we’ve chosen to vary between 90-110%.
  
• The Fan Speed multiplier was chosen to be either 30, 45, or 80%.

• Nozzle temperature ranges were chosen to be 225 C to 245 C. In Rao 2015, it was recommended not to go below 225 to prevent clogging. We have also noticed material burning above 245.

• Note: Prints 3 and 7-13 were done with a fan speed of 0% and were included in the dataset.

• Note: fan speeds are not possible between 0 to ~27%.

The R code for creating the parameter combinations and randomly ordering them:

nozzle.temp <- seq(225,245,10)                 
bed.temp <- seq(100,120,10)                    
print.speed <- seq(20,80,20)                  
extrusion.multiplier <- seq(90,110,10)
fan.speed.percent <- c(30,45,80)

# Combine them into the possible combinations
parameter.list <- list(Nozzle.Temp=nozzle.temp, Bed.Temp=bed.temp, 
                       Print.Speed=print.speed, Extrusion.Multiplier =extrusion.multiplier, 
                       Fan.Speed.Percent=fan.speed.percent)
all.options <- expand.grid(parameter.list)

# randomly order the rows
set.seed(42)
randomized.all.options <- all.options[sample(nrow(all.options),nrow(all.options)),]
# reset row numbers
rownames(randomized.all.options) <- NULL
# save these combinations as a .csv
#write.table("varied_parameters_print_order.txt", quote=F, x = randomized.all.options, col.names=FALSE)

There are 324 possible combinations:

We chose to randomly order the prints. If you look at the 3D representation of a problem surface, we could imagine the red to indicate an error while the dark blue indicates the highest fidelity. Imagine we do not randomly select the printing parameters and our problem surface is built in a straight line up the hill leading to the worst print conditions. By randomly sampling, we are hoping to approximate the surface of the problem prior to completing all 324 possible prints.

3D representation of a problem surface.

• Note: a thermocouple was used to determine the actual nozzle temperature in comparison to the pronterface read-out.
  

  	Pronterface.Temperature.Setting	Thermocouple.Reading	Difference
  1	260	257	3.0
  2	255	252	3.0
  3	250	247	3.0
  4	245	243	2.0
  5	240	238	2.0
  6	235	234	1.0
  7	230	229	1.0
  8	225	223	2.0
  9	220	219	1.0
  10	215	214	1.0
  11	210	209	1.0
  12	205	204	1.0
  13	200	200	0.0
  14	195	195.1	-0.1
  15	190	190.2	-0.2
  16	185	185.3	-0.3
  17	180	180.5	-0.5
  18	175	175.3	-0.3

These readings indicate we can rely on the pronterface reading of nozzle temperature.

2.1.4 Automatically Writing G-Code

As seen above, there are 324 possible combinations, requiring 324 gcodes. A python script was written to parse a base.gcode into the 324 required gcodes.

• Note: The base gcode needs to have a 100% extrusion multiplier.

with open("/path/to/parameter combinations") as print_order_file:
    for line in print_order_file:
        parameter_list = line.split()
        print_number = parameter_list[0]
        nozzle_temperature = parameter_list[1]
        bed_temperature = parameter_list[2]
        print_speed = parameter_list[3]
        extrusion_multiplier = parameter_list[4]
        fan_speed_percent = parameter_list[5]
        
        with open("/path/to/base gcode file") as base_gcode:
            # load in the base gcode
            gcode_data = base_gcode.readlines()
            # Edit the gcode based on the new parameters
            # Set Nozzle Temp
            gcode_data[62] = "M109 R" + nozzle_temperature + " ; wait for extruder to reach printing temp\n"
            gcode_data[570] = "M104 S" + nozzle_temperature + "\n"
            # Set Bed temp
            gcode_data[15] = "M140 S" + bed_temperature + " ; start bed heating up\n"
            gcode_data[63] = "M190 S" + bed_temperature + " ; wait for bed to reach printing temp\n"
            gcode_data[576] = "M140 S" + bed_temperature + "\n"
            # Set Print Speed and Extrusion Multiplier
            make_changes = False
            print_finished = False
            line_counter = 0
            previous_extrusion_multiplier = float(gcode_data[67][-4:-1])
            for code_line in gcode_data:
                # make sure outside of skirt
                if code_line == ";TYPE:WALL-INNER\n":
                    make_changes = True
                # make sure not at end
                if "M104 S0" in code_line:
                    make_changes = False
                # begin changing the line if it's not the skirt or after the print
                if make_changes:
                    # if the line is an extrusion and movement command
                    if all(x in code_line for x in ['G1', 'X', 'Y', 'E']):
                        split_line = code_line.split()
                        # loop through the elements of the line
                        for i in range(0,len(split_line)):
                            # if the element is setting speed
                            if "F" in split_line[i]:
                                # change the speed
                                split_line[i] = 'F' + str(int(print_speed)*60)
                            # if the element is setting the extrusion amt
                            elif "E" in split_line[i]:
                                # grab the original extrusion amount
                                extruded_amount = split_line[i][1:]
                                # calculate the new extrusion amount
                                # original*new flow/previous flow
                                new_extrude_amount = round(float(extruded_amount)*float(extrusion_multiplier)/previous_extrusion_multiplier,5)
                                # TO DO: round to 5 decimals
                                split_line[i] = 'E' + str(new_extrude_amount)
                        # rewrite the line
                        gcode_data[line_counter] = " ".join(split_line) + '\n'
                    elif " E" in code_line:
                        split_up = code_line.split()
                        for i in range(0, len(split_up)):
                            if "E" in split_up[i]:
                                extrusion = float(split_up[i][1:])
                                new_extrude_amount = round(extrusion*float(extrusion_multiplier)/previous_extrusion_multiplier,5)
                                split_up[i] = "E" + str(new_extrude_amount)
                        gcode_data[line_counter] = " ".join(split_up) + '\n'             
                line_counter += 1
            # Set Extrusion Multiplier
            gcode_data[67] = "M221 S100\n"
            # Set Fan Speed Percent
            gcode_data[70] = "M106 S" + str(int(255*int(fan_speed_percent)/100)) + '\n'
            # Add the parameter settings to the file as comments
            print_parameters = [None] * 8
            print_parameters[0] = "; GCode auto generated from a base CURA file\n"
            print_parameters[1] = "; Kevin Hunt\n"
            print_parameters[2] = "; Print Number " + print_number + ' of 324\n'
            print_parameters[3] = "; Nozzle Temperature " + nozzle_temperature + ' C\n'
            print_parameters[4] = "; Bed Temperature " + bed_temperature + ' C\n'
            print_parameters[5] = "; Print Speed " + print_speed + ' mm/sec\n'
            print_parameters[6] = "; Extrusion Multiplier " + extrusion_multiplier + ' %\n'
            print_parameters[7] = "; Fan Speed " + fan_speed_percent + " %\n"
            gcode_data = print_parameters + gcode_data
            # write the new gcode to a file
            new_filename = "gcode folder\\3_sides_of_a_cube_print_" + print_number + ".gcode"
            with open(new_filename, 'w') as new_gcode:
                new_gcode.writelines(gcode_data)

At the time of printing, the ambient room temperature and ambient humidity was recorded. Prints are sent to the Lulzbot Mini printer using Pronterface. Approximately 141 prints varying the printing parameters were completed in this internship.

2.3.1 Refine the .STL

Our original 3 sides of a cube .STL file has very few vertices. However, our gcode is segmented every 2 mm, thus we want the geometric properties of that exact location in the printed part. Thus, we need more vertices than OpenSCAD created for this part. More vertices can be obtained by refining the model using NetFabb.

NetFabb has an option to import a mesh, and, for our purposes, define a maximum edge length to refine the model. We used a student version of this software for this project.
A step-by-step process goes as follows:

• File -> Open Project
  
• Prepare -> Repair Part

Under the new menu select:

• Mesh Edit -> Refine Triangle Mesh

Setting the first value, Max. Edge Length, to 0.1mm (or the desired max).

It’s important afterwards to click “Apply Repair” in the bottom left and remove the old part.

Now, go to:

• File -> Export Path -> as STL
  
  
  

Refining the model with NetFabb with 2mm maximum edge length (not 0.1mm).

The reason we have chosen 0.1mm for the maximum edge length is to ensure that every gcode segment maps almost exactly to an STL vertex. Once refined to a max edge length of 0.1mm, there were around ~3,000,000 vertices in our STL file (Hence the need for code optimizations later on).

2.3.2 Use Libigl

In order to obtain geometric information for each vertex of the .STL, a software package known as Libigl was used; it was designed specifically to obtain geometric data from meshes. They provide an in-depth tutorial and installation instructions and an official GitHub page.

• Note: We did not install this on Windows for obvious reasons - and also Libigl uses CMake to create executable files, so Windows required a Visual Studio compiler which is both a pain to set up and potentially not free. We recommend running on MacOS/Linux for this part.

Within Libigl, we looked specifically at the vertex normals, Gaussian Curvature, and Discrete Mean Curvature (along with X, Y, Z coordinates). While Libigl is written in C++, it has bindings in Python which made working these values into the training data much easier later on.

Libigl has done a fantastic job at minimizing the amount of code you need to write. In addition, there are Python tutorials on their Github for each of the three types of geometric data we’ve used.

For example, code to calculate and view the vertex normals looks as follows:

import sys, os  
  
# Add the igl library to the module's search path
sys.path.insert(0, '/Users/person/.../libigl/python/') # Add own path 
import pyigl as igl # Libigl's python bindings  
  
stl_path = '/path/to/STL file'  
  
# Libigl uses the eigen class:
# More info found here: https://eigen.tuxfamily.org/dox/group__TutorialMatrixClass.html
OV = igl.eigen.MatrixXd()
OF = igl.eigen.MatrixXi()  
  
def compute_vertex_normals():
    # Syntax for computing normals
    igl.per_vertex_normals(V, F, igl.PER_VERTEX_NORMALS_WEIGHTING_TYPE_AREA, N_vertices)  
      
    return N_vertices
  
if __name__ == "__main__":  
  
    V = igl.eigen.MatrixXd()
    F = igl.eigen.MatrixXi()
    SVI = igl.eigen.MatrixXi()
    SVJ = igl.eigen.MatrixXi()  
    N_vertices = igl.eigen.MatrixXd()
  
    # Read in mesh and remove duplicate vertices, since multiple triangles
    # can contain the same vertex - reducing redundant training data
    igl.read_triangle_mesh(stl_path, OV, OF)
    igl.remove_duplicate_vertices(OV, OF, 0.0, V, SVI, SVJ, F)  
  
    N_vertices = compute_vertex_normals()
  
    # Vertex normals will be stored in N_vertices
    print(N_vertices)

Code for all the geometric data is also below:

import sys, os  
  
# Add the igl library to the modules search path
sys.path.insert(0, '/Users/person/.../libigl/python/') # Add own path
import pyigl as igl  
  
stl_path = '/path/to/STL file'   
  
OV = igl.eigen.MatrixXd()
OF = igl.eigen.MatrixXi()    
  
# Compute vertex normals  
def compute_vertex_normals():
    igl.per_vertex_normals(V, F, igl.PER_VERTEX_NORMALS_WEIGHTING_TYPE_AREA, N_vertices)
    
    return N_vertices
  
# Compute Gaussian Curvature  
def compute_gaussian():
    K = igl.eigen.MatrixXd()
    igl.gaussian_curvature(V, F, K)  
  
    return K  
# Compute Discrete Mean Curvature - in Libigl's tutorial
def compute_dmc():  
  
    HN = igl.eigen.MatrixXd()
    L = igl.eigen.SparseMatrixd()
    M = igl.eigen.SparseMatrixd()
    Minv = igl.eigen.SparseMatrixd()  
  
    igl.cotmatrix(V, F, L)
    igl.massmatrix(V, F, igl.MASSMATRIX_TYPE_VORONOI, M)  
  
    igl.invert_diag(M, Minv)  
  
    # Laplace-Beltrami of position
    HN = -Minv * (L * V)  
  
    # Extract magnitude as mean curvature
    H = HN.rowwiseNorm()  
  
    # Compute curvature directions via quadric fitting
    PD1 = igl.eigen.MatrixXd()
    PD2 = igl.eigen.MatrixXd()  
  
    PV1 = igl.eigen.MatrixXd()
    PV2 = igl.eigen.MatrixXd()  
  
    igl.principal_curvature(V, F, PD1, PD2, PV1, PV2)  
  
    # Mean curvature
    H = 0.5 * (PV1 + PV2)  
  
    return H  
  
  
if __name__ == "__main__":  
  
    V = igl.eigen.MatrixXd()
    F = igl.eigen.MatrixXi()
    SVI = igl.eigen.MatrixXi()
    SVJ = igl.eigen.MatrixXi()  
    N_vertices = igl.eigen.MatrixXd()   
    
    igl.read_triangle_mesh(stl_path, OV, OF)
    igl.remove_duplicate_vertices(OV, OF, 0.0, V, SVI, SVJ, F)  
  
    N_vertices = compute_vertex_normals()
    K = compute_gaussian() # K stores Gaussian Curvature for each vertex
    H = compute_dmc() # H stores DMC values for each vertex
  
    C = igl.eigen.MatrixXd()
    igl.jet(H, True, C) # Convert into RGB values  
  
    # Libigl's viewer
    viewer = igl.glfw.Viewer()
    viewer.data().set_mesh(V, F)
    viewer.data().set_colors(C)
    viewer.launch()

Set up to view the Discrete Mean Curvature values, what appears should look something like this:

Discrete Mean Curvature Visualization (on 2mm edge length mesh)

2.3.3 Other Geometric Information

In addition to the geometric data obtained from Libigl, we felt angles between the Z axis, as well as print direction might also be important to include. Angles between the Z axis might be able to better tell us about overhangs, and print direction (defined as a vector facing out of the front of the printer) might give us information about how orientation affects quality.

The angles found between the Z axis can be seen in the following picture:

We can use the vertex normals found through Libigl to find both \(\theta\) and \(\phi\) and write functions that return these angles for every vertex in the STL file.

The print angle was found just by using the standard angle between two vectors equation:

\(\alpha = cos^{-1}(\frac{\textbf{u} \cdot \textbf{v}}{|\textbf{u}| |\textbf{v}|})\)

Where \(v\) is the vector \((0, 1, 0)\) - facing out of the print bed. In the code, compute_theta_vert() computes \(\theta\) for each vertex, and compute_phi_vert() computes \(\phi\) for each vertex.

The code for this looks as follows, and can also be visualized using Libigl’s viewer:

• Note: N_vertices is still found from Libigl, so running this code separate from the above Libigl code will result in an error. We have all the geometric data calculation functions in one file.

import numpy as np
import math  
  
# Add the igl library to the modules search path
sys.path.insert(0, '/Users/person/.../libigl/python/') # Add own path 
import pyigl as igl
from iglhelpers import e2p  
  
# Compute theta from vertical (shown in picture)
def compute_theta_vert():    
    
    theta_vert = igl.eigen.MatrixXd(len(numpy_vertices), 1)  
  
    # Iterate over entire array calculating theta for each vertex
    for i in range(len(numpy_vertices)):
        if numpy_vertices[i][1] == 0:
            theta_vert[i] = 0
        else:
            theta_vert[i] = math.atan2(numpy_vertices[i][1], numpy_vertices[i][0])  
  
    return theta_vert  
  
# Compute phi from vertical (shown in picture)
def compute_phi_vert():  
  
    phi_vert = igl.eigen.MatrixXd(len(numpy_vertices), 1)  
    for i in range(len(numpy_vertices)):
        phi_vert[i] = math.acos(numpy_vertices[i][2] / np.linalg.norm(numpy_vertices[i]))  
  
    return phi_vert  
  
# Compute print angle according to above equation
def compute_angle_print():  
  
    print_angle = np.array([0, 1, 0])
    angle_print = igl.eigen.MatrixXd(len(numpy_vertices), 1)  
  
    for i in range(len(numpy_vertices)):
        angle_print[i] = math.acos(np.dot(print_angle, numpy_vertices[i]) / np.linalg.norm(numpy_vertices[i]))
  
    return angle_print  
      
if __name__ == "__main__":  
  
    numpy_vertices = e2p(N_vertices) # Libigl's helper function to convert eigen matrices into numpy arrays
    
    tv = compute_theta_vert()
    pv = compute_phi_vert()
    ap = compute_angle_print()  
      
    C = igl.eigen.MatrixXd()
    igl.jet(tv, True, C) # Convert into RGB values  
      
    viewer = igl.glfw.Viewer()
    viewer.data().set_mesh(V, F)
    viewer.data().set_colors(C)
    viewer.launch()

Theta Vertex Visualization (on 2mm edge length mesh)

Obtaining geometric data for all of the vertices took around 15 minutes. However, here we have gathered geometric data about the vertices of the .STL file and what we need to train our model is the geometric information of the gcode points. Thus, we need an efficient way to convert STL vertex information into Gcode segment information. This involved several steps talked about below.

2.3.4 Align Gcode and STL Points

Our first challenge can be perfectly summed up in the picture below.

Printer toolpath in blue. .STL points in black. .STL points need to be scaled to match the toolpath.

Here the gcode toolpath is in light blue and 2 mm refined STL vertices are the black dots. The STL vertices on the top and right side go beyond that of the toolpath. This is an issue as Gcode segments that lie on the corner might not receive the correct ‘corner-like’ information from Libigl.

In CURA, the toolpath is adjusted to take the nozzle width of 0.5 mm and material deposition into account. So the dimensions of the toolpath are 49.5mm x 49.5mm x 49.4mm. In addition, the walls of our original STL file were 1mm wide, however Cura creates a toolpath along the walls of only 0.5mm.

We expect this will be less of an issue for less sharp and thin models. While our ‘3 sides of a cube’ is simple and easy to print, it’s somewhat of a double-edged sword as it also means in order to perfectly capture the geometries, we need to exactly align the STL file and Gcode segments.

Initially, we scaled all of the X, Y, Z coordinates to meet the Gcode toolpath, however this did not appropriately scale the wall widths. This approach should be fine for most models as the wall width will not be exposed. We remedied this issue for our part by simply adjusting the CAD model to reflect the Gcode path. We realize this isn’t the most generalizable option, however, again, we feel this issue is exacerbated in our model. In the future, it might be worth looking into better ways to match the toolpath with the STL file.

Now the dimensions of .STL and gcode points line up, but they do not have the same values for X,Y,Z. This was corrected so we could directly associate geometric information from the STL file with Gcode segments. This was done by finding the minimum X, Y, Z coordinates of both the STL and Gcode files and finding the difference between each of them. Then, we added that difference to the STL file to obtain the same coordinate system.

• Note: If you rotate the model within CURA, you need to export it both as Gcode as well as STL and realign coordinate systems based on this new STL file. The ‘add an offset’ method will not work if you have an old STL file with a rotated Gcode toolpath. More information about this code can be found in the database section below.

2.3.5 Map Gcode to Nearest STL

Now we have the gcode and stl points so they correctly overlap, now we need to map the gcode points to the closest STL point and grab the geometric information of that STL point. Here is a visual representation of what we are doing. We have colored the “gcode” by it’s associated vertex. This indicates that our mapping algorithm works as expected.

• Note: To produce this image, the gcode was cut up into 0.1 mm segments and mapped to the nearest vertex refined to 2 mm.

G-code toolpath colored by the matching .STL point. Evidence the .STL to G-code mapping works.

There’s no guarantee that every STL point will directly correspond to every other Gcode point. This is the case because we’re segmenting the STL file only every 0.1mm, so there will always be that slight misalignment. In fact, only a few STL points directly corresponded to Gcode points, so we needed a way to find the closest vertex to each Gcode segment.

This was done through what’s known as a KD Tree. Scikit Learn - a fantastic package for all things Machine Learning in Python - uses a KD Tree for their K-Nearest Neighbors algorithm. We were able to take the functionality they provide and use it for our purposes(implementation/examples)

The code below uses the KD Tree to associate STL vertices with Gcode segments. The remove_duplicates() function, well, removes the duplicate vertices in the STL file (as a result of multiple trianges containing the same vertex). The reorient() function adds the offset described in the ‘Correctly Align the Gcode and STL Files’ section above. Finally the KD Tree is created in the get_KDTree() function.

# Given an STL file and spliced Gcode segments, this program matches each of
# the given GCode segments to the closest vertex in the STL file.
import sys, os
sys.path.insert(0, '/Users/person/.../libigl/python/') # Add own path  
  
import pandas as pd
from sklearn.neighbors import KDTree
import _pickle as cPickle
import pyigl as igl
from iglhelpers import e2p
import time  
  
IMPORT_TREE = True  
  
stl_path = '/path/to/STL file'
gcode_path = '/path/to/segmented gcode file'  
tree_path = '/path/to/pre-built tree' # Only used after running once  
    
def remove_duplicates():  
  
    # Use libigl's function to read in STL file
    OV = igl.eigen.MatrixXd()
    OF = igl.eigen.MatrixXi()
    igl.read_triangle_mesh(stl_path, OV, OF)  
  
    V = igl.eigen.MatrixXd()
    F = igl.eigen.MatrixXi()
    SVI = igl.eigen.MatrixXi()
    SVJ = igl.eigen.MatrixXi()  
  
    igl.remove_duplicate_vertices(OV, OF, 0.0, V, SVI, SVJ, F);
    print("\nFinished removing duplicates.")
    print("The mesh now has " + str(V.rows()) + " triangles.")  
  
    return V  
  
def reorient(V):  
  
    vertices = e2p(V) # Libigl helper function to create numpy array from MatrixXd
    df_vertices = pd.DataFrame(data = vertices, columns=['X', 'Y', 'Z'])  
  
    df_vertices['X'] = df_vertices['X'].astype(float)
    df_vertices['Y'] = df_vertices['Y'].astype(float)
    df_vertices['Z'] = df_vertices['Z'].astype(float)
    min_stl_x = df_vertices['X'].min()
    min_stl_y = df_vertices['Y'].min()
    min_stl_z = df_vertices['Z'].min()  
  
    df_gcode_arr = pd.read_csv(gcode_path, names=['X', 'Y', 'Z'])
    min_gcode_x = df_gcode_arr['X'].min()
    min_gcode_y = df_gcode_arr['Y'].min()
    min_gcode_z = df_gcode_arr['Z'].min()  
  
    # Calculate distance between min gcode values and min stl values
    dist_x = min_gcode_x - min_stl_x
    dist_y = min_gcode_y - min_stl_y
    dist_z = min_gcode_z - min_stl_z  
  
    # Add that distance (above) as offset to each STL to have the STL file
    # and G-code be in the same coordinate system
    df_vertices['X'] += dist_x
    df_vertices['Y'] += dist_y
    df_vertices['Z'] += dist_z  
  
    print("Finished reorienting.")
    return df_vertices, df_gcode_arr  
  
def get_KDTree(df_vertices, df_gcode_arr):  
  
    if IMPORT_TREE == True:
        tree = pd.read_pickle(tree_path)
    else:
        tree = KDTree(df_vertices)
        print("Saving tree as pickle file...")
        with open(tree_path, "wb") as output_file:
            cPickle.dump(tree, output_file)
        print(" done.")  
  
    df_gcode = pd.DataFrame(columns=['index', 'gx', 'gy', 'gz', 'vx', 'vy', 'vz', 'distance'])  
  
    # Query KD Tree to find nearest vertex for every gcode segment
    total_dist = 0
    for i in range(len(df_gcode_arr.index)):
        dist, ind = tree.query([[df_gcode_arr['X'][i].astype(float),
                                df_gcode_arr['Y'][i].astype(float),
                                df_gcode_arr['Z'][i].astype(float)]], k = 1)  
  
        df_gcode.loc[i] = [i+1,
                           df_gcode_arr['X'][i],
                           df_gcode_arr['Y'][i],
                           df_gcode_arr['Z'][i],
                           df_vertices.ix[ind[0][0], 'X'],
                           df_vertices.ix[ind[0][0], 'Y'],
                           df_vertices.ix[ind[0][0], 'Z'],
                           dist[0][0]
                          ]  
  
        total_dist += dist[0][0]  
  
    average_dist = round(total_dist / len(df_gcode_arr.index), 5)
    print("Average distance from vertex: ",  average_dist)  
  
    return df_gcode  
  
if __name__ == '__main__':  
  
    t0 = time.clock()  
  
    V = remove_duplicates()
    df_vertices, df_gcode_arr = reorient(V)
    df_gcode = get_KDTree(df_vertices, df_gcode_arr)  
  
    df_gcode.to_pickle("/Users/.../df_gcode.pickle")
    df_gcode.to_csv("/Users/.../df_gcode.csv")  
  
    t1 = time.clock()
    print("Run time: " + str(round(t1 - t0, 4)) + "\n")

• Note: Creating the KD Tree for ~3,000,000 vertices takes some time (around 15 minutes), however once created, you can save that tree and load it in - as seen in the first if/else statement of the get_KDTree() function. The actual act of assigning vertices to Gcode segments is highly optimized through the use of a KD Tree, so this part of the code above only takes around 20 seconds.

Including all segments, the average distance to the STL vertex is 0.04mm, which lines up with the fact that we segmented the STL file so that no edge is longer than 0.1mm.

Also, we can plot the geometric data of each gcode point to ensure we have correctly aligned/scaled properly:

G-code points colored by geometric properties from the associated .STL points. Evidence the mapping works as expected.

Now we have a data set where the first 10 elements look like:

• Note: The final training data will ultimately drop the vertex information and distances between Gcode segments and STL vertices.

2.3.6 Database Optimizations

Quick aside: Originally, we were segmenting the STL file every 2mm, however realized that we were often associating vertex information multiple layers above/below that gcode segment with that Gcode segment. This is obviously not ideal, as we really wanted each Gcode segment to be mapped to a vertex on the same level as the previous segment.

This is where we decided to further refine our model to a 0.1mm max edge length resulting in the ~3,000,000 vertices (from ~14,000). Dealing with this much data meant rethinking how we were storing the data; moving away from arrays to a more database-like approach with the Pandas DataFrame

This gave us a massive performance boost in creating our training data as now we were using data structure made for incredibly large amounts of data like our own. You can see the database structure in action in the code above.

2.4.1 Training Image Classifier

Image Classification is an approach that lends itself very well to transfer learning, a technique that modifies a neural network pre-trained to classify a much larger data set. The early layers of the neural network have already learned how to pull out basic features from an image (e.g. lines and textures) and the last few layers do the work of classification. This allows us to replace the last few layers with new starting layers that will be trained to classify our images with our classes. Here, the MATLAB tutorial for using AlexNet for transfer learning was implemented to classify images with Blob, Delamination, None or Warp.

• Note: The classifications of each error are sometimes very subjective. Here are some examples of what we classified as Blob, Warp, Delamination or None. There is a folder of all of our training data images for Image Classification (with subfolders of each error). Our Image Classification model was trained on ~430 segmented pictures of each error. Looking at these pictures will give you the best idea of what we felt constituted each error.

Examples of the Error types used for training the Image Classification Model

To begin, training data for the neural network to learn from must be created. Because we want high refinement of the location of errors in our print, images of printed walls were segmented into > 600 smaller squares. These smaller images were hand classified as Blob, Delamination, None or Warp by placing the image into a folder of that name. Not all of the error types occur at the same rate: ~95% of the images were None. When the network was trained on an unbalanced dataset, the network learned to classify every image as None because it was acheiving 95% accuracy. Hence, the classes needed to be balanced. The dataset was limited by delamination, starting out there were ~200 training images of delamination, so only ~1000 images total could be used to train the network. However, even with this limitation the AlexNet trained by Sean and Evan could achieve 93.8 % training accuracy and 88% test accuracy.

Sean and Evan came up with what we are calling “pseudo-semi-supervised learning”. Semi-supervised learning is the approach of using a model to classify and then using the results of that as new training data. What we are calling pseudo-semi-supervised is manually checking the results and reclassifying where the model was wrong. We were able to scour all of the prints so far and come up with more than twice as many delaminations. As we did this, we also reclassified many of the training images. There were many dark images, blurry images, blue printing material, and prints with white instead of black backgrounds. This training data was not representative of the images coming from our dataset. This is due to the fact that the training data used to train AlexNet was created before most of the pictures were taken so the technique of picture taking changed. Thus, despite the model trained having a decent test/training accuracy, it’s accuracy was poor when the model was used on our dataset. This is an important point; despite the high test accuracy, if the training data was not representative of the real world dataset, the model will have less accuracy in use. This came up a few times after our first round of modifying the training data. The network trained with >98% test accuracy, but when used to classify our full dataset, it did not catch a lot of blobs, or it put a lot of nones in delamination. Adding some more smaller blobs to the training set of Blobs and replacing too washed out or too dark Nones allowed us to achieve a model with 98.26% test accuracy that also showed high accuracy when manually checking the classifications of our dataset.

Other optimizations were added to the training done by us. We were able to see that delaminations mostly occured on slow prints with high nozzle temperature so we started printing these print conditions to try and come up with more delamination images. Training of AlexNet was very slow (~4 hrs) when we adopted it from Sean and Evan, so we added parallelization to the training which dropped the training time down to ~ 1 hour. We also tried various techniques: using ADAM optimizer instead of SGD, when using SGD increasing momentum from 0 to the default of 0.9, and others. The best results (by a hair) were seen with ADAM optimizer and ADAM optimizer trained very quickly ~20 minutes. The final model is created by the code in “transferLearning.m” and saved as “98.26_after_modifications.mat”:

The Image Classification Training curve acheiving 98.26% validation accuracy.

The Image Classification Training confusion matrix with 100% test accuracy

There are some remaining concerns with this model. First, k-fold cross-validation was not used when training and testing this model. K-fold cross-validation adds more confidence in the generalizability of the test accuracy to a new dataset. Secondly, classification depends on what we deemed the errors were. Areas that we deemed as warping, might not be warping. Finally, there are other error types that could be classified (e.g. stringiness, missing material, burning). The images could be manualy moved to new folders of these error types and the MATLAB AlexNet Model slightly modified and retrained to pick up on these errors. Currently, some of these error types are possibly captured by the current approach (i.e stringness and missing material are possibly caught by blobs). This depends on how the human classification treated these errors when manually assigning an error type.

2.4.2 Taking Images

The Image classifier will need to use an edge detection method to find the printed part within the image. To aid the edge detection, a black background a chosen for maximum contrast against the natural ABS. Pictures of the printed part were taken in the order of the bottom with label, then front (XZ) then left (YZ). This order is very important for using a script to automatically prepare the data for classification. Great care was taken to ensure the images were in focus and bright.

2.4.3 Organizing Images

In order for the MATLAB script to know which print number and which face it is working on, the images need to be organized in a specific way. Each Print Number needs it’s own folder named with the Print Number. Inside the folder should be the two faces named PrintNumber_XZ, PrintNumber_YZ.

The faces imaged of the 3-sides-of-a-cube printed parts.

We wrote a python script that could loop over one folder of images direct from the camera and rename and organize the images in this way. The first step is to ensure that the images from the camera are in a single folder and in the order of label, front (XZ), left (YZ), next label, etc. At the beginning of the python script, manually create a list of the print numbers in the order they are in the folder. Upon running the script, it will loop through folder of images and use the list and the order to rename the images as PrintNumber_Label, PrintNumber_XZ, PrintNumber_YZ. All the label images are moved into a single folder together. The two face images are placed into their own folder named by the PrintNumber. With this organization, these images are ready to feed into the Image Classification Network.

"""
Renames images from the camera with their appropriate print number and face plane and moves them to a new file structure in a target folder.
First: Place images in a root folder. Make sure they are in the order of "Label image, XZ face, YZ face"
Second: Create a list below with the print numbers in the appropriate order
Third: Set the holding folder and the target folder.
"""
import os, shutil
root = "path/to/folder of images in order of label, front, left, label"
new_location = "path/to/outfolder"
first_print = 79
last_print = 85
printList =  list(range(first_print,last_print+1))
print(printList)
pictureNames = []
for printNumber in printList:
    pictureNames.append(str(printNumber)+"_Label.jpg")
    pictureNames.append(str(printNumber)+"_XZ.jpg")
    pictureNames.append(str(printNumber)+"_YZ.jpg")
print(pictureNames)
for prints in printList:
    folderName = os.path.join(new_location, str(prints))
    if not os.path.exists(folderName):
        os.mkdir(folderName)
    if not os.path.exists(os.path.join(new_location, "Labels")):
        os.mkdir(os.path.join(new_location, "Labels"))
for dirs,subdirs,images in os.walk(root):
    for i in range(len(images)):
        if ("Label" in pictureNames[i]):
            new_loc = os.path.join(new_location, "Labels") 
        else:
            new_loc = os.path.join(new_location, pictureNames[i].split("_")[0])
        os.rename(os.path.join(root, images[i]), os.path.join(new_loc, pictureNames[i]))

2.4.4 Classifying Point Errors

So far we have discussed how the Image Classifier was trained and how to organize the images for classification. Next, a MATLAB program, “examplescript2.m” will prompt for the input folder of images organized as discussed above and an output folder to save the segment images for future training data. The program will take in images of printed part faces and:

1. Use edge detection to identify the pixel location of the printed part in the image
   
2. Segment the face into >600 square images
   
3. Use the Image Classifier
   
4. Save these classified images in a folder named for the error type for later use as training data
   
5. Use the edge detection to identify the gcode points contained in image, and
   
6. Write these gcode points, print numbers, and error type probabilites to a text file.

The results of running this image classifier on print 11 can be see here:

The locations of the blobs in print 11 and the results of the Image Classification of print 11

The pictures above are visual confirmations that the classification of and mapping to the gcode points works. These error classifications are then taken to be the true error at each segment when we train our prediction model next, but one should note that there is possibility of error in this classification data.

4.1.1 Getting the Predictions

For each of the following sections, we will need all of the possible 10,149 points with all of the possible 432 parameter combinations (there are more than 324 because fan speed 0 was included in the training data). A python script was used to load in the points, do a cartesian join with the possible parameter combinations, normalize the data, and send it all through the model. All of these predictions were then saved to be loaded in for the following tasks.

import numpy as np
import pickle
import pandas
from sklearn import preprocessing
import time
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import matplotlib.cm as cm
# load the image classification data from pickle
data = pandas.read_pickle("/path/to/training data") # Will be pickle file
# drop prints 1 2 3 from the dataframe because they have em 121%
data = data[~data['Print Number'].isin(["1", "2", "3", "167"])].reset_index()
# drop the print number column
data = data.loc[:, data.columns != 'Print Number']
# grab only the points information
points = data[['gx', 'gy', 'gz', 'Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature']].drop_duplicates()
print("There are ", str(points.shape[0]), " points in this object")
# get the unique values of printings parameters
# 'Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed'
possible_NozzleTemps = data['Nozzle Temperature'].unique()
possible_bedTemp = data['Bed Temperature'].unique()
possible_printSpeed = data['Print Speed'].unique()
possible_extrusionMultiplier = data['Extrusion Ratio'].unique()
possible_fanSpeed = data['Fan Speed'].unique()
# all possible combos of the changing parameters
potential_param_combos = pandas.DataFrame([[a,b,c,d,e] for a in possible_NozzleTemps for b in possible_bedTemp\
                                            for c in possible_printSpeed for d in possible_extrusionMultiplier for e in possible_fanSpeed])
print("There are ", str(potential_param_combos.shape[0]), " possible printing parameter combinations")
potential_param_combos.columns = ['Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed']
# merge the points with all the combos
print("...Merging the points with the combinations...")
points['key'] = 0
potential_param_combos['key'] = 0
points_to_check = pandas.merge(points, potential_param_combos, on='key').drop('key', axis=1)
# normalize the data
column_names = points_to_check.columns
min_max_scaler = preprocessing.MinMaxScaler()
np_scaled = min_max_scaler.fit_transform(points_to_check)
send_to_model = pandas.DataFrame(np_scaled)
send_to_model.columns = column_names
# reorder the columns to match original order
# 'Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed'
send_to_model = send_to_model[['gx', 'gy', 'gz', 'Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature', 'Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed']]
# load the model from disk
print("...Loading the model...")
filename = '/path/to/gradientBoostingClassifier_model_dropEM121.sav' # Will be found in code folder
loaded_model = pickle.load(open(filename, 'rb'))
# run the points to check through the model
print("...Making point predictions...")
results = loaded_model.predict_proba(send_to_model)
print("...Done making predictions")
# add the probability predictions and a classification column to the non normalized data
results = pandas.DataFrame(results)
results.columns = ['Blob', 'Delamination', 'None', 'Warp']
results_class = results.idxmax(axis=1)
points_to_check['Predicted Error'] = results_class
points_to_check = pandas.concat([points_to_check.reset_index(drop=True), results], axis=1)
# save this dataframe as AllPossiblePredictions
points_to_check.to_pickle("AllPossiblePredictions.pickle")

Transposed for viewing.

4.1.2 Predicting Future Prints

In picking the prints we wanted to test, we chose 4 of the prints at random that we hadn’t yet done - as to remain as objective as possible. The Error Predictions for these prints were obtained from our model.

import pickle
import pandas
# params for print
print_number = 138
nt = 245
bt = 100
ps = 60
em = 90
fs = 80
# load in the prediction data
points_to_check = pandas.read_pickle("AllPossiblePredictions.pickle")
# points_to_check = points_to_check.loc[0:500000,:]
# grab the row indexes for points of print of interest
print_of_interest = points_to_check[(points_to_check['Fan Speed']==fs) \
                            & (points_to_check['Bed Temperature']==bt)\
                            & (points_to_check['Nozzle Temperature']==nt)\
                            & (points_to_check['Print Speed']==ps)\
                            & (points_to_check['Extrusion Ratio']==em)]
             
# save the predictions
print_of_interest.to_pickle('Predictions of Print ' + str(print_number) + '.pickle')

After printing, we photographed each of them, and ran them back through the Image Classifier to obtain the “true” errors at each segment. This could then be compared with our Error Prediction model’s output for accuracy.

The prints chosen at random were:

• Print 161:
  • Nozzle Temperature: 235 °C
  • Bed Temperature: 100 °C
  • Print Speed: 60 mm/s
  • Extrusion Multiplier: 110%
  • Fan Speed: 80%
• Print 173
  • Nozzle Temperature: 245 °C
  • Bed Temperature: 120 °C
  • Print Speed: 60 mm/s
  • Extrusion Multiplier: 100%
  • Fan Speed: 30%
• Print 239
  • Nozzle Temperature: 245 °C
  • Bed Temperature: 100 °C
  • Print Speed: 40 mm/s
  • Extrusion Multiplier: 100%
  • Fan Speed: 30%
• Print 274
  • Nozzle Temperature: 245 °C
  • Bed Temperature: 100 °C
  • Print Speed: 20 mm/s
  • Extrusion Multiplier: 90%
  • Fan Speed: 45%

The accuracy score for the prediction was determined by calculating how many of the XYZ errors were accurately predicted compared to the Image Classification.

import glob
import pandas as pd
import numpy as np
from sklearn.metrics import accuracy_score  
  
list_dir = glob.glob("/path/to/folder of predictions/*.pickle")
df_alex = pd.read_csv("/path/to/defects predictions.txt", header = None)
df_alex.columns = ['gx', 'gy', 'gz', 'None', 'Blob', 'Delamination', 'Warp', 'Print Number']
y_class = df_alex[['Warp', 'None', 'Blob', 'Delamination']].idxmax(axis = 1)
y_class = pd.DataFrame(data = y_class, columns = ['Classification'])
df_alex = pd.concat([df_alex, y_class], axis = 1)  
  
avg_list = []  
  
for element in list_dir:
    df = pd.read_pickle(element).reset_index()
    string = element.strip().split(".")
    print_number = string[0][-3:]  
  
    df_subset = df_alex.loc[df_alex['Print Number'] == print_number].reset_index()
    accuracy = accuracy_score(df['Predicted Error'], df_subset['Classification'])
    print(print_number, accuracy)  
  
    avg_list.append(accuracy)  
  
print("AVERAGE ACCURACY: ", np.mean(avg_list))  

More information about these results can be found in the Results section.

4.1.3 Optimizing for Quality, Speed or Errors

With our trained model, we would like to predict the printing parameters to optimize for speed, quality or errors. In all of the optimizations below, we chose to hold Bed Temperature and Nozzle Temperature constant. Holding bed temperature constant was decided upon due to global effect of varying the bed temperature throughout the print. Holding the Nozzle Temperature constant is due to the physical constraints of the printer.

The results of a temperature ramp while printing in the worst conditions for nozzle heating.

Here we see, under the worst printing conditions for nozzle heating, it takes ~42 seconds to ramp up 10 degrees Celsius. At our slowest printing conditions, 20 mm/sec, one face layer is completed in 2.5 seconds. Waiting from segment to segment for even a second drastically affects the print quality. Thus, it is outside the current printers capabilities to change the nozzle temperature within a single layer.

Holding the nozzle temperature constant throughout the print does not seem to lower the probability of None Error and in the case of the optimizing for speed, there does not seem to be a loss in the print speeds predicted.

Allowing the temperature to change (top) vs holding the temperature constant (bottom) while optimizing for quality does not affect the predicted probability of None Type Error.

Allowing the temperature to change (top) vs holding the temperature constant (bottom) while optimizing for speed does not affect the predicted probability of None Type Error nor the Print Speeds chosen

In these optimizations, we have done a very simple approach for choosing the best parameters. In the future, a more advanced optimization approach may provide a better print; printer capabilities for changing extrusion ratio and nozzle temperature and/or desired structural properties could be included when deciding on the best local printing parameters.

4.1.3.1 Quality

In order to optimize for quality, the best overall Bed Temperature and Nozzle Temperature was found. “Best” was determined by finding the Bed Temperature and Nozzle Temperature combination that contained the most None Predictions of the 10,149 points. This Bed Temperature and Nozzle Temperature was kept constant throughout the decisions.

Points that could never be None error at this Bed Temperature and Nozzle Temperature had their Print Speed, Extrusion Ratio and Fan Speed selected for the lowest probability of Blob. Points that could be None at this Bed Temperature and Nozzle Temperature had their Print Speed, Extrusion Ratio and Fan Speed selected for the highest probability of None.

import numpy as np
import pickle
import pandas
from sklearn import preprocessing
import time
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import matplotlib.cm as cm  
  
# load in the prediction data
points_to_check = pandas.read_pickle("AllPossiblePredictions.pickle")
# points_to_check = points_to_check.loc[0:500000,:]  
  
# get the total number of points in the objects
points = points_to_check[['gx', 'gy', 'gz', 'Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature']].drop_duplicates()
print("There are ", str(points.shape[0]), " points in this object")  
  
# create a boolean for if there error is a None or not
points_to_check['isNoneError'] = points_to_check['Predicted Error'] == 'None'  
  
# how many points are never none 
bad_points = points_to_check.groupby(['gx','gy','gz']).filter(lambda s: s.isNoneError.sum()==0)[['gx', 'gy', 'gz']].drop_duplicates()
number_bad_points = bad_points.shape[0]
potential_printing_accuracy = round(float((points.shape[0] - number_bad_points)/points.shape[0]*100),2)
print("There are ", str(number_bad_points), " points in this object that are predicted to always print poorly.")
print("The potential printing accuracy is ", str(potential_printing_accuracy), "%")  
  
# these columns aren't needed
good_points = points_to_check.drop(columns=['Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature'])  
  
# now pick the best bed temp and nozzle temp
best_points = good_points[good_points.groupby(['gx', 'gy', 'gz', 'Bed Temperature', 'Nozzle Temperature'])['None'].transform(max) == good_points['None']]
bed_temp_prediction_probs = best_points.groupby(['Bed Temperature', 'Nozzle Temperature'])['None'].agg(["mean", "std"])  
  
# first get all the points at these parameters
good_points = good_points[(good_points['Bed Temperature']==bed_temp_prediction_probs['mean'].idxmax()[0]) & (good_points['Nozzle Temperature']==bed_temp_prediction_probs['mean'].idxmax()[1])]  
  
# next grab the points that are never none
bad_points = good_points.groupby(['gx','gy','gz']).filter(lambda s: s.isNoneError.sum()==0)  
  
# now create a new data frame of the points that can be none
remaining_points = pandas.merge(good_points, bad_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)  
  
# for the bad points, we need to grab the parameters that give the lowest blob probability
best_bad_points = bad_points[bad_points.groupby(['gx', 'gy', 'gz'])['Blob'].transform(min) == bad_points['Blob']]  
  
# now grab the highest none for the good_points
best_points = remaining_points[remaining_points.groupby(['gx', 'gy', 'gz'])['None'].transform(max) == remaining_points['None']]  
  
# now combine the decisions for the bad points and the good points into one df
best_points = best_points.append(best_bad_points, ignore_index=True)  
  
# plot the decisions
param_list = ['Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed', 'None']
i=1
fig = plt.figure(figsize=(8.5,10))
for parameter in param_list:
    ax = fig.add_subplot(3,2,i, projection = '3d')
    ax.view_init(azim=225)
    line = ax.scatter(best_points['gx'].astype(float).tolist(), best_points['gy'].astype(float).tolist(), best_points['gz'].astype(float).tolist(), c = best_points[parameter].astype(float).tolist(), cmap=cm.get_cmap('viridis', len(best_points[parameter].unique())))
    if parameter == 'None':
        cb = plt.colorbar(line)
        plt.title('Probability of None Type Error')
    else:
        cb = plt.colorbar(line, ticks=best_points[parameter].unique())
        plt.title('Changes in ' + parameter)
    i+=1  
  
plt.show()

Predicted parameters and probability of none type error when optimizing for quality.

4.1.3.2 Speed

In order to optimize for speed, the best overall bed temperature and nozzle temperature with a print speed of 80 mm/sec was found. “Best” was determined by finding the bed temperature and nozzle temperature combination that contained the most None Predictions of the 10,149 points. This bed temperature and nozzle temperature was kept constant throughout the decisions.

First, all points that could have None error at 80 mm/sec print speed at this Bed Temperature and Nozzle Temperature were grabbed. Fan Speed and Extrusion Ratio for each point were selected to maximize the probability of None error.

Points that could not be none at 80 mm/sec then underwent the same procedure at 60 mm/sec. From there, remaining points were grabbed that could produce None error at 40 mm/sec and finally at 20 mm/sec.

Points that could not be None under any of these conditions had their Print Speed, Extrusion Ratio and Fan Speed selected for the highest probability of None at first. However, we later found that instead selecting for the lowest probability of Blob gave better results.

'''
This is very Kludgey code. It should have been functionalized for DRYness. When first getting comfortable with this code, start with optimizing for quality.
'''
import numpy as np
import pickle
import pandas
from sklearn import preprocessing
import time
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import matplotlib.cm as cm  
  
# start timer
t0 = time.clock()  
  
# load in the prediction data
points_to_check = pandas.read_pickle("AllPossiblePredictions.pickle")
# points_to_check = points_to_check.loc[0:500000,:]  
  
# get the total number of points in the objects
points = points_to_check[['gx', 'gy', 'gz', 'Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature']].drop_duplicates()
print("There are ", str(points.shape[0]), " points in this object")  
  
# create a boolean for if there error is a None or not
points_to_check['isNoneError'] = points_to_check['Predicted Error'] == 'None'  
  
# how many points are never none
bad_points = points_to_check.groupby(['gx','gy','gz']).filter(lambda s: s.isNoneError.sum()==0)[['gx', 'gy', 'gz']].drop_duplicates()
number_bad_points = bad_points.shape[0]
potential_printing_accuracy = round(float((points.shape[0] - number_bad_points)/points.shape[0]*100),2)
print("There are ", str(number_bad_points), " points in this object that are predicted to always print poorly.")
print("The potential printing accuracy is ", str(potential_printing_accuracy), "%")  
  
# get the number of good points are within each subgroup that cannot easily be varied (exclude fan speed and extrusion multiplier)
# grab all the rows with "None" error, groupby bed temp, nozzle temp and print speed, but only count each point once
parameter_quality = points_to_check[points_to_check['Predicted Error']=='None'].drop_duplicates(['Bed Temperature', 'Nozzle Temperature', 'Print Speed', 'gx','gy','gz'])
parameter_quality = parameter_quality.groupby(['Bed Temperature', 'Nozzle Temperature', 'Print Speed']).size().to_frame(name='count').reset_index()
parameter_quality['percent'] = parameter_quality['count'].apply(lambda x: round(float(x/points.shape[0]*100),2))  
  
# print(parameter_quality.sort_values(['Print Speed', 'percent'], ascending=[False, False]))  
  
print("The best parameter combination at 80 print speed is: ")
best_combo = parameter_quality[parameter_quality['Print Speed']==80].ix[parameter_quality[parameter_quality['Print Speed']==80]['percent'].idxmax()]
print(best_combo)  
  
# using the top parameter combo as a guide, what are the good points
good_points = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==best_combo['Print Speed'])]
good_points = good_points.drop_duplicates(['gx','gy','gz'])[['gx', 'gy', 'gz']]  
  
print(good_points.shape[0])  
  
# we now have the good points and the bad points, where are the potential gain points
remaining_points = points_to_check[['gx', 'gy', 'gz']].drop_duplicates(['gx', 'gy', 'gz'])
# drop the points that are good
remaining_points = pandas.merge(remaining_points, good_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
# drop the points that will always be bad
remaining_points = pandas.merge(remaining_points, bad_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
  
print("By locally varying the print parameters, we can gain ", str(remaining_points.shape[0]), " more None type points")  
  
# select the printing parameter combinations that these points are good at
def good_param_combos():
    remaining_point_data = pandas.merge(points_to_check, remaining_points, on=['gx', 'gy', 'gz'], how='outer', indicator=True).query('_merge=="both"').drop('_merge', axis=1)
    remaining_points_constant_bed_temp = remaining_point_data[(remaining_point_data['Bed Temperature']==best_combo['Bed Temperature']) & (remaining_point_data['Predicted Error']=='None')]
    num_remaining_points_constant_bed_temp = remaining_points_constant_bed_temp.drop_duplicates(['gx','gy','gz']).shape[0]  
  
    remaining_points_quality = remaining_points_constant_bed_temp[remaining_points_constant_bed_temp['Predicted Error']=='None'].drop_duplicates(['Bed Temperature', 'Nozzle Temperature', 'Print Speed', 'gx','gy','gz'])
    remaining_points_quality = remaining_points_quality.groupby(['Bed Temperature', 'Nozzle Temperature', 'Print Speed']).size().to_frame(name='count').reset_index()
    remaining_points_quality['percent of remaining points'] = remaining_points_quality['count'].apply(lambda x: round(float(x/remaining_points.shape[0]*100),2))  
  
    print("However, we cannot change bed temperature. Without varying bed temperature, we can gain ", str(num_remaining_points_constant_bed_temp), " more None type points")
    print("These points are good at the following print parameters:")
    print(remaining_points_quality.sort_values('percent of remaining points', ascending=False).sort_values('Print Speed', ascending=False))  
  
good_param_combos()  
  
print("Of these we can gain a decent amount of points by not changing Nozzle Temperature and dropping speed down to 60.")  
  
good_points_half = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==60)]
good_points_half = good_points_half.drop_duplicates(['gx','gy','gz'])[['gx', 'gy', 'gz']]
#drop points that were already in good_points
good_points_half = pandas.merge(good_points_half, good_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("We have gained ", str(good_points_half.shape[0]), " good points")  
  
# So now how many points do we have left to gain
# drop the points that are now good
remaining_points = pandas.merge(remaining_points, good_points_half, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("There are still ", str(remaining_points.shape[0]), " points that can be gained by varying parameters.")  
  
print("Of these we can gain a decent amount of points by not changing Nozzle Temperature and dropping speed down to 40.")  
  
good_points2 = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==40)]
good_points2 = good_points2.drop_duplicates(['gx','gy','gz'])[['gx', 'gy', 'gz']]
#drop points that were already in good_points and good_points_half
good_points2 = pandas.merge(good_points2, good_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
good_points2 = pandas.merge(good_points2, good_points_half, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("We have gained ", str(good_points2.shape[0]), " good points")  
  
# So now how many points do we have left to gain
# drop the points that are now good
remaining_points = pandas.merge(remaining_points, good_points2, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("There are still ", str(remaining_points.shape[0]), " points that can be gained by varying parameters.")  
  
# select the printing parameter combinations that these points are good at
good_param_combos()  
  
print("We can gain a decent amount of the remaining points while still not changing bed temperature or nozzle temperature if we drop print speed down to 20")  
  
good_points3 = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==20)]
good_points3 = good_points3.drop_duplicates(['gx','gy','gz'])[['gx', 'gy', 'gz']]
# drop points that were already in good_points and good_points2
good_points3 = pandas.merge(good_points3, good_points, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
good_points3 = pandas.merge(good_points3, good_points_half, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
good_points3 = pandas.merge(good_points3, good_points2, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("We have gained ", str(good_points3.shape[0]), " good points")  
  
# So now how many points do we have left to gain
# drop the points that are now good
remaining_points = pandas.merge(remaining_points, good_points3, how = 'outer', indicator=True).query('_merge=="left_only"').drop('_merge', axis=1)
print("There are still ", str(remaining_points.shape[0]), " points that can be gained by varying parameters.")  
  
# select the printing parameter combinations that these points are good at
good_param_combos()  
  
# for the remaining points, let's get the highest probability of no error without changing the nozzle temp
points_to_fudge = bad_points.append(remaining_points)
print("Because ther are no more at the same nozzle temperature, we will not grab these points.")
print("This leaves a total of ", str(points_to_fudge.shape[0]), " points to fudge.")
print("We will use the probabilites of have a none error to pick the print conditions that give the highest likelihood of no error with a constant Nozzle and Bed Temperature.")  
  
# for the remaining points let's grab the results of points without changing the nozzle temp/bed temp
bad_param_data = points_to_check[(points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])]
points_to_fudge = points_to_fudge.merge(bad_param_data, how='left', on=['gx', 'gy', 'gz'])  
  
# now grab the best parameters for the smallest risk of Blobs
best_bad_points = points_to_fudge[points_to_fudge.groupby(['gx', 'gy', 'gz'])['Blob'].transform(min) == points_to_fudge['Blob']]  
  
# where are these bad points chosen
param_list = ['Print Speed', 'Extrusion Ratio', 'Fan Speed', 'None']
i=1
fig = plt.figure(figsize=(10,8))
for parameter in param_list:
    ax = fig.add_subplot(2,2,i, projection = '3d')
    ax.view_init(azim=225)
    line = ax.scatter(best_bad_points['gx'].astype(float).tolist(), best_bad_points['gy'].astype(float).tolist(), best_bad_points['gz'].astype(float).tolist(), c = best_bad_points[parameter].astype(float).tolist(), cmap=cm.get_cmap('viridis', len(best_bad_points[parameter].unique())))
    if parameter == 'None':
        cb = plt.colorbar(line)
        plt.title('Probability of None Type Error')
    else:
        cb = plt.colorbar(line, ticks=best_bad_points[parameter].unique())
        plt.title('Changes in ' + parameter)
    i+=1  
  
plt.show()  
  
# great, now we need to pick the best extrusion multipliers and fan speeds for the good points above
# so first grab the prediction data that matches the points and parameter conditions chosen above
good_param_data = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==best_combo['Print Speed'])]
good_points = good_points.merge(good_param_data, how='left', on=['gx', 'gy', 'gz'])  
  
good_half_param_data = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==60)]
good_points_half = good_points_half.merge(good_half_param_data, how='left', on=['gx', 'gy', 'gz'])  
  
good2_param_data = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==40)]
good_points2= good_points2.merge(good2_param_data, how='left', on=['gx', 'gy', 'gz'])  
  
good3_param_data = points_to_check[(points_to_check['Predicted Error']=='None') \
                            & (points_to_check['Bed Temperature']==best_combo['Bed Temperature'])\
                            & (points_to_check['Nozzle Temperature']==best_combo['Nozzle Temperature'])\
                            & (points_to_check['Print Speed']==20)]
good_points3= good_points3.merge(good3_param_data, how='left', on=['gx', 'gy', 'gz'])  
  
# now they are ready to be combined
good_points = good_points.append(good_points_half, ignore_index=True)
good_points = good_points.append(good_points2, ignore_index=True)
good_points = good_points.append(good_points3, ignore_index=True)  
  
# now grab the parameters that give the highest probability of None
best_points = good_points[good_points.groupby(['gx', 'gy', 'gz'])['None'].transform(max) == good_points['None']]  
  
best_points = best_points.append(best_bad_points, ignore_index=True)
# these columns aren't needed
best_points = best_points.drop(columns=['Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature'])
print(best_points)  
  
# save this as a pickle
best_points.to_pickle("SmoothingForSpeed-1 temp.pickle") # This will be loaded later for smoothing. Look for the same filename two code chunks down.  
  
# end timer
t1 = time.clock()
print(str(round(t1-t0,2)) + " seconds")  
  
param_list = ['Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed', 'None']
i=1
fig = plt.figure(figsize=(9,11))
for parameter in param_list:
    ax = fig.add_subplot(3,2,i, projection = '3d')
    ax.view_init(azim=225)
    line = ax.scatter(best_points['gx'].astype(float).tolist(), best_points['gy'].astype(float).tolist(), best_points['gz'].astype(float).tolist(), c = best_points[parameter].astype(float).tolist(), cmap=cm.get_cmap('viridis', len(best_points[parameter].unique())))
    if parameter == 'None':
        cb = plt.colorbar(line)
        plt.title('Probability of None Type Error')
    else:
        cb = plt.colorbar(line, ticks=best_points[parameter].unique())
        plt.title('Changes in ' + parameter)
    i+=1  
  
plt.show()  

Here are the locations of the points that could never be None at this specific Bed Temperature and Nozzle Temperature and parameters were chosen to minimize the probability of Blob:

The locations and chosen parameters to reduce predicted Blob probability of points that never have None Error predictions

Here are the parameters chosen for each point in order to optimize speed with the bad points chosen to minimize the probability of Blob:

Predicted parameters and probability of none type error when optimizing for speed.

4.1.3.3 Forcing Errors

In order to force errors, instead of selecting highest None probabilities, as is the case when optimizing for quality, we are selecting for the highest probabilites of Blobs, Warps, or Delaminations.

This is the code for forcing Blobs. To force Warps or Delaminations just change all instances of Blob to Warp or Delamination.

import numpy as np
import pickle
import pandas
from sklearn import preprocessing
import time
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import matplotlib.cm as cm  
  
# load in the prediction data
points_to_check = pandas.read_pickle("AllPossiblePredictions.pickle")
# points_to_check = points_to_check.loc[0:500000,:]  
  
# get the total number of points in the objects
points = points_to_check[['gx', 'gy', 'gz', 'Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature']].drop_duplicates()
print("There are ", str(points.shape[0]), " points in this object") 
  
# create a boolean for if there error is a Blob or not
points_to_check['isBlobError'] = points_to_check['Predicted Error'] == 'Blob'  
  
# how many points are never Blob
bad_points = points_to_check.groupby(['gx','gy','gz']).filter(lambda s: s.isBlobError.sum()==0)[['gx', 'gy', 'gz']].drop_duplicates()
number_bad_points = bad_points.shape[0]
potential_printing_accuracy = round(float((points.shape[0] - number_bad_points)/points.shape[0]*100),2)
print("There are ", str(number_bad_points), " points in this object that are predicted to always print poorly.")
print("The potential printing accuracy is ", str(potential_printing_accuracy), "%")  
  
# these columns aren't needed
good_points = points_to_check.drop(columns=['Vert Norm X', 'Vert Norm Y', 'Vert Norm Z', 'Theta Vert', 'Phi Vert', 'Print Angle', 'Discrete Mean Curve', 'Gaussian Curvature'])  
  
# now grab the highest probability Blobs parameters for each grouop of Bed Temp and Nozzle Temp
best_points = good_points[good_points.groupby(['gx', 'gy', 'gz', 'Bed Temperature', 'Nozzle Temperature'])['Blob'].transform(max) == good_points['Blob']]  
  
# now pick Bed Temp and Nozzle Temp group had the highest Blob Prob  
bed_temp_prediction_probs = best_points.groupby(['Bed Temperature', 'Nozzle Temperature'])['Blob'].agg(["mean", "std"])  
  
# grab that group of data for the parameters
best_points = best_points[(best_points['Bed Temperature']==bed_temp_prediction_probs['mean'].idxmax()[0]) & (best_points['Nozzle Temperature']==bed_temp_prediction_probs['mean'].idxmax()[1])]  
  
# plot
param_list = ['Nozzle Temperature', 'Bed Temperature', 'Print Speed', 'Extrusion Ratio', 'Fan Speed', 'Blob']
i=1
fig = plt.figure(figsize=(8.5,10))
for parameter in param_list:
    ax = fig.add_subplot(3,2,i, projection = '3d')
    ax.view_init(azim=225)
    line = ax.scatter(best_points['gx'].astype(float).tolist(), best_points['gy'].astype(float).tolist(), best_points['gz'].astype(float).tolist(), c = best_points[parameter].astype(float).tolist(), cmap=cm.get_cmap('viridis', len(best_points[parameter].unique())))
    if parameter == 'Blob':
        cb = plt.colorbar(line)
        plt.title('Probability of Blob Type Error')
    else:
        cb = plt.colorbar(line, ticks=best_points[parameter].unique())
        plt.title('Changes in ' + parameter)
    i+=1  
  
plt.show()  
  

Predicted parameters and probability of Blob type error when optimizing for Blobs.

4.1.4 Smoothing

Randomly varying Extrusion Ratio and Print Speed every 2 mm revealed physical limitations of the system.

Randomizing Print Speed (left) and Extrusion Multiplier (right) every 2 mm

We see that changing the Print Speed comprimises the printing quality and we believe, due to the compressible fluid nature of ABS, changes in the extrusion multiplier are not immediately put in to effect. We also believe, this poor quality of changing the Print Speed is due in part to the same issue. Although extrusion multiplier is not changed in the varying print speed part, when the speed increases the rate of extrusion must increase as well. We think this extrusion rate cannot keep pace with the print speed due to the compressible fluid nature of the material. There may also be some mechanical jerking from the rapid change in speed affecting the print quality. We see a similar issue in one of our optimization tests:

The effect of changing Print Speed and Extrusion Multiplier during the print on print quality

Here we see very clear demarcations where Print Speed and Extrusion Multiplier have changed during the print. It is also interesting to note that Print Speed and Extrusion Multiplier tend to change together. The printing path is from the left to the right in the picture of the printed part above. So we see as print speed and extrusion multiplier are supposed to quickly drop, there is an overextrusion of material for roughly 1 cm (5 XYZ segments). Near the right edge of this printed face, Print Speed and Extrusion Multiplier are supposed to increase. Here we see an underextrusion of material for rougly 1 cm.

In order to overcome these issues due to a quick change in Print Speed and/or Extrusion Multiplier, we created a very rough technique for smoothing. This code will look at a window of segments for a single parameter and if the last value is different than the second to last value, the whole window takes incremental steps from the first value of the window to the last. So for a window size of 5:

There are some issues with this approach and future work could be done to improve the smoothing algorithm. For now, this algoithm does a nice job.

Window Size 7 Smoothing

This smoothing is applied to Print Speed, Extrusion Multiplier and Fan Speed.

• Note: The smoothing algorithm will set fan speeds between 0-30, all of these will be realized as Fan Speed 0 due to the physical limitations of the fan.

import numpy as np
import pickle
import pandas as pd
from sklearn import preprocessing
import time
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
import matplotlib.cm as cm  
  
# load the image classification data from pickle
data = pd.read_pickle("./mapped_SmoothingForSpeed-1 temp.pickle") # Output from the code 2 chunks above  
  
# function to perform smoothing
def smooth(df):
    windowsize = 5
    for i in range(1,df.shape[0]):
        if i < windowsize-1:
            continue
        else:
            df_subset = df.iloc[i-windowsize+1:i+1]
            if df_subset.iloc[-1].item() != df_subset.iloc[-2].item():
                new_array = np.linspace(df_subset.iloc[0].item(), df_subset.iloc[-1].item(), num=len(df_subset))
                df.iloc[i-windowsize+1:i+1] = new_array
    return(df)  
  
# function to plot the results of the smoothing
def plots(best_points):
    param_list = ['Print Speed', 'Extrusion Ratio', 'Fan Speed']
    for parameter in param_list:
        fig = plt.figure()
        ax = fig.add_subplot(111, projection = '3d')
        line = ax.scatter(best_points['gx'].astype(float).tolist(), best_points['gy'].astype(float).tolist(), best_points['gz'].astype(float).tolist(), c = best_points[parameter].astype(float).tolist(), cmap=cm.get_cmap('viridis', len(best_points[parameter].unique())))
        cb = plt.colorbar(line)
        plt.title('Changes in ' + parameter)
    plt.show()  
  
# do the smoothing
data['Extrusion Ratio'] = smooth(data['Extrusion Ratio'])
data['Print Speed'] = smooth(data['Print Speed'])
data['Fan Speed'] = smooth(data['Fan Speed'])  
  
# save the results
data.to_pickle("/.mapped_SmoothingForSpeed-1 temp_SMOOTH.pickle")  
  
# plot the results
plots(data)  

Here are the results of applying the smoothing algorithm to the predicted parameters for optimizing speed and quality:

Smoothing the Predictions for Optimizing for Quality (Before: Top, After: Bottom) with a Window Size of 10

Smoothing the Predictions for Optimizing for Speed (Before: Top, After: Bottom) with a Window Size of 7

4.2.1 Outer points to all gcode points

Since we were only able to capture errors of the outer walls (due to the limitations we had in photographing), it was necessary to communicate the optimal printing parameters found throughout our model on the outer walls to the rest of the print (inner walls and bottom).

To associate the outer wall segments with the inner wall segments, we used another KD Tree to find the nearest outer wall segment to each inner wall segment. We then used Pandas to merge information about printing parameters with the associated Gcode segments.

For the entire bottom (including 3 layers), we looked at the predictions of our model for the bottom 3 layers of the outer wall and generalized the parameters that were seen most often there. Thus the bottom three layers of our print were always constant.

Code for this can be seen below:

from sklearn.neighbors import KDTree
import pandas as pd
import numpy as np    
  
BOTTOM_COUNT = 3  
  
def gcode_tree():  
  
    best_print = pd.read_pickle("/.mapped_SmoothingForSpeed-1 temp.pickle")  
  
    # Get the Z value of the BOTTOM_COUNT (3rd) layer
    df_z_bottom = best_print['gz'].unique().astype(float)
    df_z_bottom = pd.DataFrame(data = df_z_bottom, columns = ['Z'])
    df_z_bottom = df_z_bottom.sort_values(['Z']).reset_index()
    max_bottom = df_z_bottom['Z'].nsmallest(BOTTOM_COUNT)[BOTTOM_COUNT - 1]  
  
    all_segs = pd.read_csv("path/to/segmented_gcode_points.txt", header = None, names = ['gx', 'gy', 'gz'])  
  
    outer_segs = best_print[['gx', 'gy', 'gz']]  
  
    # Create a KD Tree of just the outer points and query every point segment of gcode
    gcode_tree = KDTree(outer_segs)
    dist, ind = gcode_tree.query(all_segs, k = 1)
    ind = ind.ravel()  
  
    df_dist = pd.DataFrame(data = dist).reset_index().drop(['index'], axis = 1)
    df_nn = pd.DataFrame(data = outer_segs.iloc[ind]).reset_index().drop(['index'], axis = 1)
    df_dist.columns = ['Distance']
    df_nn.columns = ['nx', 'ny', 'nz']  
  
    # Make a dataframe of each gcode segment and it's nearest outer wall point
    # Eventually removing the bottom 3 layers from this
    gcode_map = pd.DataFrame(data = all_segs, columns = ['gx', 'gy', 'gz'])
    gcode_map = pd.concat([gcode_map, df_nn, df_dist], axis = 1)
    # print(gcode_map) -- might be useful here  
  
    print_params = best_print[['gx', 'gy', 'gz',
                               'Nozzle Temperature', 'Bed Temperature', 'Print Speed',
                               'Extrusion Ratio', 'Fan Speed']].reset_index().drop(['index'], axis = 1)  
  
    print_params = print_params.rename(index = str, columns = {"gx": "nx", "gy": "ny", "gz": "nz"})  
  
    # Merge the gcode mapping with the associated print parameters of the outer wall
    df_with_bottom = pd.merge(gcode_map, print_params, how = 'left', on = ['nx', 'ny', 'nz'])  
  
    # Cut out predictions for bottom and apply blanket prediction
    df_without_bottom = df_with_bottom.loc[df_with_bottom['gz'] > max_bottom]
    df_bottom = df_with_bottom.loc[df_with_bottom['gz'] <= max_bottom]  
  
    df_bottom['Nozzle Temperature'] = df_bottom['Nozzle Temperature'].value_counts().idxmax().astype(float)
    df_bottom['Bed Temperature'] = df_bottom['Bed Temperature'].value_counts().idxmax().astype(float)
    df_bottom['Print Speed'] = df_bottom['Print Speed'].value_counts().idxmax().astype(float)
    df_bottom['Extrusion Ratio'] = df_bottom['Extrusion Ratio'].value_counts().idxmax().astype(float)
    df_bottom['Fan Speed'] = df_bottom['Fan Speed'].value_counts().idxmax().astype(float)  
  
    df_final = pd.concat([df_bottom, df_without_bottom])
    df_final = df_final.drop(['nx', 'ny', 'nz', 'Distance'], axis = 1)  
  
    return df_final  
  
if __name__ == "__main__":
    mapped_gcode = gcode_tree()
    print(mapped_gcode)
    mapped_gcode.to_pickle("./mapped_gcode_SmoothingForSpeed.pickle")

Before and after pictures of this code can also be seen here:

Mapping external parameters to internal walls and applying the most frequent bottom parameter combination to the entire bottom.

• Note: The right picture contains all the gcode segments in our print.

These pictures show that the algorithm did indeed work as intended, and provided a decent way to generalize the predictions of our model on the outer wall to every segment seen in the Gcode.

In merging with Pandas the way we did, we preserved the original order of Gcode segments in the print, thus communicating our parameter combinations to Gcode was made much simpler.

4.2.2 Communicate to gcode

Next, these printing parameters for each segment need to be communicated to the printer via gcode. Here is a script that takes in the original gcode, cuts up the lines, and rewrites the segmented gcode with the new printing parameters for each gcode segment. In order to run this, place any mapped Gcode segments from the above sections into a folder and give the code the path to that folder. It will write the gcode for all of the pickles in that file. The code is very similar to the code used to section up the gcode, it just uses the new predicted parameters to rewrite the gcode instead of the original.

import math
import numpy as np
import pickle
import pandas
import glob, os  
  
segment_size = 2 #mm
offset = 8 # lines
os.chdir("/path/to/Done Mapping Pickles") # (folder)
# load the point parameter info from pickle
for filename in glob.glob('*.pickle'):
    data = pandas.read_pickle(filename)  
      
    bed_temperature = str(data['Bed Temperature'].unique()[0])
    nozzle_temperature = str(data['Nozzle Temperature'].unique()[0])  
  
    with open(".\\3_sides_of_a_cube_print_24.gcode") as base_gcode:
        # load in the base gcode
        gcode_data = base_gcode.readlines()
        gcode_to_write = []  
  
        # Set the overall bed temperature
        gcode_data[15+offset] = "M140 S" + bed_temperature + " ; start bed heating up\n"
        gcode_data[63+offset] = "M190 S" + bed_temperature + " ; wait for bed to reach printing temp\n"
        gcode_data[576+offset] = "M140 S" + bed_temperature + "\n"  
  
        # Set the first nozzle temp
        gcode_data[62+offset] = "M109 R" + nozzle_temperature + " ; wait for extruder to reach printing temp\n"  
  
        # Set First Extrusion Multiplier
        gcode_data[67+offset] = "M221 S100\n"  
  
        # Set First Fan Speed Percent to be off
        gcode_data[71+offset] = "M106 S0\n"  
  
        # Set Print Speed and Extrusion Multiplier
        make_changes = False
        print_finished = False
        line_counter = 0
        previous_extrusion_multiplier = 100.0
        X_prev = 0
        Y_prev = 0
        current_X = 0
        current_Y = 0
        current_Z = 0
        gcode_split_points = []
        for code_line in gcode_data:
            line_to_write = code_line
            # make sure outside of skirt
            # TO DO: changed this, check works
            if ";TYPE:" in code_line and code_line != ";TYPE:SKIRT\n":
                make_changes = True
            # make sure not at end
            if "M104 S0" in code_line:
                make_changes = False
            # store the current X and Y
            if "X" in code_line:
                current_X = float([x[1:] for x in code_line.split() if "X" in x][0])
            if " Y" in code_line:
                current_Y = float([y[1:] for y in code_line.split() if "Y" in y][0])
            if " Z" in code_line:
                current_Z = float([y[1:] for y in code_line.split() if "Z" in y][0])
            #TO DO: prevent infill split up
            # begin changing the line if it's not the skirt or after the print
            if make_changes:
                # if the line is an extrusion and movement command
                if all(x in code_line for x in ['G1', 'X', 'Y', 'E']):
                    # grab the amount of extrudate
                    extrusion = float([e[1:] for e in code_line.split() if "E" in e][0])  
  
                    # start segmenting the line up
                    # a^2 + b^2 = c^2
                    extrusion_distance = math.sqrt((current_X-X_prev)**2+(current_Y-Y_prev)**2)  
  
                    # based on the segment sizes, how many whole segments will be needed
                    number_of_new_segments = max(1,round(extrusion_distance/segment_size))  
  
                    # based on that many segments, how far do we have to step in the X and Y direction each time
                    X_step_distance = (current_X-X_prev)/number_of_new_segments
                    Y_step_distance = (current_Y-Y_prev)/number_of_new_segments
                    # store the new coordinates for the segments
                    new_X_coords = []
                    new_Y_coords = []
                    for i in range(0,number_of_new_segments-1):
                        new_X_coords.append(X_prev + X_step_distance*(i+1))
                        new_Y_coords.append(Y_prev + Y_step_distance*(i+1))
                    new_X_coords.append(current_X)
                    new_Y_coords.append(current_Y)  
  
                    # segment the extrusion amount up
                    new_extrusion_amount = round(extrusion/number_of_new_segments,5) 
  
                    # write these new instructions out
                    line_to_write = []  
                      
                    import numpy as np
                    for i in range(0,number_of_new_segments):
                        gx = new_X_coords[i]
                        gy = new_Y_coords[i]
                        gz = current_Z  
  
                        # grab the row pertaining to the point from the data
                        current_row = data[(np.isclose(data['gx'],gx)) & (np.isclose(data['gy'],gy)) & (np.isclose(data['gz'],gz))]  
  
                        # grab the parameters for this point
                        nozzle_temperature = str(current_row['Nozzle Temperature'].item())
                        print_speed = current_row['Print Speed'].item()
                        extrusion_multiplier = current_row['Extrusion Ratio'].item()
                        fan_speed_percent = current_row['Fan Speed'].item()  
  
                        # Set Fan Speed Percent in the gcode
                        fan_speed_set = "M106 S" + str(int(255*fan_speed_percent/100)) + ' ; Set Fan Speed Percent\n'  
  
                        # Set Hot end temp in the gcode
                        hot_end = "M104 S" + nozzle_temperature + " ; Set Nozzle Temp and continue without waiting\n"  
                          
                        # calculate the new extrusion amount
                        # MAKE SURE original extrusion amount was 100% !!! This depends on the base gcode used
                        new_extrude_amount = round(new_extrusion_amount*extrusion_multiplier/float(previous_extrusion_multiplier),5)  
  
                        # write the gcode printing command
                        gcode_command = "G1 F" + str(int(print_speed)*60) + " X" + str(gx) \
                        + " Y" + str(gy) + " E" + str(new_extrude_amount) + "\n"  
                          
                        line_to_write.extend([fan_speed_set, hot_end, gcode_command])
                        # if you need to get the new split points unccmment this line and the one writing the file below
                        # gcode_split_points.append(str(new_X_coords[i]) + "," + str(new_Y_coords[i]) + "," + str(current_Z) + "\n")  
  
            # store the previous X and Y
            X_prev = current_X
            Y_prev = current_Y  
  
            gcode_to_write.extend(line_to_write)  
  
        # write the new gcode to a file
        new_filename = filename.split('.')[0] + ".gcode"
        with open(new_filename, 'w') as new_gcode:
            new_gcode.writelines(gcode_to_write)  
  
        # # use to write out the new gcode points
        # with open("segmented_gcode_points.txt", 'w') as points:
        #     points.writelines(gcode_split_points)  

And now we have come full circle, CAD –> .STL –> gcode –> segment match with geometric properites –> run these points through model with various printing parameters –> communicate optimum parameters back into gcode.

5.2.1 Quality

The results of optimizing for quality can be seen below:

Prediction Model error predictions (before smoothing) vs Print vs Image Classification of Errors of the Optimized for Quality Print

We printed this part twice, and it terms of quality, they ranked 6th and 12th out of 146 prints according to our Image Classification Model. This print took 24 minutes 13 seconds to complete. Comparing this with each time for baseline prints, we get:

• 20 mm/s: 45:27
• 40 mm/s: 23:31
• 60 mm/s: 16:03
• 80 mm/s: 12:32

• Quality Optimizations: 24:13

• Note: As mentioned previously, the AlexNet Image Classification Model was dependent on our own interpretation errors. Outputs of these rankings should only be used to roughly get an idea of where the prints fall.

In addition, we can look at the average percentage of each print classified (by our Image Classification model) to not have an error. This can be compared with our optimized-for-quality print below:

• 20 mm/s: 74.21%
• 40 mm/s: 76.04%
• 60 mm/s: 72.20%
• 80 mm/s: 69.59%
• Quality optimizations: 86.41%

The print time of the different Print Speeds and our Optimized for Quality print and their Image Classification Qualities.

This is fantastic news, as this means the use of machine learning in 3D printing is an extremely valid tactic!

5.2.2 Speed

In optimizing for speed, our goal, broadly, was to keep as high of speed as possible while also decreasing errors. Results of our speed optimizations can be seen below:

Prediction Model error predictions (before smoothing) vs Print vs Image Classification of Errors of the Optimized for Speed Print

Our print with local speed optimizations took 16 minutes 56 seconds to complete. This can be compared to baseline prints - each with speed and time:

20 mm/s: 45:27 40 mm/s: 23:31 60 mm/s: 16:03 80 mm/s: 12:32 *Our print: 16:56

Reiterating the average percentages of no error at each print speed from above, and including our optimized-for-speed print, we get:

20 mm/s: 74.21% 40 mm/s: 76.04% 60 mm/s: 72.20% 80 mm/s: 69.59% *Speed optimizations: 76.61%

The print time of the different Print Speeds, Optimized for Quality print, and Optimized for Speed print and their Image Classification Qualities.

The constant print speed percentages fall in line with what we expect - that there’s a sweet spot in speed that isn’t too slow, but isn’t too fast.

More importantly, our print optimizing speed still outranks the average quality of each constant print speed.

Again, this is good news for machine learning in 3D printing!

5.2.3 Errors

Forcing each error type is also a great measure of our model’s predictive power. As a reminder, here is our model’s attempt at forcing blobs to occur:

Predicted printing parameters and probability of Blob Type Error to optimize for blobs.

In creating these localized optimizations, we’ve modified previous segment information in a way that our model can’t account for. What might work best for constant parameters is not guaranteed to work if the parameters are constantly being varied. For example, going quickly from low print speeds to high print speeds around the corners actually had a negative impact on forcing blobs due to fact that quickly speeding up had the effect of thinning out the material.

5.2.4 Comments

Firstly, this is very intriguing as it’s getting into an aspect of 3D printing not yet explored. In addition, our quick and dirty smoothing algorithm significantly improved the performance in this regard. Looking into more advanced smoothing techniques/optimizations will absolutely improve performance. Finally, and what’s beautiful about machine learning in this regard, is that adding this print to the training data (while also adding previous segments’ parameter information) will result in a model that begins to understand these interplay between the parameters with localized optimizations.