API  4.3
For MATLAB, Python, Java, and C++ users
exampleMocoInverse.py

This is an example using the MocoInverse tool with a complex model to prescribe walking. This example also shows how to track electromyography data.

1 # -------------------------------------------------------------------------- #
2 # OpenSim Moco: exampleMocoInverse.py #
3 # -------------------------------------------------------------------------- #
4 # Copyright (c) 2020 Stanford University and the Authors #
5 # #
6 # Author(s): Christopher Dembia #
7 # #
8 # Licensed under the Apache License, Version 2.0 (the "License") you may #
9 # not use this file except in compliance with the License. You may obtain a #
10 # copy of the License at http://www.apache.org/licenses/LICENSE-2.0 #
11 # #
12 # Unless required by applicable law or agreed to in writing, software #
13 # distributed under the License is distributed on an "AS IS" BASIS, #
14 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. #
15 # See the License for the specific language governing permissions and #
16 # limitations under the License. #
17 # -------------------------------------------------------------------------- #
18 
19 # This example shows how to use the MocoInverse tool to exactly prescribe a
20 # motion and estimate muscle behavior for walking. The first example does not
21 # rely on electromyography data, while the second example penalizes deviation
22 # from electromyography data for a subset of muscles.
23 #
24 # See the README.txt next to this file for more information.
25 
26 import opensim as osim
27 
28 def solveMocoInverse():
29 
30  # Construct the MocoInverse tool.
31  inverse = osim.MocoInverse()
32 
33  # Construct a ModelProcessor and set it on the tool. The default
34  # muscles in the model are replaced with optimization-friendly
35  # DeGrooteFregly2016Muscles, and adjustments are made to the default muscle
36  # parameters.
37  modelProcessor = osim.ModelProcessor('subject_walk_armless.osim')
38  modelProcessor.append(osim.ModOpAddExternalLoads('grf_walk.xml'))
39  modelProcessor.append(osim.ModOpIgnoreTendonCompliance())
40  modelProcessor.append(osim.ModOpReplaceMusclesWithDeGrooteFregly2016())
41  # Only valid for DeGrooteFregly2016Muscles.
42  modelProcessor.append(osim.ModOpIgnorePassiveFiberForcesDGF())
43  # Only valid for DeGrooteFregly2016Muscles.
44  modelProcessor.append(osim.ModOpScaleActiveFiberForceCurveWidthDGF(1.5))
45  modelProcessor.append(osim.ModOpAddReserves(1.0))
46  inverse.setModel(modelProcessor)
47 
48  # Construct a TableProcessor of the coordinate data and pass it to the
49  # inverse tool. TableProcessors can be used in the same way as
50  # ModelProcessors by appending TableOperators to modify the base table.
51  # A TableProcessor with no operators, as we have here, simply returns the
52  # base table.
53  inverse.setKinematics(osim.TableProcessor('coordinates.sto'))
54 
55  # Initial time, final time, and mesh interval.
56  inverse.set_initial_time(0.81)
57  inverse.set_final_time(1.79)
58  inverse.set_mesh_interval(0.02)
59 
60  # By default, Moco gives an error if the kinematics contains extra columns.
61  # Here, we tell Moco to allow (and ignore) those extra columns.
62  inverse.set_kinematics_allow_extra_columns(True)
63 
64  # Solve the problem and write the solution to a Storage file.
65  solution = inverse.solve()
66  solution.getMocoSolution().write('example3DWalking_MocoInverse_solution.sto')
67 
68  # Generate a PDF with plots for the solution trajectory.
69  model = modelProcessor.process()
70  report = osim.report.Report(model,
71  'example3DWalking_MocoInverse_solution.sto',
72  bilateral=True)
73  # The PDF is saved to the working directory.
74  report.generate()
75 
76 def solveMocoInverseWithEMG():
77 
78  # This initial block of code is identical to the code above.
79  inverse = osim.MocoInverse()
80  modelProcessor = osim.ModelProcessor('subject_walk_armless.osim')
81  modelProcessor.append(osim.ModOpAddExternalLoads('grf_walk.xml'))
82  modelProcessor.append(osim.ModOpIgnoreTendonCompliance())
83  modelProcessor.append(osim.ModOpReplaceMusclesWithDeGrooteFregly2016())
84  modelProcessor.append(osim.ModOpIgnorePassiveFiberForcesDGF())
85  modelProcessor.append(osim.ModOpScaleActiveFiberForceCurveWidthDGF(1.5))
86  modelProcessor.append(osim.ModOpAddReserves(1.0))
87  inverse.setModel(modelProcessor)
88  inverse.setKinematics(osim.TableProcessor('coordinates.sto'))
89  inverse.set_initial_time(0.81)
90  inverse.set_final_time(1.79)
91  inverse.set_mesh_interval(0.02)
92  inverse.set_kinematics_allow_extra_columns(True)
93 
94  study = inverse.initialize()
95  problem = study.updProblem()
96 
97  # Add electromyography tracking.
98  emgTracking = osim.MocoControlTrackingGoal('emg_tracking')
99  emgTracking.setWeight(50.0)
100  # Each column in electromyography.sto is normalized so the maximum value in
101  # each column is 1.0.
102  controlsRef = osim.TimeSeriesTable('electromyography.sto')
103  # Scale the tracked muscle activity based on peak levels from
104  # "Gait Analysis: Normal and Pathological Function" by
105  # Perry and Burnfield, 2010 (digitized by Carmichael Ong).
106  soleus = controlsRef.updDependentColumn('soleus')
107  gasmed = controlsRef.updDependentColumn('gastrocnemius')
108  tibant = controlsRef.updDependentColumn('tibialis_anterior')
109  for t in range(0, controlsRef.getNumRows()):
110  soleus[t] = 0.77 * soleus[t]
111  gasmed[t] = 0.87 * gasmed[t]
112  tibant[t] = 0.37 * tibant[t]
113  emgTracking.setReference(osim.TableProcessor(controlsRef))
114  # Associate actuators in the model with columns in electromyography.sto.
115  emgTracking.setReferenceLabel('/forceset/soleus_r', 'soleus')
116  emgTracking.setReferenceLabel('/forceset/gasmed_r', 'gastrocnemius')
117  emgTracking.setReferenceLabel('/forceset/gaslat_r', 'gastrocnemius')
118  emgTracking.setReferenceLabel('/forceset/tibant_r', 'tibialis_anterior')
119  problem.addGoal(emgTracking)
120 
121  # Solve the problem and write the solution to a Storage file.
122  solution = study.solve()
123  solution.write('example3DWalking_MocoInverseWithEMG_solution.sto')
124 
125  # Write the reference data in a way that's easy to compare to the solution.
126  controlsRef.removeColumn('medial_hamstrings')
127  controlsRef.removeColumn('biceps_femoris')
128  controlsRef.removeColumn('vastus_lateralis')
129  controlsRef.removeColumn('vastus_medius')
130  controlsRef.removeColumn('rectus_femoris')
131  controlsRef.removeColumn('gluteus_maximus')
132  controlsRef.removeColumn('gluteus_medius')
133  controlsRef.setColumnLabels(['/forceset/soleus_r', '/forceset/gasmed_r',
134  '/forceset/tibant_r'])
135  controlsRef.appendColumn('/forceset/gaslat_r', gasmed)
136  osim.STOFileAdapter.write(controlsRef, 'controls_reference.sto')
137 
138  # Generate a report comparing MocoInverse solutions without and with EMG
139  # tracking.
140  model = modelProcessor.process()
141  output = 'example3DWalking_MocoInverseWithEMG_report.pdf'
142  ref_files = [
143  'example3DWalking_MocoInverseWithEMG_solution.sto',
144  'controls_reference.sto']
145  report = osim.report.Report(model,
146  'example3DWalking_MocoInverse_solution.sto',
147  output=output, bilateral=True,
148  ref_files=ref_files)
149  # The PDF is saved to the working directory.
150  report.generate()
151 
152 
153 # This problem solves in about 5 minutes.
154 solveMocoInverse()
155 
156 # This problem penalizes the deviation from electromyography data for a
157 # subset of muscles, and solves in about 30 minutes.
158 solveMocoInverseWithEMG()