In [1]:

    

import sys
sys.path.insert(0, '..')
import time
import matplotlib.pyplot as plt
%matplotlib notebook 
from IPython.display import set_matplotlib_formats
set_matplotlib_formats('pdf', 'png')
plt.rcParams['savefig.dpi'] = 75

plt.rcParams['figure.autolayout'] = False
plt.rcParams['figure.figsize'] = 10, 6
plt.rcParams['axes.labelsize'] = 18
plt.rcParams['axes.titlesize'] = 20
plt.rcParams['font.size'] = 16
plt.rcParams['lines.linewidth'] = 2.0
plt.rcParams['lines.markersize'] = 8
plt.rcParams['legend.fontsize'] = 14

plt.rcParams['text.usetex'] = True
plt.rcParams['font.family'] = "serif"
plt.rcParams['font.serif'] = "cm"
plt.rcParams['text.latex.preamble'] = "\usepackage{subdepth}, \usepackage{type1cm}"


import numpy as np

from Configuration import Configuration
from MotorUnitPool import MotorUnitPool
from NeuralTract import NeuralTract
from AfferentPool import AfferentPool
from SynapsesFactory import SynapsesFactory
from jointAnkleForceTask import jointAnkleForceTask
from numba import jit, prange
import scipy as sc
from scipy.signal import *


    

In [2]:

    

conf = Configuration('confForceVariability.rmto')
conf.simDuration_ms = 5000 # Here I change simulation duration without changing the Configuration file.
t = np.arange(0.0, conf.simDuration_ms, conf.timeStep_ms)


    

In [3]:

    

FR = np.array([95.0, 110.0 ,125.0, 140.0, 155.0, 170.0, 185.0, 200.0, 250])
GammaOrder = np.array([7, 5, 4, 4, 4, 3, 2, 2, 1])
condVelS = np.linspace(44/2.0, 44*2.0, 5)
condVelFR = np.linspace(51/2.0, 51*2.0, 5)
condVelFF1 = np.linspace(52/2.0, 52*2.0, 5)
condVelFF2 = np.linspace(53/2.0,53*2, 5)


forceSOL = np.empty((len(t),len(FR), len(condVelS)))
forceMG = np.empty((len(t),len(FR), len(condVelS)))
forceLG = np.empty((len(t),len(FR), len(condVelS)))
torque = np.empty((len(t), len(FR), len(condVelS)))


for m in xrange(1,2):
    
    conf.changeConfigurationParameter('axonDelayCondVel:MG-S',condVelS[m],condVelFR[m])
    conf.changeConfigurationParameter('axonDelayCondVel:SOL-S',condVelS[m],condVelFR[m])
    conf.changeConfigurationParameter('axonDelayCondVel:LG-S',condVelS[m],condVelFR[m])
    conf.changeConfigurationParameter('axonDelayCondVel:MG-FR',condVelFR[m],condVelFF1[m])
    conf.changeConfigurationParameter('axonDelayCondVel:SOL-FR',condVelFR[m],condVelFF1[m])
    conf.changeConfigurationParameter('axonDelayCondVel:LG-FR',condVelFR[m],condVelFF1[m])
    conf.changeConfigurationParameter('axonDelayCondVel:MG-FF',condVelFF1[m],condVelFF2[m])
    conf.changeConfigurationParameter('axonDelayCondVel:SOL-FF',condVelFF1[m],condVelFF2[m])
    conf.changeConfigurationParameter('axonDelayCondVel:LG-FF',condVelFF1[m],condVelFF2[m])
    
    
    
    pools = dict()
    pools[0] = MotorUnitPool(conf, 'SOL')
    pools[1] = MotorUnitPool(conf, 'MG')
    pools[2] = MotorUnitPool(conf, 'LG')
    pools[3] = NeuralTract(conf, 'CMExt')
    ankle = jointAnkleForceTask(conf, pools)
    Syn = SynapsesFactory(conf, pools)
    del Syn



    for l in xrange(3,4):
        print('Firing Rate ' + str(FR[l]))
        tic = time.time()
        for i in xrange(0,len(t)-1): 
            pools[3].atualizePool(t[i], FR[l], GammaOrder[l])
            ankle.atualizeAnkle(t[i], 0.0)
            for j in xrange(3):
                pools[j].atualizeMotorUnitPool(t[i])
            ankle.computeTorque(t[i])

        toc = time.time() 


        forceSOL[:,l,m] = np.squeeze(pools[0].Muscle.force)
        forceMG[:,l,m] = np.squeeze(pools[1].Muscle.force)
        forceLG[:,l,m] = np.squeeze(pools[2].Muscle.force)
        torque[:,l,m] = np.squeeze(ankle.ankleTorque_Nm)


        print str(toc - tic) + ' seconds'
        for k in xrange(0, len(pools)):
            pools[k].reset()


    

    




