""" 
  This file is licensed under the Allen Institute Terms of Use, which can be found at:
http://www.alleninstitute.org/Media/policies/terms_of_use_content.html
--------------------------------------------------------------------------

Script to generate NEST data for Fig.1 (used in panel h)) 

NOTE: Data files originally generated and used for Fig. have been provided

Run as : 

$ python Fig1h_EIFCond_generate_NEST_data.py <simulator>

where <simulator> is 'neuron', 'nest', 'brian', etc
 
"""
###################################################################################
# Preamble
###################################################################################


import matplotlib as mpl
mpl.use('Agg')

# Import necessary packages:
from pyNN.utility import Timer

# Get simulator; probably use "nest"
from pyNN.nest import *

timer = Timer()

###################################################################################
# Settings and Parameters
###################################################################################

# Settings:
simtime     = 205.0     # simulation time [ms]
NE          = 10000     # number of exc neurons
dt          = 0.05      # simulation step length [ms]
FR_bg = 100             #firing rate in Hz

# single neuron parameters
tauMem      = 20.0    # neuron membrane time constant [ms]
tauSynE     = dt      # synaptic time constant [ms]
tauRef      = 0.0     # refractory time [ms]
U0          = 0.0     # resting potential [mV]
theta       = 20.0    # threshold potential [mV]
#firstbinedge = 5.086142476173294e-03
firstbinedge = 5.094321737688319e-02

#Synaptic Weight and Rate#
J = 1.0 #peak synaptic conductance(uS)

# put cell parameters into a dict
neuronClass = EIF_cond_exp_isfa_ista
cell_params = {'tau_m'      : tauMem,
               'tau_syn_E'  : tauSynE,
               'tau_refrac' : tauRef,
               'v_rest'     : U0,
               'v_reset'    : U0,
               'v_thresh'   : 19.0,    #(mV and called V_T in the paper)
               'v_spike'    : theta, 
               'a'          : 0.0,
               'b'          : 0.0,
               'delta_T'    : 2.0,     #(mV) slope factor for EIF neurons
               'tau_w'      : dt,
               'e_rev_E'    : 100.0,   #(mV)
               'cm'         : 1.0}     # (nF)

###################################################################################
# Create Network
###################################################################################

extra = {'threads' : 8}
rank = setup(timestep = dt, min_delay = dt, spike_precision = 'off_grid', **extra) 

# Start construction timer:
timer.start() 

# Create a target population and initialize the potentials
E_net = Population(NE, neuronClass, cell_params, label="E_net")
vinit_distr = RandomDistribution('uniform', [U0, firstbinedge])
E_net.initialize('v', vinit_distr)

# Create a source population
ex_spikes = Population(NE, SpikeSourcePoisson, {'rate' : FR_bg}, "ex_spikes")

# Record spike-times for target neural population:
E_net.record()

# Connect the source to the target using a user-defined connector
bg_Connector_E = OneToOneConnector(weights = J, delays = None, verbose = True)
E_input = Projection(ex_spikes, E_net, bg_Connector_E, target="excitatory")

# Read out time used for building
buildCPUTime = timer.elapsedTime()

###################################################################################
# Run Simulation
###################################################################################

# Run, measure computer time
timer.start() # start timer on construction
print "%d Running simulation for %g ms." % (rank, simtime)
run(simtime)
simCPUTime = timer.elapsedTime()

E_rate = E_net.meanSpikeCount()*1000.0/simtime

# Report run times
print("Build time         : %g s" % buildCPUTime)   
print("Simulation time    : %g s" % simCPUTime)
print("Excitatory rate    : %g Hz" % E_rate)
print("dt                 : %g ms" %dt)

# Write spike data to file
E_net.printSpikes("EIFCO_5mV100HzFin1920p05_st_nest")

end()