/* --------------------------------------------------------------
   init-octopus-simulation.hoc
   NEURON code to simulate dendritic NMDA spike in a ball and stick model.

   To run:
   1. Compile mod files in the folder into a new "special".
   2. Run Terminal and change to directory containing this file, e.g.:
   3. Execute the following:
        mod/x86_64/special init-octopus-simulation-

   Suraj Honnuraiah     Last edit: 04-Feb-25
   -------------------------------------------------------------- */

strdef modelName, loadProgram, cellName, outputFile, cellPath, loadProgram
objref trunc, secR, fi,vC, mbSec
objref sh, axonal, dendritic, dendritic_only, initZone, stimSec
objref trunc,mbSec, middleSec
objref mydendrites

load_file("stdgui.hoc")
load_proc("nrnmainmenu")
objref sh, st, axonal, dendritic, dendritic_only
objref mydendrites,tuftS,temp,mytuft, mbSec


objref ampa, nmda, gaba, s, datafile, pploc[300], inhloc[300], gly[300], ampa[300]
objref ncl, weight
objref apc
strdef name
create soma
access soma
    
    
tstop = 1000
steps_per_ms = 40
dt = 0.025
    
    // --------------------------------------------------------------
    // passive & active membrane
    // --------------------------------------------------------------
    
    ra        = 150
    global_ra = ra
    rm        = 5000
    c_m       = 0.9
    v_init    = -65
    celsius   = 37
    resting_Vm=-65
    
    
    nsyn=10
    
    
    scl=1e-1
    ghf=1
    
    
    gna_soma=0.83
    gna_dend=0.83
    
    gh_soma=0.07*scl*ghf                //0.07
    gh_dend=0.07*scl*ghf                //0.07
    gka_soma=0.07*scl               //0.07
    gka_dend=0.07*scl               //0.07
    gkht_soma=0.1875*scl               //0.1875
    gkht_dend=0.1875*scl               //0.1875
    gklt_soma=0.75*scl                //0.75
    gklt_dend=0.75*scl                //0.75
    
    ginhb_max=0.01
    gexct_max=0.005
    
    
    /*
    gh_soma=0
    gh_dend=0
    gka_soma=0
    gka_dend=0
    gkht_soma=0
    gkht_dend=0
    gklt_soma=0
    gklt_dend=0
    g_leak=0.0016666667
    */
    
    
    
    proc init_cell() {
        
        // passive
        forall {
            insert pas
            Ra = ra
            cm = c_m
            g_pas = 1/rm
            e_pas = v_init
        }
    }
    
    
    proc load_3dcell() {
        
        // $s1 filename
        
        aspiny = 0
        forall delete_section()
            xopen($s1)
            access soma
            
            dendritic = new SectionList()
            
            // make sure no compartments exceed 50 uM length
            forall {
                diam_save = diam
                n = L/50
                nseg = n + 1
                if (n3d() == 0) diam = diam_save
                    dendritic.append()
                }
                
                dendritic_only = new SectionList()
                forsec dendritic dendritic_only.append()
                    soma  dendritic_only.remove()
                    
                    //create_axon()
                    
                    init_cell()
                    
                    //if (!aspiny) add_spines(dendritic_only,spine_dens)
                }
                load_3dcell("morphology_octopus.hoc")
                
                for i=0,116 {
                    
                    pploc[i]= new SectionList ()
                    inhloc[i]=new SectionList ()
                    dend[i] pploc[i].append()
                }
                
                
                for j=0, 116{
                    
                    forsec pploc[j] {
                        
                        insert jsrnaf ena_jsrnaf = 50
                        insert klt ek_lt = -70
                        insert kht ek_kht = -70
                        insert na ena_na = 50
                        insert ka ek_ka=-70
                        insert hcnobo eh_hcnobo = -38
                        insert leak g_leak=0 erev_leak = -65
                        
                        gnabar_jsrnaf = gna_soma
                        gkhtbar_kht = gkht_soma
                        gkltbar_klt = gklt_soma
                        gkabar_ka = gka_soma
                        gbar_hcnobo = gh_soma
                        
                    }
                    
                }
                
                soma {
                    insert jsrnaf ena_jsrnaf = 50
                    insert klt ek_lt = -70
                    insert kht ek_kht = -70
                    insert na ena_na = 50
                    insert ka ek_ka=-70
                    insert hcnobo eh_hcnobo = -38
                    insert leak g_leak=0  erev_leak = -65
                    
                    gnabar_jsrnaf = gna_soma
                    gkhtbar_kht = gkht_soma
                    gkltbar_klt = gklt_soma
                    gkabar_ka = gka_soma
                    gbar_hcnobo = gh_soma
                }
                
                
                
                // distal dendrites list
                
                dend[19] inhloc[0].append
                dend[17] inhloc[1].append
                dend[51] inhloc[2].append
                dend[52] inhloc[3].append
                dend[49] inhloc[4].append
                dend[98] inhloc[5].append
                dend[65] inhloc[6].append
                dend[64] inhloc[7].append
                
                
                // proximal dendrites list
                
                dend[79] inhloc[8].append
                dend[34] inhloc[9].append
                dend[32] inhloc[10].append
                dend[69] inhloc[11].append
                dend[88] inhloc[12].append
                dend[91] inhloc[13].append
                dend[40] inhloc[14].append
                dend[93] inhloc[15].append
                dend[16] inhloc[16].append
                dend[78] inhloc[17].append
                //soma     inhloc[18].append
                
                
                for z=0,nsyn{
                    
                    /*
                    soma {gly[z]=new syn_gi(1)
                    gly[z].onset=400
                    gly[z].gmax=ginhb_max
                    gly[z].tau0=0.5
                    gly[z].tau1=5
                }
                */
                
                forsec inhloc[z] {gly[z]=new syn_gi(0)
                    gly[z].onset=400
                    gly[z].gmax=ginhb_max
                    gly[z].tau0=0.5
                    gly[z].tau1=5
                }
                
                
            }
            
            
            for e=0,nsyn{
                
                forsec inhloc[e] {ampa[e]=new syn_g_Oct(1)
                    ampa[e].onset=400
                    ampa[e].gmax=gexct_max
                }
                
                
            }
            
            
            
            
            objectvar stim,stim2,stim3,stim4,stim5,stim6,stim7,stim8,stim9,stim10
                
                amp=0
                dur=100
                soma stim = new IClamp(0.5)
                
                stim.del = 400
                stim.dur = dur
                stim.amp = amp
                
                /*
                amp=1
                dur=300
                dend[113] stim1 = new IClamp(0.5)
                
                stim1.del = 100
                stim1.dur = dur
                stim1.amp = amp
                
                
                amp=1
                dur=300
                dend[169] stim2 = new IClamp(0.5)
                
                stim2.del = 100
                stim2.dur = dur
                stim2.amp = amp
                
                
                amp=1
                dur=300
                dend[0] stim3 = new IClamp(0.5)
                
                stim3.del = 100
                stim3.dur = dur
                stim3.amp = amp
                
                
                amp=1
                dur=300
                dend[6] stim4 = new IClamp(0.5)
                
                stim4.del = 100
                stim4.dur = dur
                stim4.amp = amp
                
                amp=1
                dur=300
                dend[7] stim5 = new IClamp(0.5)
                
                stim5.del = 100
                stim5.dur = dur
                stim5.amp = amp
                
                */
                
                objref apc
                proc update_apc(){
                    soma{
                        apc=new APCount(0.5)
                        apc.thresh=0
                    }
                }
                
                
                
                proc update_current(){
                    
                    stim.amp = amp
                    
                }
                
                
                proc update_gna_soma(){
                    
                    soma{
                        gnabar_jsrnaf=gna_soma
                    }
                    
                }
                
                proc update_gna_dend(){
                    
                    for n=0, 116{
                        
                        forsec pploc[n] {
                            gbar_na = gna_dend
                        }
                    }
                    
                }
                
                
                proc update_gka_dend(){
                    
                    for n=0, 116{
                        
                        forsec pploc[n] {
                            gkabar_ka = gka_dend
                        }
                    }
                }
                
                proc update_gka_soma(){
                    
                    soma{
                        gkabar_ka = gka_soma
                    }
                    
                }
                
                proc update_gh_dend(){
                    
                    for n=0, 116{
                        
                        forsec pploc[n] {
                            gbar_hcnobo = gh_dend
                        }
                    }
                    
                }
                
                
                proc update_gh_soma(){
                    
                    soma{
                        gbar_hcnobo = gh_soma
                    }
                    
                }
                
                
                
                proc update_gkht_dend(){
                    
                    for n=0, 116{
                        
                        forsec pploc[n] {
                            gkhtbar_kht = gkht_dend
                        }
                    }
                }
                
                proc update_gkht_soma(){
                    
                    soma{
                        gkhtbar_kht = gkht_soma
                    }
                    
                }
                
                proc update_gklt_dend(){
                    
                    for n=0, 116{
                        
                        forsec pploc[n] {
                            gkltbar_klt = gklt_dend
                        }
                    }
                }
                
                proc update_gklt_soma(){
                    
                    soma{
                        gkltbar_klt = gklt_soma
                    }
                    
                }
                
                
                proc update_isyn_strength(){
                    
                    for l=0, nsyn{
                        
                        forsec inhloc[l] {
                            gly[l].gmax=ginhb_max
                        }
                        
                        
                    }
                    
                }
                
                
                proc update_esyn_strength(){
                    
                    for m=0, nsyn{
                        
                        forsec inhloc[m] {
                            ampa[m].gmax=gexct_max
                        }
                        
                        
                    }
                    
                }
                
                proc update_Rm(){
                    
                    forall{
                        g_pas = 1/rm
                    }
                    
                }
                
                
                /********************************************************************/
                
                proc update_init(){
                    finitialize(v_init)
                    fcurrent()
                }
                
                
                /*
                proc Stimulus_Frequency() {
                    update_apc()
                    b=0
                    for(i=0;i<=6;i+=1) {
                        gka_soma=0.5*i
                        update_gka_soma()
                        print"gka_soma:", gka_soma
                        b=b+1
                        file_index+=b
                        sprint(name,"%.1f",file_index)
                        wopen(name)
                        for(k=0;k<=35;k+=1){
                            amp=0.05*k
                            update_current()
                            finitialize()
                            while (t<tstop){
                                fadvance ()
                            }
                            print"",amp,apc.n
                            fprint("%f\t%f\n",amp,apc.n)
                        }
                    }
                }
                
                
                Stimulus_Frequency()
                */
                
                objectvar g[20]         // max 20 graphs
                    ngraph = 0
                    
                    proc addgraph() { local ii  // define subroutine to add a new graph
                        // addgraph("variable", minvalue, maxvalue)
                        ngraph = ngraph+1
                        ii = ngraph-1
                        g[ii] = new Graph()
                        g[ii].size(0,tstop,$2,$3)
                        g[ii].xaxis()
                        g[ii].yaxis()
                        g[ii].addvar($s1,1,0)
                        g[ii].save_name("graphList[0].")
                        graphList[0].append(g[ii])
                    }
                    
                    if(ismenu==0) {
                        nrnmainmenu()         // create main menu
                        nrncontrolmenu()      // crate control menu
                        ismenu=1
                    }
                    
                    objref stfunc,shape
                    shape=new Shape(0)
                    shape.view(-792.112, -154.611, 1363.32, 1185.94, 5, 340, 300, 260)
                    shape.show(0)
                    proc make_shape_plot(){//DRAWS THE POINTS ON THE CELL
                        shape.point_mark_remove()
                        
                        for n=0,nsyn{
                            shape.point_mark(gly[n], 2, 4, 12)
                            shape.point_mark(ampa[n], 3, 4, 12)
                        }
                        
                        
                    }///END SHAPE
                    
                    
                    xpanel("Parameters")
                    xvalue("Clamp Amplitude","amp",1,"update_current()",1,1)
                    xvalue("Rm","rm",1,"update_Rm()",1,1)
                    xvalue("iSynapse Strength","ginhb_max",1,"update_isyn_strength()",1,1)
                    xvalue("eSynapse Strength","gexct_max",1,"update_esyn_strength()",1,1)
                    xvalue("Dend_GNa","gna_dend",1,"update_gna_dend()",1,1)
                    xvalue("Soma_GNa","gna_soma",1,"update_gna_soma()",1,1)
                    xvalue("Dend_Gh","gh_dend",1,"update_gh_dend()",1,1)
                    xvalue("Soma_Gh","gh_soma",1,"update_gh_soma()",1,1)
                    xvalue("Dend_Gka","gka_dend",1,"update_gka_dend()",1,1)
                    xvalue("Soma_Gka","gka_soma",1,"update_gka_soma()",1,1)
                    xvalue("Soma_GKHT","gkht_soma",1,"update_gkht_soma()",1,1)
                    xvalue("Dend_GKHT","gkht_dend",1,"update_gkht_dend()",1,1)
                    xvalue("Soma_GKLT","gklt_soma",1,"update_gklt_soma()",1,1)
                    xvalue("Dend_GKLT","gklt_dend",1,"update_gklt_dend()",1,1)
                    
                    
                    
                    xpanel()
                    
                    
                    addgraph("soma.v(0.5)",-75,90)
                    addgraph("dend[19].v(0)",-75,90)
                    addgraph("dend[79].v(0)",-75,90)
                    addgraph("gly[0].i",-1,1)
                    addgraph("gly[8].i",-1,1)
                    //addgraph("gly[1].i",-1,1)
                    //addgraph("dend[79].v(0.5)",-75,90)
                    //addgraph("dend[0].v(0.5)",-75,90)
                    //addgraph("dend[6].v(0.5)",-75,90)
                    //addgraph("dend[7].v(0.5)",-75,90)
                    //addgraph("dend[113].v(0.5)",-75,90)
                    //addgraph("dend[169].v(0.5)",-75,90)
                    //addgraph("a3_12122222112.v(0.5)",-75,90)
                    
                    make_shape_plot()