def test_point_process_coord():
    h('synapse.loc(0.15)')
    h('access soma')
    x, y, z = cell.get_pp_coord(h.synapse)
    assert x == 1.5
    assert y == 0
    assert z == 0

def initModel(h,par,vecpar,recpar, verbose):
    '''
    Initializes the model.
    Creates the axon and so on.
    '''
    passValuesToNeuron(par,verb=verbose)
    h('{load_file("MRGnodeHFS.hoc")}')
    passValuesToNeuron(vecpar,[a for a in vecpar.keys()],verb=verbose)
	#h('{load_file("nrngui.hoc")}')
	#h('{load_proc("nrnmainmenu")}')
    h('{buildModel()}')
    if verbose:
        print("Passed parameters and built model.")
    # insert recorders and record action potential timestamps!
    segmentsToRecordV = [];
    segmentsNames = [int(h.intrinsicNode), int(h.HFSreferenceNode), int(h.axonnodes-1)] #segments to record
    for ii in segmentsNames:
       	segmentsToRecordV.append(h.node[ii](0.5))
    rec = None
    for ii in range(0,len(segmentsToRecordV)):
       	rec = insertRecorders(segmentsToRecordV[ii],{'node'+str(segmentsNames[ii]):'_ref_v'},rec)
    h.node[int(h.axonnodes-1)].push()
    apc = h.APCount(h.node[int(h.axonnodes-1)](0.5))
    apc.thresh               = 0
    rec['apc']               = h.Vector()
    apc.record(rec['apc'])
    #rec['electrodeWaveform'] = h.rec_electrode
    if verbose:
        print("Inserted recorders and APCount. Returning recorders.")
    ##########################################################################
    return rec

 def useSTDP(self, mechanism, parameters, ddf):
     """
     Set this connection to use spike-timing-dependent plasticity.
     
     `mechanism`  -- the name of an NMODL mechanism that modifies synaptic
                     weights based on the times of pre- and post-synaptic spikes.
     `parameters` -- a dictionary containing the parameters of the weight-
                     adjuster mechanism.
     `ddf`        -- dendritic delay fraction. If ddf=1, the synaptic delay
                     `d` is considered to occur entirely in the post-synaptic
                     dendrite, i.e., the weight adjuster receives the pre-
                     synaptic spike at the time of emission, and the post-
                     synaptic spike a time `d` after emission. If ddf=0, the
                     synaptic delay is considered to occur entirely in the
                     pre-synaptic axon.
     """
     self.ddf = ddf
     self.weight_adjuster = getattr(h, mechanism)(0.5)
     self.pre2wa = state.parallel_context.gid_connect(int(self.source), self.weight_adjuster)
     self.pre2wa.threshold = self.nc.threshold
     self.pre2wa.delay = self.nc.delay * (1-ddf)
     self.pre2wa.weight[0] = 1
     # directly create NetCon as wa is on the same machine as the post-synaptic cell
     self.post2wa = h.NetCon(self.target._cell.source, self.weight_adjuster,
                             sec=self.target._cell.source_section)
     self.post2wa.threshold = 1
     self.post2wa.delay = self.nc.delay * ddf
     self.post2wa.weight[0] = -1
     for name, value in parameters.items():
         setattr(self.weight_adjuster, name, value)
     # setpointer
     i = len(h.plastic_connections)
     h.plastic_connections.append(self)
     h('setpointer plastic_connections._[%d].weight_adjuster.wsyn, plastic_connections._[%d].nc.weight' % (i,i))

    def setup_sections(self):
        start = h.startsw()
        
        ###################################################
        # set up sections
        self.sections = []
        # old style, but it is need for section_name in hoc
        h(self.section_def_template % (self.name, len(self.tree)))
        for sec in h.allsec():
            self.sections.append(sec)


        ###################################################
        # connect sections
        for i,sec in enumerate(self.sections):
            parent = self.tree[i]
            #print "%d to %d" % (i, tree[i])
            if(parent != 0):
                sec.connect(self.sections[parent-1], 1, 0)

        self.num_compartment = 0
        for sec in h.allsec():
            self.num_compartment += 1

        self.setup_time += h.startsw() - start