Muscle spindle from muscle SOL built.
Motor Unit Pool SOL built
Muscle spindle from muscle MG built.
Motor Unit Pool MG built
Muscle spindle from muscle LG built.
Motor Unit Pool LG built
Descending Command CMExt built
Ankle joint for Force Task built
All the 210959 synapses were built
All the 0 synaptic noises were built
Firing Rate 140.0
18042.3814211 seconds

In [4]:

    

forceLG.shape

plt.figure()
plt.plot(t, torque[:,:,1])


    

    















    












    
Out[4]:





[<matplotlib.lines.Line2D at 0x7f4237175e10>,
 <matplotlib.lines.Line2D at 0x7f4237175f50>,
 <matplotlib.lines.Line2D at 0x7f4235c71090>,
 <matplotlib.lines.Line2D at 0x7f4235c71190>,
 <matplotlib.lines.Line2D at 0x7f4235c71290>,
 <matplotlib.lines.Line2D at 0x7f4235c71390>,
 <matplotlib.lines.Line2D at 0x7f4235c71490>,
 <matplotlib.lines.Line2D at 0x7f4235c71590>,
 <matplotlib.lines.Line2D at 0x7f4235c71690>]

In [5]:

    

plt.figure()
plt.plot(t, forceSOL)


    

    















    












    




---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-5-51d14c67d9da> in <module>()
      1 plt.figure()
----> 2 plt.plot(t, forceSOL)

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/pyplot.pyc in plot(*args, **kwargs)
   3315                       mplDeprecation)
   3316     try:
-> 3317         ret = ax.plot(*args, **kwargs)
   3318     finally:
   3319         ax._hold = washold

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/__init__.pyc in inner(ax, *args, **kwargs)
   1896                     warnings.warn(msg % (label_namer, func.__name__),
   1897                                   RuntimeWarning, stacklevel=2)
-> 1898             return func(ax, *args, **kwargs)
   1899         pre_doc = inner.__doc__
   1900         if pre_doc is None:

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/axes/_axes.pyc in plot(self, *args, **kwargs)
   1404         kwargs = cbook.normalize_kwargs(kwargs, _alias_map)
   1405 
-> 1406         for line in self._get_lines(*args, **kwargs):
   1407             self.add_line(line)
   1408             lines.append(line)

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/axes/_base.pyc in _grab_next_args(self, *args, **kwargs)
    405                 return
    406             if len(remaining) <= 3:
--> 407                 for seg in self._plot_args(remaining, kwargs):
    408                     yield seg
    409                 return

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/axes/_base.pyc in _plot_args(self, tup, kwargs)
    383             x, y = index_of(tup[-1])
    384 
--> 385         x, y = self._xy_from_xy(x, y)
    386 
    387         if self.command == 'plot':

/opt/miniconda2/lib/python2.7/site-packages/matplotlib/axes/_base.pyc in _xy_from_xy(self, x, y)
    245         if x.ndim > 2 or y.ndim > 2:
    246             raise ValueError("x and y can be no greater than 2-D, but have "
--> 247                              "shapes {} and {}".format(x.shape, y.shape))
    248 
    249         if x.ndim == 1:

ValueError: x and y can be no greater than 2-D, but have shapes (100000,) and (100000, 9, 5)

In [ ]:

    

plt.figure()
plt.plot(t, forceMG)


    

In [ ]:

    

plt.figure()
plt.plot(t, forceLG)


    

In [ ]:

    

np.save('forceSol',forceSOL,)
np.save('forceMG', forceMG)
np.save('forceLG', forceLG)
np.save('torque', torque)


    

In [ ]:

    

X=np.load('forceSol.npy')


    

In [ ]:

    




    

In [ ]: