#!/bin/env python
# -*- coding: utf-8 -*-

from brian2 import *
import random
import sys

if(len(sys.argv) != 2):
    print("Usage : ", sys.argv[0] , " shuffle(0,1) ")
    sys.exit(-1)

shuffle = int(sys.argv[1])

tau = 20 * msecond        # constant de temps 
taud = 2 * msecond
Vt = -55 * mvolt          # seuil pour spiker
Vr = -75 * mvolt          # potentiel de reset
Vs = -65 * mvolt          # potentiel de repos
dt = defaultclock.dt

g_weight = 0.07 * (Vt - Vs)
print("Poids du couplage" ,g_weight)
I0 = 40 * mV
t_init = [20*ms, 70* ms,  125 * ms]

# On fixe la décroissance du premier stimulus à -0.05 mvolt/msecond
alpha = [-0.05*mvolt/msecond] +  ([0] * (len(t_init)-1))
t_cross = 200 * msecond
for i,(t1,a1) in enumerate(zip (t_init[1:],alpha[1:])):
    alpha[i+1] = alpha[0] * (t_cross - t_init[0])/(t_cross - t1)

# On shuffle la liste des temps
if(shuffle):
    random.shuffle(t_init)
print(t_init)

tmax = 800*msecond


# Notre groupe contient 3 neurones intègre et tire 
eqs = '''
dV/dt = (-(V-Vs) + I)/tau : volt
I : volt
'''

eqsd = '''
dV/dt = -(V-Vs)/taud : volt
I : volt
'''

N = len(t_init)
Ga = NeuronGroup(N=N, model=eqs,
                 threshold='V > Vt', reset='V = Vr', method='linear')
Si = Synapses(Ga, Ga, 'w : volt', on_pre='V_post += w')
Si.connect(condition='i!=j')
Si.w = g_weight


# Notre neurone détecteur de synchronie
Gd = NeuronGroup(N=1, model=eqsd,
                 threshold='V > Vt', reset='V = Vr', method='linear')
Sd = Synapses(Ga, Gd, 'w : volt', on_pre='V_post += w')
Sd.connect()
Sd.w = 4.5*mvolt


# On monitore ses spikes et son potentiel
spikemon = SpikeMonitor(Ga)
spikemond = SpikeMonitor(Gd)
statemon_input = StateMonitor(Ga, 'I', record=True)

# On initialise les variables des neurones
Ga.V = Vs
Ga.I = 0


@check_units(t=second, result=volt)
def get_currents(t):
    I = I0 + array(alpha)*volt/second*(t - array(t_init)*second)
    I[where(I < 0)] = 0
    I[where(I > I0)] = 0
    return I

# Pour injecter le courant, on se définit
# une fonction appelée à chaque itération
@network_operation(when='start')
def myoperation():
    global Ga
    # On modifie les entrées du neurones
    Ga.I = get_currents(defaultclock.t)

# On lance la simulation
run(tmax)

figure(figsize=(10,10))
subplot(211)
for i, ts in zip(spikemon.i, spikemon.t):
    plot([ts/msecond, ts/msecond], [i, i+1], 'k')

for ts in spikemond.spike_trains()[0]:
    plot([ts/msecond, ts/msecond], [len(Ga),len(Ga)+1], 'r')

xlim([0, tmax/msecond])
ylim([-1, N+2])
title("Trains de spikes")

subplot(212)
plot(statemon_input.t, statemon_input.I.T)
xlim([0, tmax/second])
title(u'Courants d\'entrée')


if(shuffle):
    savefig('tout_connecte_shuffle.png')
else:
    savefig('tout_connecte.png')

show()