def test_get_point_processes():
    h('synapse.loc(0.15)')
    h('access soma')
    pps = cell.get_point_processes()
    assert len(pps)==1
    assert pps[0][0].same(h.synapse)
    assert pps[0][1] == 1.5

def load_morphology(filename):
    swc = h.Import3d_SWC_read()
    swc.input(filename)

    imprt = h.Import3d_GUI(swc, 0)
    h("objref this")
    imprt.instantiate(h.this)

def test_current_balance_synapse_at_section_end():

    h('synapse.loc(0)')
    cell.initialize(dt=0.025)
    t, I = cell.integrate(1)
    coord = cell.get_seg_coords()
    assert (np.abs(total_current(coord, I))<1e-6).all()

    def test_calc_potential_from_multi_dipoles01(self):
        neuron.h('forall delete_section()')
        soma = neuron.h.Section(name='soma')
        dend1 = neuron.h.Section(name='dend1')
        dend2 = neuron.h.Section(name='dend2')
        dend1.connect(soma(0.5), 0)
        dend2.connect(dend1(1.0), 0)
        morphology = neuron.h.SectionList()
        morphology.wholetree()
        radii = [300, 400, 500, 600]
        sigmas = [0.3, 1.5, 0.015, 0.3]
        electrode_locs = np.array([[0., 0., 290.],
                                   [10., 90., 300.],
                                   [-90, 50., 400.],
                                   [110.3, -100., 500.]])
        cell = cell_w_synapse_from_sections(morphology)
        t_point = -1

        MD_4s = LFPy.FourSphereVolumeConductor(radii, sigmas, electrode_locs)
        dipoles, dipole_locs = cell.get_multi_current_dipole_moments()
        p = dipoles[:,t_point,:]
        Np = p.shape[0]
        Nt = 1
        Ne = electrode_locs.shape[0]
        pot_MD = MD_4s.calc_potential_from_multi_dipoles(cell)[:,t_point]
        pot_sum = np.zeros((Ne, Nt))
        for i in range(Np):
            dip = np.array([p[i]])
            dip_loc = dipole_locs[i]
            fs = LFPy.FourSphereVolumeConductor(radii, sigmas, electrode_locs)
            pot = fs.calc_potential(dip, dip_loc)
            pot_sum += pot
        pot_sum = pot_sum.reshape(4)
        np.testing.assert_almost_equal(pot_MD, pot_sum)
        np.testing.assert_allclose(pot_MD, pot_sum, rtol=1E-4)

    def create_head(self, neck, head_vol, big_spine):
        """Create the head of the spine and populate it with the right channels"""
        name_sec = self.id + "_head"
        h("create " + name_sec)
        head = getattr(h, name_sec)
        
        
        if big_spine:
            head.L = 1
            head.diam = 1.175
            r = head.diam/2.
            self.head_vol = math.pi * r * r * head.L
        else:
            head.L = 0.5
            head.diam = math.sqrt(head_vol / (head.L * math.pi) ) * 2
        
        self.Ra = 150.0
        head.nseg = 1
        head.connect(neck)
        
        #head.insert("pas")
        head.insert("kir")
        head.insert("can")
        head.insert("caq")
        head.insert("car")
        head.insert("skkca")

        h.factors_caltrack()
        head.insert("caltrack")
        
        h.factors_catrack()
        head.insert("catrack")
        
        return head

