Providing a GUI for the spline caculations
Fixed bug with import
Added functionality for when e is pressed, points exported to csv
Added functionality for editing points
Added functionality for adding points
Change-Id: Ied5e6155a8ebfbbd6d2eecba34532245197609a8
diff --git a/y2018/control_loops/python/spline_generate.py b/y2018/control_loops/python/spline_generate.py
new file mode 100644
index 0000000..4376f37
--- /dev/null
+++ b/y2018/control_loops/python/spline_generate.py
@@ -0,0 +1,366 @@
+#!/usr/bin/python
+import numpy as np
+import matplotlib.pyplot as plt
+from frc971.control_loops.python import drivetrain
+
+# used to define properties of the drivetrain, changes depending on robot
+# see yXXXX/control_loops/python/drivetrain.py for the current values
+
+kDrivetrain = drivetrain.DrivetrainParams(
+ J = 6.0,
+ mass=68.0,
+ robot_radius=0.616 / 2.0,
+ wheel_radius=0.127 / 2.0 * 120.0 / 118.0,
+ G_low=46.0 / 60.0 * 20.0 / 48.0 * 14.0 / 62.0,
+ G_high=62.0 / 44.0 * 20.0 / 48.0 * 14.0 / 62.0,
+ q_pos_low=0.12,
+ q_pos_high=0.14,
+ q_vel_low=1.0,
+ q_vel_high=0.95,
+ efficiency=0.70,
+ has_imu=True,
+ force=True,
+ kf_q_voltage=13.0,
+ controller_poles=[0.82, 0.82],
+)
+
+drivetrain = drivetrain.Drivetrain(kDrivetrain)
+# set up coefficients for Hemrite basis function evaluation
+coeffs = np.array([[1, 0, 0, -10, 15, -6], [0, 1, 0, -6, 8, -3], [0, 0, 0.5, -1.5, 1.5, -0.5], [0, 0, 0, 0.5, -1, 0.5], [0, 0, 0, -4, 7, -3], [0, 0, 0, 10, -15, 6]])
+coeffs_prime = np.empty_like(coeffs)
+for ii in range(0, len(coeffs)):
+ for jj in range(0, len(coeffs[ii]) - 1):
+ coeffs_prime[ii][jj] = (jj + 1) * coeffs[ii][jj]
+
+def RungeKutta(f, x, dt):
+ """4th order RungeKutta integration of F starting at X."""
+ a = f(x)
+ b = f(x + dt / 2.0 * a)
+ c = f(x + dt / 2.0 * b)
+ d = f(x + dt * c)
+
+ return x + dt * (a + 2.0 * b + 2.0 * c + d) / 6.0
+
+def normalize(v):
+ norm = np.linalg.norm(v)
+ return v / norm
+
+def theta(v):
+ return np.arctan2(v[1], v[0])
+
+# evaluate Nth hermite basis function at t
+def nth_H(N, t):
+ return coeffs[N][0] + coeffs[N][1]*t + coeffs[N][2]*t**2 + coeffs[N][3]*t**3 + coeffs[N][4]*t**4 + coeffs[N][5]*t**5
+
+def nth_H_prime(N, t):
+ return coeffs[N][0] + coeffs[N][1]*t + coeffs[N][2]*t**2 + coeffs[N][3]*t**3 + coeffs[N][4]*t**4
+
+# class defining a quintic Hermite spline, with utilities for modification and plotting
+class Hermite_Spline:
+ # init method given known parameters, ie savefile loading(if necessary)
+ def __init__(self, start, control1, control2, end, resolution = 200):
+ self.start = start
+ self.end = end
+ self.control1 = control1
+ self.control2 = control2
+
+ self.points = np.array([])
+ self.velocities = []
+ self.accelerations = []
+ self.arc_lengths = []
+ self.thetas = []
+ self.omegas = []
+ self.curvatures = []
+
+ self.shifted_points = []
+
+ self.Ks = []
+ self.dKs = []
+
+ # coefficients are po, v0, a0, a1, v1, p1
+ self.coeffs = np.array([])
+ self.compute_coefficients()
+ self.resolution = resolution
+ self.setup()
+
+ # take paramters and compute coeffcicents for Hermite basis functions, to be called every time he change control points
+ def compute_coefficients(self):
+ self.coeffs = np.append(self.coeffs, np.array(self.start))
+ self.coeffs = np.append(self.coeffs, np.array(self.control1) - np.array(self.start))
+ self.coeffs = np.append(self.coeffs, [0,0])
+ self.coeffs = np.append(self.coeffs, [0,0])
+ self.coeffs = np.append(self.coeffs, np.array(self.end) - np.array(self.control2))
+ self.coeffs = np.append(self.coeffs, np.array(self.end))
+
+ self.coeffs = np.reshape(self.coeffs, newshape = (6, 2))
+
+ # setters for control points, set coefficients
+ def set_positions(self, p1 = None, p2 = None):
+ if p1 != None:
+ self.start = p1
+ if p2 != None:
+ self.end = p2
+
+ self.compute_coefficients()
+
+ def set_controls(self, c1 = None, c2 = None):
+ if c1 != None:
+ self.control1 = c1
+ if c2 != None:
+ self.control2 = c2
+
+ self.compute_coefficients()
+
+ def set_velocities(self, v1 = None, v2 = None):
+ if v1 != None:
+ self.control1 = self.start + v1
+ if v2 != None:
+ self.control2 = self.end + v2
+
+ self.compute_coefficients()
+
+ def get_smoothness(self):
+ K = self.get_curvature()
+ return np.sum(np.abs(np.gradient(K)))
+
+ # given Basis functions and controls compute coordinate given t
+ def spline_eval_hermite(self, t):
+ return np.array(self.coeffs[0]*nth_H(0, t) + self.coeffs[1]*nth_H(1, t)+ self.coeffs[2]*nth_H(2, t) + self.coeffs[3]*nth_H(3, t) + self.coeffs[4]* nth_H(4, t)+ self.coeffs[5]*nth_H(5, t))
+
+ # given Basis functions and controls compute velocity given t
+ def spline_eval_hermite_v(self, t):
+ return normalize(np.array(self.coeffs[0]*nth_H_prime(0, t) + self.coeffs[1]*nth_H_prime(1, t)+ self.coeffs[2]*nth_H_prime(2, t) + self.coeffs[3]*nth_H_prime(3, t) + self.coeffs[4]* nth_H_prime(4, t)+ self.coeffs[5]*nth_H_prime(5, t)))
+
+ # take coefficients and compute spline points/properties
+ def setup(self, resolution_multiplier = 10, dt = .000001):
+ points = []
+ velocities = []
+ accelerations = []
+ s = []
+ thetas = []
+ omegas = []
+ curvatures = []
+
+ last_point = self.spline_eval_hermite(0)
+ distance = 0
+
+ # iterate through interim points and compute pos_vectors, and at predefined points arc length,
+ # velocity, and acceleration vectors and store them at their associated index
+ for i in range(0, self.resolution * resolution_multiplier):
+ t = i / (1.0 * self.resolution * resolution_multiplier)
+
+ current_point = self.spline_eval_hermite(t)
+ current_point_dt = self.spline_eval_hermite(t + dt)
+ current_s = np.linalg.norm(current_point - last_point)
+
+ ds = np.linalg.norm(current_point_dt - current_point)
+
+ distance = current_s + distance
+ # at important points compute important values and store
+ if i % resolution_multiplier == 0:
+ s.append(distance)
+ points.append(current_point)
+
+ v = self.spline_eval_hermite_v(t)
+ v_dt = self.spline_eval_hermite_v(t + dt)
+ theta_t = theta(v)
+ theta_dt = theta(v_dt)
+
+ a = (v_dt - v) / ds
+ omega = (theta_dt - theta_t) / ds
+ if np.linalg.norm(v) == 0:
+ curvature = 0
+ else:
+ curvature = np.linalg.det(np.column_stack((v, a)) / (np.linalg.norm(v)**(3/2)))
+
+ velocities.append(v)
+ accelerations.append(a)
+ thetas.append(theta_t)
+ omegas.append(omega)
+ if curvature == 0:
+ curvatures.append(0.0001)
+ else:
+ curvatures.append(curvature)
+
+ last_point = current_point
+
+ self.arc_lengths = np.array(s)
+ self.points = np.reshape(points, newshape = (-1, 2))
+ self.velocities = np.reshape(velocities, newshape = (-1, 2))
+ self.accelerations = np.reshape(accelerations, newshape = (-1, 2))
+ self.thetas = np.array(thetas)
+ self.omegas = np.array(omegas)
+ self.curvatures = np.array(curvatures)
+
+
+ def plot_diagnostics(self):
+ plt.figure("Spline")
+ plt.title('Spline')
+ plt.plot(self.points[:, 0], self.points[:, 1])
+ # plt.scatter(self.points[:, 0], self.points[:, 1])
+
+ plt.figure("Diagnostics")
+
+ plt.subplot(2, 2, 1)
+ plt.title('theta')
+ plt.xlabel('arc_length')
+ plt.ylabel('theta')
+ theta, = plt.plot(self.arc_lengths, self.thetas, label = 'theta')
+ plt.legend(handles = [theta])
+
+ plt.subplot(2, 2, 2)
+ plt.title('omegas')
+ plt.xlabel('arc_length')
+ plt.ylabel('omega')
+ omega, = plt.plot(self.arc_lengths, self.omegas, label = 'omega')
+ plt.legend(handles = [omega])
+
+ plt.subplot(2, 2, 3)
+ plt.title('Velocities')
+ plt.xlabel('arc_length')
+ plt.ylabel('velocity')
+ dxds, = plt.plot(self.arc_lengths, self.velocities[:, 0], label = 'dx/ds')
+ dyds, = plt.plot(self.arc_lengths, self.velocities[:, 1], label = 'dy/ds')
+ plt.legend(handles = [dxds, dyds])
+
+ plt.subplot(2, 2, 4)
+ plt.title('Accelerations')
+ plt.xlabel('arc_length')
+ plt.ylabel('acceleration')
+ dx2ds2, = plt.plot(self.arc_lengths, self.accelerations[:, 0], label = 'd^2x/ds^2')
+ dy2ds2, = plt.plot(self.arc_lengths, self.accelerations[:, 1], label = 'd^2x/ds^2')
+ plt.legend(handles = [dx2ds2, dy2ds2])
+
+# class defining a number of splines with convinience methods
+class Path:
+ def __init__(self):
+ self.splines = []
+ self.knot_accels = []
+
+ def add_spline(self, spline):
+ self.splines.append(spline)
+
+ def get_K(self):
+ curvatures = []
+ for spline in self.splines:
+ curvatures.append(spline.curvatures)
+ return np.array(curvatures).flatten()
+
+ def get_S(self):
+ arc_lengths = []
+ for spline in self.splines:
+ arc_lengths.append(spline.arc_lengths)
+ return np.array(arc_lengths).flatten()
+
+ def get_points(self):
+ points = []
+ for spline in self.splines:
+ points.append(spline.points)
+ return points
+
+ def get_velocities(self, i):
+ velocities = []
+ for spline in self.splines:
+ velocities.append(spline.points)
+ return velocities
+
+ def remove_spine(self, i):
+ if i < len(self.splines):
+ self.splines.pop(i)
+ else:
+ print("index %f out of bounds, no spline of that index" % i)
+
+ def join(self, first_priority = False):
+ for i in range(0, len(self.splines)):
+ if first_priority & i != len(self.splines):
+ print("unfinished")
+
+
+# class which takes a Path object along with constraints and reparamterizes it with respect to time
+class Trajectory:
+ def __init__(self, path, max_angular_accel=3, max_voltage=11, max_normal_accel = .2):
+ self.path = path
+ self.A = drivetrain.A_continuous
+ self.B = drivetrain.B_continuous
+ self.robot_radius = drivetrain.robot_radius
+ self.Kv = 100
+ self.robot_radius = 3
+ self.max_angular_accel = max_angular_accel
+ self.max_voltage = max_voltage
+ self.max_normal_accel = max_normal_accel
+
+ self.max_velocities_adhering_to_normal_accel = []
+ self.max_velocities_adhering_to_voltage = []
+ self.path.splines[0].setup(resolution_multiplier = 100)
+
+ self.set_max_v_adhering_to_normal_accel()
+ self.max_voltageK_pass()
+
+ def set_max_v_adhering_to_normal_accel(self):
+ Ks = self.path.get_K()
+ accels = np.full_like(Ks, fill_value = self.max_normal_accel)
+ max_velocities = np.sqrt(np.abs(accels / Ks))
+ self.max_velocities_adhering_to_normal_accel = max_velocities
+
+ def max_voltageK_pass(self):
+ max_ds_dt = []
+ Ks = self.path.get_K()
+ turning_radii = np.full_like(Ks, fill_value = 1) / np.abs(Ks)
+
+
+
+ # compute max steady-state velocity given voltage constraints
+ for i in range(0, len(Ks)):
+ v_ratio = (turning_radii[i] + self.robot_radius) / (turning_radii[i] - self.robot_radius)
+ matrix = np.array([[self.A[1, 1], self.A[1, 3], self.B[1, 1]], [self.A[3, 1] - 1, self.A[3, 3], self.B[3, 1]], [-1, v_ratio, 0]])
+ sols = np.array([-1 * self.max_voltage * self.B[1, 0], -1 * self.max_voltage * self.B[3, 0], 0])
+ Vs = np.dot(np.linalg.inv(matrix), sols)
+ max_ds_dt.append((Vs[0] + Vs[1]) / 2)
+
+ self.max_velocities_adhering_to_voltage = max_ds_dt
+ # compute the maximum acceleration we can ask for given voltage and, ya know, staying on the path.
+
+
+ '''
+ These methods use the continuous form of our drivetrain state equation
+ in order to compute the maximum acceleration which adheres to the path
+ and voltage constraints, as well as any arbitary set of constraints
+ on velocity as a function of arc_length
+ '''
+
+ def forward_accel_pass(self):
+ points = self.path.get_points()
+ velocities = self.path.get_velocities()
+ curvatures = self.path.get_K()
+ arc_lenghts = self.path.get_S()
+
+ for i in range(0, len(points)):
+ Xn1 =
+
+
+ def backward_accelaration_pass(self):
+
+ print("max backward accel pass")
+
+
+ def plot_diagnostics(self, i = 0):
+
+ plt.figure('max velocity')
+ plt.title('max_v_normal_accel')
+ plt.xlabel('arc_length')
+ plt.ylabel('max V')
+ max_v_normal = plt.plot(self.path.get_S(), self.max_velocities_adhering_to_normal_accel, label = 'ds/dt (normal)')# , label = 'ds/dt')
+ curvature = plt.plot(self.path.get_S(), 1000 * np.abs(self.path.get_K()), label = 'K')
+ max_v_K_V = plt.plot(self.path.get_S(), self.max_velocities_adhering_to_voltage, label = 'ds/dt (voltage)')
+ plt.legend(handles = [max_v_normal[0], curvature[0], max_v_K_V[0]])
+
+def main():
+ A = Hermite_Spline(np.array([0,0]), np.array([0,400]), np.array([200,300]), np.array([200,200]), resolution = 200)
+ A.plot_diagnostics()
+ path = Path()
+ path.add_spline(A)
+ trajectory = Trajectory(path, 0)
+ trajectory.plot_diagnostics()
+ plt.show()
+
+main()