From e9f4f691e554000f9decf84da56568d6150badb2 Mon Sep 17 00:00:00 2001
From: Tobias Klauser <tklauser@distanz.ch>
Date: Thu, 22 Nov 2012 19:22:42 +0100
Subject: Initial import (based on Leg6DoF9Muscles model)

---
 LocomotorPrimitives.cpp | 113 ++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 113 insertions(+)
 create mode 100644 LocomotorPrimitives.cpp

(limited to 'LocomotorPrimitives.cpp')

diff --git a/LocomotorPrimitives.cpp b/LocomotorPrimitives.cpp
new file mode 100644
index 0000000..7ceecb4
--- /dev/null
+++ b/LocomotorPrimitives.cpp
@@ -0,0 +1,113 @@
+/**
+ * Simulation for the Locomotor Primitives (Dominici et. al 2012) in a Human
+ * Leg model with 14 Muscles.
+ *
+ * Author: Tobias Klauser <tobias.klauser@uzh.ch>
+ */
+
+#include <OpenSim/OpenSim.h>
+
+#include "LocomotorPrimitivesManager.h"
+#include "LocomotorPrimitivesController.h"
+
+#define NO_SIM		0
+
+static const std::string MODEL_NAME = "LocomotorPrimitives";
+
+// Define the initial and final simulation time
+static const double initialTime = 0.0;
+static const double finalTime = 5.0;
+
+static void constructModel(OpenSim::Model &model)
+{
+	std::cout << "Constructing " << MODEL_NAME << " model" << std::endl;
+
+	// Define controller for the model
+	LocomotorPrimitivesController *control = new LocomotorPrimitivesController(0.01);
+	control->setActuators(model.updActuators());
+	model.addController(control);
+
+	// Set default activation and fiber length on all muscles of the model
+	const OpenSim::Set<OpenSim::Actuator> &actSet = model.getActuators();
+	int sz = actSet.getSize();
+	std::cout << " + Defining initial muscle states for " << sz << " actuators:" << std::endl;
+
+	// Set default activation and fiber length of all muscles
+	for (int i = 0; i < sz; i++) {
+		OpenSim::ActivationFiberLengthMuscle *m = dynamic_cast<OpenSim::ActivationFiberLengthMuscle *>(&actSet.get(i));
+		std::cout << "   " << m->getName() << std::endl;
+		m->setDefaultActivation(0.01);
+		m->setDefaultFiberLength(0.1);
+	}
+
+	// Create the force reporter for obtaining the forces applied to the model
+	// during a forward simulation
+	OpenSim::ForceReporter *reporter = new OpenSim::ForceReporter(&model);
+	model.addAnalysis(reporter);
+
+	std::cout << "Finished constructing model" << std::endl;
+}
+
+void simulateModel(OpenSim::Model &model)
+{
+	if (NO_SIM) {
+		std::cout << "Skipping simulation as per NO_SIM=" << NO_SIM << std::endl;
+		return;
+	}
+
+	std::cout << "Simulating model " << MODEL_NAME << std::endl;
+	std::cout << " + Initializing system" << std::endl;
+	SimTK::State &si = model.initSystem();
+
+	std::cout << " + Creating integrator" << std::endl;
+	// Create the integrator and manager for the simulation
+	SimTK::RungeKuttaMersonIntegrator integrator(model.getMultibodySystem());
+	integrator.setAccuracy(1.0e-4);
+	std::cout << " + Creating simulation manager" << std::endl;
+	
+	LocomotorPrimitivesManager manager(model, integrator);
+
+	// Integrate from initial time to final time
+	manager.setInitialTime(initialTime);
+	manager.setFinalTime(finalTime);
+	std::cout << " + Integrating from " << std::fixed << initialTime << " to " << std::fixed << finalTime << std::endl;
+	manager.integrate(si);
+
+	// Save the model states from forward integration
+	OpenSim::Storage statesDegrees(manager.getStateStorage());
+	statesDegrees.print(MODEL_NAME + "_states.sto");
+
+	// Save the forces
+	OpenSim::Analysis &reporter = model.getAnalysisSet()[0];
+	((OpenSim::ForceReporter &)reporter).getForceStorage().print(MODEL_NAME + "_forces.mot");
+}
+
+int main(void)
+{
+	clock_t ts_start = clock();
+
+	try {
+		// Create an OpenSim model and set its name
+		OpenSim::Model osimModel("../../leg6dof9musc.osim");
+		osimModel.setName(MODEL_NAME);
+		constructModel(osimModel);
+
+		osimModel.setUseVisualizer(false);
+		simulateModel(osimModel);
+	} catch (OpenSim::Exception ex) {
+		std::cout << ex.getMessage() << std::endl;
+		return 1;
+	} catch (std::exception ex) {
+		std::cout << ex.what() << std::endl;
+		return 1;
+	} catch (...) {
+		std::cout << "UNRECOGNIZED EXCEPTION" << std::endl;
+		return 1;
+	}
+
+	std::cout << "main() routine time = " << 1.e3*(clock()-ts_start)/CLOCKS_PER_SEC << "ms" << std::endl;
+
+    std::cout << "Simulation " + MODEL_NAME << " completed successfully." << std::endl;
+	std::cin.get();
+	return 0;
+}
\ No newline at end of file
-- 
cgit v1.2.3-54-g00ecf