def test_current_balance_synapse_in_segment():
    h('synapse.loc(0.1)')
    h('soma {insert pas}')
    cell.initialize(dt=0.025)
    t, I = cell.integrate(1)
    coord = cell.get_seg_coords()
    assert (np.abs(total_current(coord, I))<1e-6).all()

    def test_get_multi_dipole_potential00(self):
        neuron.h('forall delete_section()')
        soma = neuron.h.Section(name='soma')
        dend1 = neuron.h.Section(name='dend1')
        dend2 = neuron.h.Section(name='dend2')
        dend3 = neuron.h.Section(name='dend3')
        dend4 = neuron.h.Section(name='dend4')
        dend5 = neuron.h.Section(name='dend5')
        dend1.connect(soma(0.5), 0)
        dend2.connect(dend1(1.0), 0)
        dend3.connect(dend2(1.0), 0)
        dend4.connect(dend3(1.0), 0)
        dend5.connect(dend4(1.0), 0)
        morphology = neuron.h.SectionList()
        morphology.wholetree()
        electrode_locs = np.array([[0., 0., 10000.]])
        cell, electrode = cell_w_synapse_from_sections_w_electrode(morphology, electrode_locs)
        sigma = 0.3
        t_point = 0

        MD_inf = LFPy.InfiniteVolumeConductor(sigma)
        pot_MD = MD_inf.get_multi_dipole_potential(cell, electrode_locs)
        pot_cb = electrode.LFP

        np.testing.assert_almost_equal(pot_MD, pot_cb)
        np.testing.assert_allclose(pot_MD, pot_cb, rtol=1E-4)

 def _do_callback(self):
     if callable(self._callable):
         self._callable()
     else:
         h(self._callable)
     if self._thread is not None:
         self.start()

