|
2 | 2 | # -*- coding: utf-8 -*-
|
3 | 3 |
|
4 | 4 | import numpy as np
|
| 5 | +import matplotlib.pyplot as plt |
5 | 6 |
|
6 | 7 | from dtk import process
|
7 | 8 |
|
@@ -55,3 +56,126 @@ def find_constant_speed(time, speed, plot=False):
|
55 | 56 | fig.show()
|
56 | 57 |
|
57 | 58 | return len(time) - (new_indice), time[len(time) - new_indice]
|
| 59 | + |
| 60 | + |
| 61 | +def gait_landmarks_from_grf(mot_file, |
| 62 | + right_grfy_column_name='ground_force_vy', |
| 63 | + left_grfy_column_name='1_ground_force_vy', |
| 64 | + threshold=1e-5, do_plot=False, min_time=None, |
| 65 | + max_time=None): |
| 66 | + """ |
| 67 | + Obtain gait landmarks (right and left foot strike & toe-off) from ground |
| 68 | + reaction force (GRF) time series data. |
| 69 | +
|
| 70 | + Parameters |
| 71 | + ---------- |
| 72 | + mot_file : str |
| 73 | + Name of *.mot (OpenSim Storage) file containing GRF data. |
| 74 | + right_grfy_column_name : str, optional |
| 75 | + Name of column in `mot_file` containing the y (vertical) component |
| 76 | + of GRF data for the right leg. |
| 77 | + left_grfy_column_name : str, optional |
| 78 | + Same as above, but for the left leg. |
| 79 | + threshold : float, optional |
| 80 | + Below this value, the force is considered to be zero (and the |
| 81 | + corresponding foot is not touching the ground). |
| 82 | + do_plot : bool, optional (default: False) |
| 83 | + Create plots of the detected gait landmarks on top of the vertical |
| 84 | + ground reaction forces. |
| 85 | + min_time : float, optional |
| 86 | + If set, only consider times greater than `min_time`. |
| 87 | + max_time : float, optional |
| 88 | + If set, only consider times greater than `max_time`. |
| 89 | +
|
| 90 | + Returns |
| 91 | + ------- |
| 92 | + right_foot_strikes : np.array |
| 93 | + All times at which right_grfy is non-zero and it was 0 at the |
| 94 | + preceding time index. |
| 95 | + left_foot_strikes : np.array |
| 96 | + Same as above, but for the left foot. |
| 97 | + right_toe_offs : np.array |
| 98 | + All times at which left_grfy is 0 and it was non-zero at the |
| 99 | + preceding time index. |
| 100 | + left_toe_offs : np.array |
| 101 | + Same as above, but for the left foot. |
| 102 | +
|
| 103 | + Notes |
| 104 | + ----- |
| 105 | + Source modifed from: |
| 106 | +
|
| 107 | + https://github.com/fitze/epimysium/blob/master/epimysium/postprocessing.py |
| 108 | +
|
| 109 | + """ |
| 110 | + data = dataman.storage2numpy(mot_file) |
| 111 | + |
| 112 | + time = data['time'] |
| 113 | + right_grfy = data[right_grfy_column_name] |
| 114 | + left_grfy = data[left_grfy_column_name] |
| 115 | + |
| 116 | + # Time range to consider. |
| 117 | + if max_time == None: max_idx = len(time) |
| 118 | + else: max_idx = nearest_index(time, max_time) |
| 119 | + |
| 120 | + if min_time == None: min_idx = 1 |
| 121 | + else: min_idx = max(1, nearest_index(time, min_time)) |
| 122 | + |
| 123 | + index_range = range(min_idx, max_idx) |
| 124 | + |
| 125 | + # Helper functions |
| 126 | + # ---------------- |
| 127 | + def zero(number): |
| 128 | + return abs(number) < threshold |
| 129 | + |
| 130 | + def birth_times(ordinate): |
| 131 | + births = list() |
| 132 | + for i in index_range: |
| 133 | + # 'Skip' first value because we're going to peak back at previous |
| 134 | + # index. |
| 135 | + if zero(ordinate[i - 1]) and (not zero(ordinate[i])): |
| 136 | + births.append(time[i]) |
| 137 | + return np.array(births) |
| 138 | + |
| 139 | + def death_times(ordinate): |
| 140 | + deaths = list() |
| 141 | + for i in index_range: |
| 142 | + if (not zero(ordinate[i - 1])) and zero(ordinate[i]): |
| 143 | + deaths.append(time[i]) |
| 144 | + return np.array(deaths) |
| 145 | + |
| 146 | + right_foot_strikes = birth_times(right_grfy) |
| 147 | + left_foot_strikes = birth_times(left_grfy) |
| 148 | + right_toe_offs = death_times(right_grfy) |
| 149 | + left_toe_offs = death_times(left_grfy) |
| 150 | + |
| 151 | + if do_plot: |
| 152 | + |
| 153 | + #pl.figure(figsize=(4, 8)) |
| 154 | + ones = np.array([1, 1]) |
| 155 | + |
| 156 | + def myplot(index, label, ordinate, foot_strikes, toe_offs): |
| 157 | + ax = plt.subplot(2, 1, index) |
| 158 | + plt.plot(time[min_idx:max_idx], ordinate[min_idx:max_idx], 'k') |
| 159 | + plt.ylabel('vertical ground reaction force (N)') |
| 160 | + plt.title('%s (%i foot strikes, %i toe-offs)' % ( |
| 161 | + label, len(foot_strikes), len(toe_offs))) |
| 162 | + |
| 163 | + for i, strike in enumerate(foot_strikes): |
| 164 | + if i == 0: kwargs = {'label': 'foot strikes'} |
| 165 | + else: kwargs = dict() |
| 166 | + plt.plot(strike * ones, ax.get_ylim(), 'r', **kwargs) |
| 167 | + |
| 168 | + for i, off in enumerate(toe_offs): |
| 169 | + if i == 0: kwargs = {'label': 'toe-offs'} |
| 170 | + else: kwargs = dict() |
| 171 | + plt.plot(off * ones, ax.get_ylim(), 'b', **kwargs) |
| 172 | + |
| 173 | + |
| 174 | + myplot(1, 'left foot', left_grfy, left_foot_strikes, left_toe_offs) |
| 175 | + plt.legend(loc='best') |
| 176 | + |
| 177 | + myplot(2, 'right foot', right_grfy, right_foot_strikes, right_toe_offs) |
| 178 | + |
| 179 | + plt.xlabel('time (s)') |
| 180 | + |
| 181 | + return right_foot_strikes, left_foot_strikes, right_toe_offs, left_toe_offs |
0 commit comments