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

This is an example that uses the MocoInverse tool and EMG data to create an EMG-driven simulation of walking.

1 
2 import opensim as osim
3 import exampleEMGTracking_helpers as helpers
4 import os
5 import numpy as np
6 
7 
10 
11 # Part 1a: Load a 19 degree-of-freedom model with 18 lower-limb,
12 # sagittal-plane muscles and a torque-actuated torso. This includes a set of
13 # ground reaction forces applied to the model via ExternalLoads, which is
14 # necessary for the muscle redundancy problem. See the function definition
15 # at the bottom of this file to see how the model is loaded and constructed.
16 model = helpers.getWalkingModel()
17 
18 # Part 1b: Create the MocoInverse tool and set the Model.
19 inverse = osim.MocoInverse()
20 inverse.setModel(model)
21 
22 # Part 1c: Create a TableProcessor using the coordinates file from inverse
23 # kinematics.
24 coordinates = osim.TableProcessor('coordinates.mot')
25 coordinates.append(osim.TabOpLowPassFilter(6))
26 coordinates.append(osim.TabOpUseAbsoluteStateNames())
27 
28 # Part 1d: Set the kinematics reference for MocoInverse using the
29 # TableProcessor we just created.
30 inverse.setKinematics(coordinates)
31 inverse.set_kinematics_allow_extra_columns(True)
32 
33 # Part 1e: Provide the solver settings: initial and final time, the mesh
34 # interval, and the constraint and convergence tolerances.
35 inverse.set_initial_time(0.83)
36 inverse.set_final_time(2.0)
37 inverse.set_mesh_interval(0.04)
38 inverse.set_constraint_tolerance(1e-3)
39 inverse.set_convergence_tolerance(1e-3)
40 
41 if not os.path.isfile('effortSolution.sto'):
42  # Part 1f: Solve the problem!
43  inverseSolution = inverse.solve()
44  solution = inverseSolution.getMocoSolution()
45  solution.write('effortSolution.sto')
46 
47 
52 emgReference = osim.TimeSeriesTable('emg.sto')
53 helpers.compareSolutionToEMG(emgReference, 'effortSolution.sto')
54 
55 
58 
59 # Part 3a: Call initialize() to get access to the MocoStudy contained within
60 # the MocoInverse instance. This will allow us to make additional
61 # modifications to the problem not provided by MocoInverse.
62 study = inverse.initialize()
63 problem = study.updProblem()
64 
65 # Part 3b: Create a MocoControlTrackingGoal, set its weight, and provide
66 # the EMG data as the tracking reference. We also need to specify the
67 # reference labels for the four muscles whose EMG we will track.
68 tracking = osim.MocoControlTrackingGoal('emg_tracking')
69 tracking.setWeight(5)
70 tracking.setReference(osim.TableProcessor(emgReference))
71 tracking.setReferenceLabel('/forceset/gasmed_l', 'gastrocnemius')
72 tracking.setReferenceLabel('/forceset/tibant_l', 'tibialis_anterior')
73 tracking.setReferenceLabel('/forceset/bfsh_l', 'biceps_femoris')
74 tracking.setReferenceLabel('/forceset/glmax2_l', 'gluteus')
75 
76 # Part 3c: The EMG signals in the tracking are all normalized to have
77 # a maximum value of 1, but the magnitudes of the excitations from the
78 # effort minimization solution suggest that these signals should be
79 # rescaled. Use addScaleFactor() to add a MocoParameter to the problem that
80 # will scale the reference data for the muscles in the tracking cost.
81 tracking.addScaleFactor('gastroc_factor', '/forceset/gasmed_l', [0.01, 1.0])
82 tracking.addScaleFactor('tibant_factor', '/forceset/tibant_l', [0.01, 1.0])
83 tracking.addScaleFactor('bifem_factor', '/forceset/bfsh_l', [0.01, 1.0])
84 tracking.addScaleFactor('gluteus_factor', '/forceset/glmax2_l', [0.01, 1.0])
85 
86 # Part 3d: Add the tracking goal to the problem.
87 problem.addGoal(tracking)
88 
89 # Part 3e: Update the MocoCasADiSolver with the updated MocoProblem using
90 # resetProblem().
91 solver = osim.MocoCasADiSolver.safeDownCast(study.updSolver())
92 solver.resetProblem(problem)
93 
94 # Part 3f: Tell MocoCasADiSolver that the MocoParameters we added to the
95 # problem via addScaleFactor() above do not require initSystem() calls on
96 # the model. This provides a large speed-up.
97 solver.set_parameters_require_initsystem(False)
98 
99 if not os.path.isfile('trackingSolution.sto'):
100  # Part 3g: Solve the problem!
101  solution = study.solve()
102  solution.write('trackingSolution.sto')
103 
104 # Part 3h: Get the values of the optimized scale factors.
105 trackingSolution = osim.MocoTrajectory('trackingSolution.sto')
106 gastroc_factor = trackingSolution.getParameter('gastroc_factor')
107 tibant_factor = trackingSolution.getParameter('tibant_factor')
108 bifem_factor = trackingSolution.getParameter('bifem_factor')
109 gluteus_factor = trackingSolution.getParameter('gluteus_factor')
110 
111 
113 print('\nOptimized scale factor values:')
114 print('------------------------------')
115 print('gastrocnemius = ' + str(gastroc_factor))
116 print('tibialis anterior = ' + str(tibant_factor))
117 print('biceps femoris short head = ' + str(bifem_factor))
118 print('gluteus = ' + str(gluteus_factor))
119 
120 # Part 4b: Re-scale the reference data using the optimized scale factors.
121 gastroc = emgReference.updDependentColumn('gastrocnemius')
122 tibant = emgReference.updDependentColumn('tibialis_anterior')
123 bifem = emgReference.updDependentColumn('biceps_femoris')
124 gluteus = emgReference.updDependentColumn('gluteus')
125 for t in np.arange(emgReference.getNumRows()):
126  t = int(t) # Convert to Python built-in int type for indexing
127  gastroc[t] = gastroc_factor * gastroc[t]
128  tibant[t] = tibant_factor * tibant[t]
129  bifem[t] = bifem_factor * bifem[t]
130  gluteus[t] = gluteus_factor * gluteus[t]
131 
132 # Part 4c: Generate the plots. Compare results to the effort minimization
133 # solution.
134 helpers.compareSolutionToEMG(emgReference, 'effortSolution.sto',
135  'trackingSolution.sto')