def passive_soma(quad, show=False):
  """
  Creates the model with basic pyramidal passive properties.
  """
  # Load the hoc into neuron
  h('xopen(%s)' %quad.hocfile)
  h.load_file('stdrun.hoc')
  seclist = list(h.allsec())
  for sec in seclist:
    sec.insert('pas')
    sec.Ra = 200.
  
  # Current injection into soma or tip
  soma_sec = return_soma_seg(quad, h)
  stim_loc = 0.
  stim = h.IClamp(stim_loc, sec=soma_sec)
  stim.delay = 1 # ms
  stim.dur = 1 # ms
  stim.amp = 20 # nA
  
  # Run sim and record data
  (v, labels) = ez_record(h) # PyNeuron to record all compartments
  t, I = h.Vector(), h.Vector()
  t.record(h._ref_t)
  I.record(stim._ref_i)
  h.init()
  h.tstop = 10 # s?
  h.run()
  v = ez_convert(v) # Convert v to numpy 2D array
  
  # If show, plot, else just return v
  if show:

    def alloc_synapse_ff(self,r,post_syn,cellind,k,gidn,i):

        NCELL=self.NCELL
        SIZE=self.SIZE
        COMM = self.COMM
        RANK=self.RANK
        #from neuron import h
        h=self.h  
        pc=h.ParallelContext()
        polarity = 0        
        polarity=int(h.Cell[int(cellind)].polarity)
        if polarity==1:
            #TODO pickle load the graphs here instead of making them manually.
            self.ecm[i][gidn] = self.ecm[i][gidn] + 1
            self.ecg.add_edge(i,gidn,weight=r/0.4)
            assert np.sum(self.ecm)!=0
        else:
            self.icm[i][gidn] = self.icm[i][gidn] + 1
            self.icg.add_edge(i,gidn,weight=r/0.4)
            assert np.sum(self.icm)!=0                
            #TODO Add other edge attributes like secnames etc.
        print post_syn
        h('objref syn_')   
        h(post_syn)
        syn_=h.syn_
        h.syn_.cid=i
        h.Cell[cellind].ampalist.append(h.syn_)
        h.Cell[cellind].div.append(k['gid'])
        h.Cell[cellind].gvpre.append(k['gid'])
        nc=pc.gid_connect(k['gid'],syn_)                                        
        nc.threshold = -20
        nc.delay=1+r/0.4
        nc.weight[0]=r/0.4    
        self.nclist.append(nc)

    def test_Network_04(self):
        cellParameters = dict(
            morphology=os.path.join(LFPy.__path__[0], 'test', 'ball_and_sticks_w_lists.hoc'),
            templatefile=os.path.join(LFPy.__path__[0], 'test', 'ball_and_stick_template.hoc'),
            templatename='ball_and_stick_template',
            templateargs=None,
            passive=True,
            dt=2**-3,
            tstop=100,
            delete_sections=False,
        )
        
        synapseParameters = dict(idx=0, syntype='Exp2Syn', weight=0.002,
                                 tau1=0.1, tau2=0.1, e=0)
        
        populationParameters = dict(
            CWD=None,
            CELLPATH=None,
            Cell=LFPy.NetworkCell,
            cell_args = cellParameters,
            pop_args = dict(
                radius=100,
                loc=0.,
                scale=20.),
            rotation_args = dict(x=0, y=0),
            POP_SIZE = 1,
            name = 'test',
        )
        networkParameters = dict(
            dt=2**-3,
            tstart=0.,
            tstop=100.,
            v_init=-70.,
            celsius=6.3,
            OUTPUTPATH='tmp_testNetworkPopulation'
            )
        # set up
        network = LFPy.Network(**networkParameters)
        network.create_population(**populationParameters)
        
        cell = network.populations['test'].cells[0]
        
        # create synapses
        synlist = []
        numsynapses = 2
        for i in range(numsynapses):
            synlist.append(LFPy.Synapse(cell=cell, **synapseParameters))
            synlist[-1].set_spike_times(np.array([10+(i*10)]))
        
        network.simulate()
        
        # test that the input results in the correct amount of PSPs
        np.testing.assert_equal(ss.argrelextrema(cell.somav, np.greater)[0].size, numsynapses)

        # clean up
        network.pc.gid_clear()
        os.system('rm -r tmp_testNetworkPopulation')
        neuron.h('forall delete_section()')

 def alloc_synapse(self,r,h,sec,seg,cellind,secnames,k,i,gidn):
     '''
     Allocate a synaptic cleft from exhuastive collision detection.
     '''
     NCELL=self.NCELL
     SIZE=self.SIZE
     COMM = self.COMM
     RANK=self.RANK
     from neuron import h
     pc=h.ParallelContext()
     h=self.h
     self.visited[i][gidn] = self.visited[i][gidn] + 1              
     if r < 2.5: #2.5 micro metres.
         polarity = 0        
         polarity=int(h.Cell[int(cellind)].polarity)
         h('objref syn_')        
         if int(polarity) == int(0):
             post_syn = secnames + ' ' + 'syn_ = new FastInhib(' + str(seg.x) + ')'
             #post_syn = secnames + ' ' + 'syn_ = new GABAa(' + str(seg.x) + ')'
             self.icm[i][gidn] = self.icm[i][gidn] + 1
             self.icg.add_edge(i,gidn,weight=r/0.4)
             self.icg[i][gidn]['post_loc']=secnames
             self.icg[i][gidn]['pre_loc']=k['secnames']
             assert np.sum(self.icm)!=0
     
         else:
             if (k['gid']%2==0):
                 #TODO Find standard open source brain affiliated code for NMDA synapse
                 post_syn = secnames + ' ' + 'syn_ = new AmpaNmda(' + str(seg.x) + ')'                       
                 self.ecm[i][gidn] = self.ecm[i][gidn] + 1
                 self.ecg.add_edge(i,gidn,weight=r/0.4)
                 self.ecg[i][gidn]['post_loc']=secnames
                 self.ecg[i][gidn]['pre_loc']=k['secnames']
                 self.seclists.append(secnames)
                 assert np.sum(self.ecm)!=0
             else:
                 #TODO Find standard open source brain affiliated code for NMDA synapse
                 post_syn = secnames + ' ' + 'syn_ = new ExpSid(' + str(seg.x) + ')'                       
                 self.ecm[i][gidn] = self.ecm[i][gidn] + 1
                 self.ecg.add_edge(i,gidn,weight=r/0.4)
                 self.ecg[i][gidn]['post_loc']=secnames
                 self.ecg[i][gidn]['pre_loc']=k['secnames']
                 self.seclists.append(secnames)
                 assert np.sum(self.ecm)!=0
     
         h(post_syn)
         print post_syn
         self.synapse_list.append((r,post_syn,cellind,k,gidn,i))
         syn_=h.syn_
         h.syn_.cid=i
         h.Cell[cellind].ampalist.append(h.syn_)
         h.Cell[cellind].div.append(k['gid'])
         h.Cell[cellind].gvpre.append(k['gid'])
         nc=pc.gid_connect(k['gid'],syn_)                                        
         nc.threshold = -20
         nc.delay=1+r/0.4
         nc.weight[0]=r/0.4    
         self.nclist.append(nc)

    def destroy(self):
        """
        Takes care that neuron objects and python objects reffering to the cell
        are destroyed.

        At this point it actually destroys all the sections.
        """
        neuron.h('objref tuft')
        super(self.__class__, self).destroy()

 def __init__(self):
     """
     Create an `FinitializeHandler` object in Hoc, which will call the
     `_initialize()` method when NEURON is initialized.
     """
     h('objref initializer')
     h.initializer = self
     self.fih = h.FInitializeHandler(1, "initializer._initialize()")
     self.clear()

    def _setup_plasticity(self, synapse_type, parameters):
        """
        Set this connection to use spike-timing-dependent plasticity.

        `mechanism`  -- the name of an NMODL mechanism that modifies synaptic
                        weights based on the times of pre- and post-synaptic spikes.
        `parameters` -- a dictionary containing the parameters of the weight-
                        adjuster mechanism.
        """
        mechanism = synapse_type.model
        self.weight_adjuster = getattr(h, mechanism)(0.5)
        if synapse_type.postsynaptic_variable == "spikes":
            parameters["allow_update_on_post"] = int(False)  # for compatibility with NEST
            self.ddf = parameters.pop("dendritic_delay_fraction")
            # If ddf=1, the synaptic delay
            # `d` is considered to occur entirely in the post-synaptic
            # dendrite, i.e., the weight adjuster receives the pre-
            # synaptic spike at the time of emission, and the post-
            # synaptic spike a time `d` after emission. If ddf=0, the
            # synaptic delay is considered to occur entirely in the
            # pre-synaptic axon.
        elif synapse_type.postsynaptic_variable is None:
            self.ddf = 0
        else:
            raise NotImplementedError("Only post-synaptic-spike-dependent mechanisms available for now.")
        self.pre2wa = state.parallel_context.gid_connect(int(self.presynaptic_cell), self.weight_adjuster)
        self.pre2wa.threshold = self.nc.threshold
        self.pre2wa.delay = self.nc.delay * (1 - self.ddf)
        if self.pre2wa.delay > 1e-9:
            self.pre2wa.delay -= (
                1e-9
            )  # we subtract a small value so that the synaptic weight gets updated before it is used.
        if synapse_type.postsynaptic_variable == "spikes":
            # directly create NetCon as wa is on the same machine as the post-synaptic cell
            self.post2wa = h.NetCon(
                self.postsynaptic_cell._cell.source,
                self.weight_adjuster,
                sec=self.postsynaptic_cell._cell.source_section,
            )
            self.post2wa.threshold = 1
            self.post2wa.delay = self.nc.delay * self.ddf
            self.post2wa.weight[0] = -1
            self.pre2wa.weight[0] = 1
        else:
            self.pre2wa.weight[0] = self.nc.weight[0]
        parameters.pop("x", None)  # for the Tsodyks-Markram model
        parameters.pop("y", None)  # would be better to actually use these initial values
        for name, value in parameters.items():
            setattr(self.weight_adjuster, name, value)
        if (
            mechanism == "TsodyksMarkramWA"
        ):  # or could assume that any weight_adjuster parameter called "tau_syn" should be set like this
            self.weight_adjuster.tau_syn = self.nc.syn().tau
        # setpointer
        i = len(h.plastic_connections)
        h.plastic_connections.append(self)
        h("setpointer plastic_connections._[%d].weight_adjuster.wsyn, plastic_connections._[%d].nc.weight" % (i, i))

 def insert_glutamate_stim(cell, section="apic[63]", site=0.5):
     gmaxS = 10
     neuron.h("access %s" % section)
     glut_syn = neuron.h.glutamate(site)
     glut_syn.delay = 50
     glut_syn.ntar = 1.3
     glut_syn.gmax = gmaxS
     glut_syn.Nspike = 3
     glut_syn.Tspike = 20
     return glut_syn