source: project/release/4/nemo/trunk/examples/Golgi/Golgi_DeSouza10.nemo @ 29485

Last change on this file since 29485 was 29485, checked in by Ivan Raikov, 7 years ago

nemo: added ability to load sxslt files

File size: 75.5 KB
Line 
1
2;; Simoes de Souza FM and De Schutter E (2011) Robustness effect of
3;; gap junctions between Golgi cells on cerebellar cortex
4;; oscillations.  Neural Systems & Circuits 1:7.
5
6(nemo-model Golgi
7
8  (
9   (input (v (unit mV))
10          (celsius (unit degC)))
11
12   (const F = 96485.309)
13   (const R = 8.31342)
14
15   (defun ktf (celsius) (1000.0 * R * (celsius + 273.15) / F ))
16
17   (defun nernst (celsius ci co z)
18     (if (ci <= 0)
19         then 1e6
20         else (if (co <= 0)
21                  then -1e6
22                  else (ktf (celsius) / z * (ln (co / ci))))
23         ))
24   
25
26   (component (type defaults)
27      (const V_t     =  -35  (unit mV))
28      (const celsius =   23  (unit degC))
29      (const Ra      = 0.01  (unit ohm.cm))
30      (output celsius V_t Ra)
31      )
32
33
34   (component (type membrane-capacitance)
35         (const c = 1 (unit uf/cm2))
36         (output c )
37         )
38
39
40   (component (type geometry) (name soma)
41         (const L = 27    (unit um))
42         (const diam = 27 (unit um))
43         (output L diam ))
44
45
46   (component (type post-synaptic-conductance) (name AMPA)
47
48        (input (w from event (unit uS)))
49
50        (const tauA = 0.03 (unit ms)) ;; rise time
51        (const tauB = 0.5  (unit ms)) ;; decay time
52
53        (const  e = 0 (unit mV))
54
55        (const tau1 = (if ((tauA / tauB) > .9999) then (0.9999 * tauB) else tauA) (unit ms))
56        (const tau2 = tauB (unit ms))
57
58        (const tp = ((tau1 * tau2) / (tau2 - tau1) * ln (tau2 / tau1)))
59        (const scale_factor  = (1 / (neg (exp(neg (tp) / tau1)) + exp (neg (tp) / tau2))))
60
61        (transient (A) = (neg (A) / tau1) (onevent (A + (w * scale_factor))) (initial 0) (unit uS))
62        (transient (B) = (neg (B) / tau2) (onevent (B + (w * scale_factor))) (initial 0) (unit uS))
63
64        (g =  (B - A) (unit uS))
65
66        (output g e scale_factor)
67
68      )
69
70   (component (type decaying-pool) (name ca)
71
72      (input (ica from ion-currents))
73
74      (const  d       = 0.2  (unit um))
75      (const  cao     = 2.0  (unit mM))
76      (const  cai0    = 5e-5 (unit mM))
77      (const  beta    = 1.3  (unit /ms))
78     
79      (d (ca) =  ((neg (ica) / (2 * F * d)) * 1e4 -
80                  (beta * (ca - cai0)))
81         (initial cai0))
82     
83      (output ca cao)
84      )
85
86
87   (component (type decaying-pool) (name ca2)
88
89      (input (ica2 from ion-currents))
90
91      (const valence  = 2)
92      (const  d       = 0.2  (unit um))
93      (const  ca2o    = 2.0  (unit mM))
94      (const  ca2i0   = 5e-5 (unit mM))
95      (const  beta    = 1.3  (unit /ms))
96     
97      (d (ca2) =  ((neg (ica2) / (2 * F * d)) * 1e4 -
98                   (beta * (ca2 - ca2i0)))
99         (initial ca2i0) (unit mM))
100     
101      (output ca2 ca2o valence)
102      )
103
104
105   (component (type decaying-pool) (name na)
106
107      (input (ina from ion-currents))
108
109      (const  d       = 0.2   (unit um))
110      (const  nao0    = 140.0 (unit mM))
111      (const  nai0    = 5.0   (unit mM))
112      (const  beta    = 0.075 (unit /ms))
113     
114      (d (nai) =  ((neg (ina) / (2 * F * d)) * 1e4 - (beta * (nai - nai0)))
115         (initial nai0) (unit mM))
116     
117      (d (nao) =  ((ina / (2 * F * d)) * 1e4 - (beta * (nao - nao0)))
118         (initial nao0) (unit mM))
119     
120      (output nai nao)
121      )
122
123
124   (component (type decaying-pool) (name k)
125
126      (input (ik from ion-currents))
127
128      (const  d       = 0.2   (unit um))
129      (const  ko0     = 5.0   (unit mM))
130      (const  ki0     = 140.0 (unit mM))
131      (const  beta    = 0.075 (unit /ms))
132     
133      (d (ki) =  ((neg (ik) / (2 * F * d)) * 1e4 - (beta * (ki - ki0)))
134         (initial ki0) (unit mM))
135     
136      (d (ko) =  ((ik / (2 * F * d)) * 1e4 - (beta * (ko - ko0)))
137         (initial ko0) (unit mM))
138     
139      (output ki ko)
140      )
141
142
143   (component (type ionic-current) (name CaHVA )
144             
145      (input
146       (cai from ion-pools)
147       (cao from ion-pools))
148             
149      (component (type gate)
150
151                 ;; rate constants
152                 
153                 (Q10 = (pow (3 ((celsius - 20) / 10))))
154                 
155                 (const Aalpha_s  = 0.04944 (unit /ms))
156                 (const Kalpha_s  = 15.87301587302 (unit mV))
157                 (const V0alpha_s = -29.06 (unit mV))
158                 
159                 (const Abeta_s  = 0.08298 (unit /ms))
160                 (const Kbeta_s  = -25.641 (unit mV))
161                 (const V0beta_s = -18.66  (unit mV))
162                 
163                 (const Aalpha_u  = 0.0013  (unit /ms))
164                 (const Kalpha_u  = -18.183 (unit mV))
165                 (const V0alpha_u = -48     (unit mV))
166                 
167                 (const Abeta_u   = 0.0013  (unit /ms))
168                 (const Kbeta_u   = 83.33   (unit mV))
169                 (const V0beta_u  = -48     (unit mV))
170                 
171                 ;; rate functions
172                 
173                 (defun alpha_s (v Q10)
174                   (Q10 * Aalpha_s * exp ((v - V0alpha_s) / Kalpha_s)))
175                 
176                 (defun beta_s (v Q10)
177                   (Q10 * Abeta_s * exp ((v - V0beta_s) / Kbeta_s)))
178                 
179                 (defun alpha_u (v Q10)
180                   (Q10 * Aalpha_u * exp ((v - V0alpha_u) / Kalpha_u)))
181                 
182                 (defun beta_u (v Q10)
183                   (Q10 * Abeta_u * exp ((v - V0beta_u) / Kbeta_u)))
184
185                 (s_inf = ((alpha_s (v Q10)) / (alpha_s (v Q10) + beta_s (v Q10))))
186                 (tau_s = (1 / (alpha_s (v Q10) + beta_s (v Q10)) ) (unit ms))
187
188                 (u_inf = ((alpha_u (v Q10)) / (alpha_u (v Q10) + beta_u (v Q10)) ))
189                 (tau_u = (1 / (alpha_u (v Q10) + beta_u (v Q10)) ) (unit ms))
190
191                 (hh-ionic-gate
192                  (CaHVA  ;; ion name: exported variables will be of the form {ion}_{id}
193                   (initial-m  s_inf)
194                   (initial-h  u_inf)
195                   (m-power    2)
196                   (h-power    1)
197                   (m-inf      s_inf)
198                   (m-tau      tau_s)
199                   (h-inf      u_inf)
200                   (h-tau      tau_u)
201                   ))
202                 
203                 )
204             
205              (component (type pore)
206                         (const  gbar  = 460e-6 (unit S/cm2))
207                         (output gbar ))
208             
209              (component (type permeating-ion) (name ca)
210                         (e = (nernst (celsius cai cao 2)) (unit mV))
211                         (output e)
212                         )
213             
214              ) ;; end CaHVA current
215
216
217   (component (type voltage-clamp) (name CaHVA)
218             
219              (const vchold   = -71)
220              (const vcbase   = -69)
221              (const vcinc    = 10)
222              (const vcsteps  = 8)
223              (const vchdur   = 30)
224              (const vcbdur   = 100)
225
226              (component (type default-concentration) (name ca)
227                      (const cn       = 5e-5)
228                      (const cnout    = 2)
229                      (output cn cnout))
230             
231              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
232              )
233
234   
235   (component (type ionic-current) (name CaLVA )
236             
237      (input
238       (ca2i from ion-pools)
239       (ca2o from ion-pools))
240
241      (component (type gate)
242
243                 (const shift   = 2 (unit mV)) ; screening charge for Ca_o = 2 mM
244                 
245                 (const v0_m_inf = -50  (unit mV))
246                 (const v0_h_inf = -78  (unit mV))
247                 (const k_m_inf  = -7.4 (unit mV))
248                 (const k_h_inf  = 5.0  (unit mV))
249       
250                 (const C_tau_m   = 3)
251                 (const A_tau_m   = 1.0)
252                 (const v0_tau_m1 = -25  (unit mV))
253                 (const v0_tau_m2 = -100 (unit mV))
254                 (const k_tau_m1  = 10   (unit mV))
255                 (const k_tau_m2 = -15   (unit mV))
256                 
257                 (const C_tau_h   = 85)
258                 (const A_tau_h   = 1.0  (unit mV))
259                 (const v0_tau_h1 = -46  (unit mV))
260                 (const v0_tau_h2 = -405 (unit mV))
261                 (const k_tau_h1  = 4    (unit mV))
262                 (const k_tau_h2  = -50  (unit mV))
263                         
264                 ;; rate functions
265                 
266                 (phi_m = (pow (5.0 ((celsius - 24) / 10))))
267                 (phi_h = (pow (3.0 ((celsius - 24) / 10))))
268                 
269                 (m_inf = (1.0 / ( 1 + exp ((v + shift - v0_m_inf) / k_m_inf)) ))
270                 (h_inf = (1.0 / ( 1 + exp ((v + shift - v0_h_inf) / k_h_inf)) ))
271                 
272                 (tau_m = ( (C_tau_m + A_tau_m / ( exp ((v + shift - v0_tau_m1) / k_tau_m1) +
273                                                       exp ((v + shift - v0_tau_m2) / k_tau_m2) ) ) / phi_m)
274                        (unit ms))
275                 
276                 (tau_h = ( (C_tau_h + A_tau_h / ( exp ((v + shift - v0_tau_h1 ) / k_tau_h1) +
277                                                       exp ((v + shift - v0_tau_h2) / k_tau_h2) ) ) / phi_h)
278                        (unit ms))
279                 
280
281                 (hh-ionic-gate
282                  (CaLVA  ;; ion name: exported variables will be of the form {ion}_{id}
283                   (initial-m  m_inf)
284                   (initial-h  h_inf)
285                   (m-power    2)
286                   (h-power    1)
287                   (m-inf      m_inf)
288                   (m-tau      tau_m)
289                   (h-inf      h_inf)
290                   (h-tau      tau_h)
291                   ))
292                 )
293     
294     
295      (component (type pore)
296                 (const  gbar  = 2.5e-4 (unit S/cm2))
297                 (output gbar ))
298
299             
300      (component (type permeating-ion) (name ca2)
301                 (const valence = 2)
302                 (e = (nernst (celsius ca2i ca2o valence)) (unit mV))
303                 (output e valence))
304
305             
306      ) ;; end CaLVA current
307
308
309
310   (component (type voltage-clamp) (name CaLVA)
311             
312              (const vchold   = -71)
313              (const vcbase   = -69)
314              (const vcinc    = 10)
315              (const vcsteps  = 8)
316              (const vchdur   = 200)
317              (const vcbdur   = 30)
318
319              (component (type default-concentration) (name ca2)
320                      (const cn       = 5e-5)
321                      (const cnout    = 2)
322                      (output cn cnout))
323             
324              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
325              )
326
327
328     
329   (component (type ionic-current) (name HCN1)
330
331              (const Ehalf = -72.49 (unit mV))
332              (const c     = 0.11305 (unit /mV))
333             
334              (const rA = 0.002096 (unit /mV))
335              (const rB = 0.97596 )
336             
337              (defun r (potential)
338                ((rA * potential) + rB))
339             
340              (defun tau (potential t1 t2 t3)
341                (exp (((t1 * potential) - t2) * t3)))
342             
343              (defun o_inf (potential Ehalf c)
344                (1 / (1 + exp ((potential - Ehalf) * c))))
345             
346
347              (component (type gate)
348
349                         ;; rate constants
350       
351                         (const tCs = 0.01451)
352                         (const tDs = -4.056 (unit mV))
353                         (const tEs = 2.302585092 (unit /mV))
354
355                         (o_slow_inf = ((1 - r (v)) * o_inf (v Ehalf c)))
356                       
357                         (tau_s =  (tau (v tCs tDs tEs)) )
358
359                         (d (o_slow) = ((o_slow_inf - o_slow) / tau_s)
360                            (initial o_slow_inf))
361                         
362                         (output o_slow)
363
364                         )
365             
366              (component (type gate)
367                         
368
369                         ;; rate constants
370       
371                         (const tCf = 0.01371)
372                         (const tDf = -3.368      (unit mV))
373                         (const tEf = 2.302585092 (unit /mV))
374
375                         (o_fast_inf = (r (v) * o_inf (v Ehalf c)))
376                       
377                         (tau_f =  (tau (v tCf tDf tEf)) (unit ms))
378
379                         (d (o_fast) = ((o_fast_inf - o_fast) / tau_f)
380                            (initial o_fast_inf))
381                         
382                         (output o_fast)
383
384                         )
385             
386              (component (type pore)
387                         (const  gbar  = 5e-5 (unit S/cm2))
388                         (output gbar))
389             
390              (component (type permeating-ion) (name non-specific)
391                         (const e = -20 (unit mV))
392                         (output e ))
393             
394              ) ;; end HCN1 current
395
396
397             
398   (component (type voltage-clamp) (name HCN1)
399             
400              (const vchold   = -71)
401              (const vcbase   = -69)
402              (const vcinc    = 10)
403              (const vcsteps  = 8)
404              (const vchdur   = 30)
405              (const vcbdur   = 100)
406             
407              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
408              )
409
410
411     
412   (component (type ionic-current) (name HCN2)
413                         
414
415              (const Ehalf = -81.95 (unit mV))
416              (const c     = 0.1661 (unit /mV))
417             
418              ;; rate constants
419
420              (const rA = -0.0227 (unit /mV))
421              (const rB = -1.4694 )
422
423              (defun r (potential r1 r2)
424                (if (potential >= -64.70)
425                    then  0.0
426                    else (if (potential <= -108.70)
427                             then 1.0
428                             else (r1 * potential) + r2)))
429
430              (defun o_inf (potential Ehalf c)
431                (1 / (1 + exp((potential - Ehalf) * c))))
432             
433              (component (type gate)
434                         
435                         ;; rate constants
436                         
437                         (const tCs = 0.0152 )
438                         (const tDs = -5.2944 (unit mV))
439                         (const tEs = 2.3026  (unit /mV))
440                         
441                         (defun tau_slow (potential t1 t2 t3)
442                           (exp (t3 * ((t1 * potential) - t2))))
443
444                         (o_slow_inf = ((1 - r (v rA rB)) * o_inf (v Ehalf c)))
445
446                         (tau_s =  (tau_slow (v tCs tDs tEs)) (unit ms))
447
448                         (d (o_slow) = ((o_slow_inf - o_slow) / tau_s)
449                            (initial o_slow_inf))
450
451                         (output o_slow)
452
453                         )
454             
455              (component (type gate)
456                         
457
458                         ;; rate constants
459
460                         (const tCf = 0.0269)
461                         (const tDf = -5.6111 (unit mV))
462                         (const tEf = 2.3026  (unit /mV))
463                         
464                         (defun tau_fast (potential t1 t2 t3)
465                           (exp (t3 * ((t1 * potential) - t2))))
466
467                         (o_fast_inf = (r (v rA rB) * o_inf (v Ehalf c)))
468
469                         (tau_f =  (tau_fast (v tCf tDf tEf)) (unit ms))
470
471                         (d (o_fast) = ((o_fast_inf - o_fast) / tau_f)
472                            (initial o_fast_inf))
473
474                         (output o_fast)
475
476                         )
477             
478              (component (type pore)
479                         (const  gbar  = 8e-5 (unit S/cm2))
480                         (output gbar))
481             
482              (component (type permeating-ion) (name non-specific)
483                         (const e = -20 (unit mV))
484                         (output e ))
485             
486              ) ;; end HCN2 current
487
488
489
490   (component (type voltage-clamp) (name HCN2)
491             
492              (const vchold   = -71)
493              (const vcbase   = -69)
494              (const vcinc    = 10)
495              (const vcsteps  = 8)
496              (const vchdur   = 30)
497              (const vcbdur   = 100)
498             
499              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
500              )
501
502   
503
504   (component (type ionic-current) (name KA )
505
506              (input
507               (ki from ion-pools)
508               (ko from ion-pools)
509               )
510             
511             
512
513              (defun sigm (x y) (1 / (exp (x / y) + 1)))
514             
515              (component (type gate)
516                         
517                         ;; rate constants
518                         
519                         (Q10 = (pow (3 ((celsius - 25.5) / 10))))
520
521                         (const Aalpha_a  = 0.8147    (unit /ms))
522                         (const Kalpha_a  = -23.32708 (unit mV))
523                         (const V0alpha_a = -9.17203  (unit mV))
524                         (const Abeta_a   = 0.1655    (unit /ms))
525                         (const Kbeta_a   = 19.47175  (unit mV))
526                         (const V0beta_a  = -18.27914 (unit mV))
527
528                         (const Aalpha_b  = 0.0368    (unit /ms))
529                         (const Kalpha_b  = 12.8433   (unit mV))
530                         (const V0alpha_b = -111.33209 (unit mV))
531                         (const Abeta_b   = 0.0345    (unit /ms))
532                         (const Kbeta_b   = -8.90123  (unit mV))
533                         (const V0beta_b  = -49.9537  (unit mV))
534
535                         (const V0_ainf = -38 (unit mV))
536                         (const  K_ainf = -17 (unit mV))
537
538                         (const V0_binf   = -78.8 (unit mV))
539                         (const K_binf    = 8.4   (unit mV))
540
541                         
542                         ;; rate functions
543
544                         (defun alpha_a (v Q10)
545                           (Q10 * Aalpha_a * sigm((v - V0alpha_a) Kalpha_a)))
546
547                         (defun beta_a (v Q10)
548                           (Q10 * (Abeta_a / exp((v - V0beta_a) / Kbeta_a))))
549
550                         (defun alpha_b (v Q10)
551                           (Q10 * Aalpha_b * sigm((v - V0alpha_b) Kalpha_b)))
552
553                         (defun beta_b (v Q10)
554                           (Q10 * Abeta_b * sigm((v - V0beta_b) Kbeta_b)))
555
556
557
558                         (a_inf = (1 / (1 + exp ((v - V0_ainf) / K_ainf))))
559                         (tau_a = (1 / (alpha_a (v Q10) + beta_a (v Q10)) ))
560                         (b_inf = (1 / (1 + exp ((v - V0_binf) / K_binf))))
561                         (tau_b = (1 / (alpha_b (v Q10) + beta_b (v Q10)) ))
562
563
564                         (hh-ionic-gate
565                          (KA  ;; ion name: exported variables will be of the form {ion}_{id}
566                           (initial-m  (a_inf))
567                           (initial-h  (b_inf))
568                           (m-power    3)
569                           (h-power    1)
570                           (m-inf      a_inf)
571                           (m-tau      tau_a)
572                           (h-inf      b_inf)
573                           (h-tau      tau_b)
574                           ))
575
576                         )
577
578                     
579              (component (type pore)
580                         (const  gbar  = 0.008 (unit S/cm2))
581                         (output gbar ))
582
583             
584              (component (type permeating-ion) (name k)
585                         (e = (nernst (celsius ki ko 1)) (unit mV))
586                         (output e ))
587             
588              ) ;; end KA current
589             
590
591
592   (component (type voltage-clamp) (name KA)
593             
594              (const vchold   = -71)
595              (const vcbase   = -69)
596              (const vcinc    = 10)
597              (const vcsteps  = 8)
598              (const vchdur   = 30)
599              (const vcbdur   = 100)
600
601              (component (type default-concentration) (name k)
602                      (const cn    = 140)
603                      (const cnout = 5)
604                      (output cn cnout))
605             
606              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
607              )
608
609   
610
611   (component (type ionic-current) (name KCa )
612
613      (input
614       (cai from ion-pools)
615       (ki from ion-pools)
616       (ko from ion-pools)
617       )
618
619      (component (type gate)
620
621                 ;; rate constants
622                 (Q10 = (pow (3 ((celsius - 30) / 10))))
623                 
624                 (const Aalpha_c = 7 (unit /ms))
625                 (const Balpha_c = 1.5e-3 (unit mM))
626                 
627                 (const Kalpha_c =  -11.765 (unit mV))
628                 
629                 (const Abeta_c = 1.9 (unit /ms))
630                 (const Bbeta_c = 0.15e-3 (unit mM))
631
632                 (const Kbeta_c = -11.765 (unit mV))
633                 
634                 ;; rate functions
635                 (defun alpha_c (v cai Q10)
636                   (Q10 * Aalpha_c / (1 + (Balpha_c * exp(v / Kalpha_c) / cai))))
637                 
638                 (defun beta_c (v cai Q10)
639                   (Q10 * Abeta_c / (1 + (cai / (Bbeta_c * exp (v / Kbeta_c))) )))
640             
641                 (c_inf = ((alpha_c (v cai Q10)) / (alpha_c (v cai Q10) + beta_c (v cai Q10)) ))
642                 (tau_c = (1 / (alpha_c (v cai Q10) + beta_c (v cai Q10))) (unit ms))
643
644                 (hh-ionic-gate
645                  (KCa  ;; ion name: exported variables will be of the form {ion}_{id}
646                   (initial-m  c_inf)
647                   (m-power    1)
648                   (h-power    0)
649                   (m-inf      c_inf)
650                   (m-tau      tau_c)
651                   ))
652                 
653                 )
654     
655      (component (type pore)
656                 (const  gbar  = 0.003 (unit S/cm2))
657                 (output gbar ))
658     
659     
660      (component (type permeating-ion) (name k)
661                 (e = (nernst (celsius ki ko 1)) (unit mV))
662                 (output e ))
663     
664      (component (type modulating-ion) (name ca) )
665     
666      ) ;; end KCa current
667
668   
669   (component (type voltage-clamp) (name KCa)
670             
671              (const vchold   = -71)
672              (const vcbase   = -69)
673              (const vcinc    = 10)
674              (const vcsteps  = 8)
675              (const vchdur   = 30)
676              (const vcbdur   = 100)
677
678              (component (type default-concentration) (name k)
679                      (const kcn    = 140)
680                      (const kcnout = 5)
681                      (output kcn kcnout))
682
683              (const cnhold   = 5e-5)
684              (const cnbase   = 5e-5)
685              (const cninc    = 1e3)
686              (const cnsteps  = 1)
687              (const cnout    = 2)
688             
689              (output vchold vcbase vcinc vcsteps vchdur vcbdur
690                      cnhold cnbase cninc cnsteps cnout)
691              )
692
693
694   (component (type ionic-current) (name KM )
695             
696              (input
697               (ki from ion-pools)
698               (ko from ion-pools)
699               )
700             
701              (component (type gate)
702
703                         ;; rate constants
704                         (const Aalpha_n = 0.0033 (unit /ms))
705
706                         (const Kalpha_n  = 40  (unit mV))
707                         (const V0alpha_n = -30 (unit mV))
708                         (const Abeta_n   = 0.0033 (unit /ms))
709
710                         (const Kbeta_n  = -20 (unit mV))
711                         (const V0beta_n = -30 (unit mV))
712                         (const V0_ninf  = -35 (unit mV))
713                         (const B_ninf   = 6   (unit mV))
714                         
715                         (Q10 = (pow (3 ((celsius - 22) / 10))))
716                         
717                         ;; rate equations
718                         (defun alpha_n (v Q10)
719                           (Q10 * Aalpha_n * exp((v - V0alpha_n) / Kalpha_n) ))
720
721                         (defun beta_n (v Q10)
722                           (Q10 * Abeta_n * exp((v - V0beta_n) / Kbeta_n) ))
723
724                         (hh-ionic-gate
725                          (KM  ;; ion name: exported variables will be of the form {ion}_{id}
726                           (initial-m  ((1 / (1 + exp((neg (v - V0_ninf)) / B_ninf)))))
727                           (m-power    1)
728                           (h-power    0)
729                           (m-tau      ((1 / (alpha_n(v Q10) + beta_n (v Q10)) )))
730                           (m-inf      ((1 / (1 + exp((neg (v - V0_ninf)) / B_ninf)))))
731                           ))
732                         )
733             
734              (component (type pore)
735                         (const  gbar  = 0.001 (unit S/cm2))
736                         (output gbar ))
737             
738              (component (type permeating-ion) (name k)
739                         (e = (nernst (celsius ki ko 1)) (unit mV))
740                         (output e ))
741             
742              ) ;; end KM current
743
744
745
746   (component (type voltage-clamp) (name KM)
747             
748              (const vchold   = -71)
749              (const vcbase   = -69)
750              (const vcinc    = 10)
751              (const vcsteps  = 8)
752              (const vchdur   = 30)
753              (const vcbdur   = 100)
754
755              (component (type default-concentration) (name k)
756                      (const cn    = 140)
757                      (const cnout = 5)
758                      (output cn cnout))
759             
760              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
761              )
762
763   
764
765   (component (type ionic-current) (name KV )
766
767             
768              (input
769               (ki from ion-pools)
770               (ko from ion-pools)
771               )
772             
773              (defun linoid (x y)
774                (if ((abs (x / y)) < 1e-6)
775                    then (y * (1 - ((x / y) / 2)))
776                    else (x / (exp (x / y) - 1))
777                    ))
778
779
780              (component (type gate)
781
782                         ;; rate constants
783
784                         (const Aalpha_n  = -0.01 (unit /mV-ms))
785                         (const Kalpha_n  = -10   (unit mV))
786                         (const V0alpha_n = -26   (unit mV))
787                         (const Abeta_n   = 0.125 (unit /ms))
788       
789                         (const Kbeta_n  = -80 (unit mV))
790                         (const V0beta_n = -36 (unit mV))
791
792                         ;; rate functions
793                         (Q10 = (pow (3 ((celsius - 6.3) / 10))))
794
795                         (defun alpha_n (v Q10)
796                           (Q10 * Aalpha_n * linoid ((v - V0alpha_n) Kalpha_n)))
797
798                         (defun beta_n (v Q10)
799                           (Q10 * Abeta_n * exp((v - V0beta_n) / Kbeta_n) ))
800
801                         (n_inf = ((alpha_n (v Q10)) / (alpha_n (v Q10) + beta_n (v Q10)) ))
802                         (tau_n = (1 / (alpha_n (v Q10) + beta_n (v Q10))) (unit ms))
803
804                         (hh-ionic-gate
805                          (KV  ;; ion name: exported variables will be of the form {ion}_{id}
806                           (initial-m  n_inf)
807                           (m-power    4)
808                           (h-power    0)
809                           (m-inf     (n_inf))
810                           (m-tau     (tau_n))
811                           ))
812                         )
813
814              (component (type pore)
815                         (const  gbar  = 0.032 (unit S/cm2))
816                         (output gbar ))
817             
818              (component (type permeating-ion) (name k)
819                         (e = (nernst (celsius ki ko 1)) (unit mV))
820                         (output e ))
821             
822              ) ;; end KV current
823             
824
825   (component (type voltage-clamp) (name KV)
826             
827              (const vchold   = -71)
828              (const vcbase   = -69)
829              (const vcinc    = 10)
830              (const vcsteps  = 8)
831              (const vchdur   = 30)
832              (const vcbdur   = 100)
833
834              (component (type default-concentration) (name k)
835                      (const cn    = 140)
836                      (const cnout = 5)
837                      (output cn cnout))
838             
839              (output vchold vcbase vcinc vcsteps vchdur vcbdur)
840              )
841
842
843
844   
845   (component (type ionic-current) (name SK2 )
846
847      (input
848       (cai from ion-pools)
849       (ki from ion-pools)
850       (ko from ion-pools)
851       )
852
853      (component (type gate)
854
855                 (Q10 = (pow (3 ((celsius - 23) / 10))))
856                 
857                 (const diff = 3) ;; Diffusion factor
858                 
859                 ;; rates ca-independent
860                 (const invc1 = 80e-3  (unit /ms))
861                 (const invc2 = 80e-3  (unit /ms))
862                 (const invc3 = 200e-3 (unit /ms))
863
864                 (const invo1 = 1      (unit /ms))
865                 (const invo2 = 100e-3 (unit /ms))
866                 (const diro1 = 160e-3 (unit /ms))
867                 (const diro2 = 1.2    (unit /ms))
868
869                 ;; rates ca-dependent
870                 (const dirc2 = 200 (unit /ms))
871                 (const dirc3 = 160 (unit /ms))
872                 (const dirc4 = 80  (unit /ms))
873
874                 (invc1_t = (invc1 * Q10)  (unit /ms))
875                 (invc2_t = (invc2 * Q10)  (unit /ms))
876                 (invc3_t = (invc3 * Q10)  (unit /ms))
877                 (invo1_t = (invo1 * Q10)  (unit /ms))
878                 (invo2_t = (invo2 * Q10)  (unit /ms))
879                 (diro1_t = (diro1 * Q10)  (unit /ms))
880                 (diro2_t = (diro2 * Q10)  (unit /ms))
881                 (dirc2_t = (dirc2 * Q10)  (unit /ms))
882                 (dirc3_t = (dirc3 * Q10)  (unit /ms))
883                 (dirc4_t = (dirc4 * Q10)  (unit /ms))
884
885
886                 (dirc2_t_ca = (dirc2_t * (cai / diff)) )
887                 (dirc3_t_ca = (dirc3_t * (cai / diff)) )
888                 (dirc4_t_ca = (dirc4_t * (cai / diff)) )
889
890                 
891                 (reaction
892                  (SK2_z
893                   (transitions
894                    (<-> c1 c2 dirc2_t_ca invc1_t)
895                    (<-> c2 c3 dirc3_t_ca invc2_t)
896                    (<-> c3 c4 dirc4_t_ca invc3_t)
897                    (<-> c3 o1 diro1_t invo1_t)
898                    (<-> c4 o2 diro2_t invo2_t)
899                    )
900                   (conserve  ((1 = (c1 + c2 + c3 + c4 + o2 + o1))))
901                   (open o1 o2)
902                   (power 1)))
903                 
904                 (output SK2_z)
905                 
906                 )
907
908      (component (type pore)
909                 (const  gbar  = 0.038 (unit S/cm2))
910                 (output gbar ))
911     
912     
913      (component (type permeating-ion) (name k)
914                 (e = (nernst (celsius ki ko 1)) (unit mV))
915                 (output e ))
916
917      (component (type modulating-ion) (name ca) )
918     
919      ) ;; end SK2 current
920
921
922
923   (component (type voltage-clamp) (name SK2)
924           
925              (const vchold   = -71)
926              (const vcbase   = -69)
927              (const vcinc    = 10)
928              (const vcsteps  = 8)
929              (const vchdur   = 30)
930              (const vcbdur   = 100)
931
932              (component (type default-concentration) (name k)
933                      (const kcn    = 140)
934                      (const kcnout = 5)
935                      (output kcn kcnout))
936
937              (const cnhold   = 5e-5)
938              (const cnbase   = 5e-5)
939              (const cninc    = 1e3)
940              (const cnsteps  = 1)
941              (const cnout    = 2)
942             
943              (output vchold vcbase vcinc vcsteps vchdur vcbdur
944                      cnhold cnbase cninc cnsteps cnout)
945              )
946   
947   (component (type ionic-current) (name Na )
948
949              (input
950               (nai from ion-pools)
951               (nao from ion-pools)
952               )
953
954              (defun linoid (x y)
955                (if ((abs (x / y)) < 1e-6)
956                    then (y * (1 - ((x / y) / 2)))
957                    else (x / (1 - exp (x / y) ))
958                    ))
959             
960             
961              (component (type gate)
962             
963              ;; rate constants
964                         
965                         (Q10 = (pow (3 ((celsius - 20) / 10))))
966                         
967                         (const Aalpha_u  = 0.3 (unit /mV-ms))
968                         (const Kalpha_u  = -10 (unit mV))
969                         (const V0alpha_u = -25 (unit mV))
970       
971                         (const Abeta_u  = 12      (unit /ms))
972                         (const Kbeta_u  = -18.182 (unit mV))
973                         (const V0beta_u = -50     (unit mV))
974
975                         (const Aalpha_v  = 0.21 (unit /ms))
976                         (const Kalpha_v  = -3.333 (unit mV))
977                         (const V0alpha_v = -50    (unit mV))
978 
979                         (const Abeta_v  = 3   (unit /ms))
980                         (const Kbeta_v  = -5  (unit mV))
981                         (const V0beta_v = -17 (unit mV))
982
983                         ;; rate functions
984                         (defun alpha_u (v Q10)
985                           (Q10 * Aalpha_u * linoid((v - V0alpha_u) Kalpha_u) ))
986
987                         (defun beta_u (v Q10)
988                           (Q10 * Abeta_u * exp((v - V0beta_u) / Kbeta_u) ))
989
990                         (defun alpha_v (v Q10)
991                           (Q10 * Aalpha_v * exp((v - V0alpha_v) / Kalpha_v) ))
992                               
993                         (defun beta_v (v Q10)
994                           (Q10 * Abeta_v / (1 + exp((v - V0beta_v) / Kbeta_v) )))
995
996                         (u_inf = ((alpha_u (v Q10)) / (alpha_u (v Q10) + beta_u (v Q10)) ))
997                         (tau_u = (1 / (alpha_u (v Q10) + beta_u (v Q10))))
998
999                         (v_inf = ((alpha_v (v Q10)) / (alpha_v (v Q10) + beta_v (v Q10)) ))
1000                         (tau_v = (1 / (alpha_v (v Q10) + beta_v (v Q10)) ))
1001                         
1002                         (hh-ionic-gate
1003                          (Na  ;; ion name: exported variables will be of the form {ion}_{id}
1004                           (initial-m  (u_inf))
1005                           (initial-h  (v_inf))
1006                           (m-power    3)
1007                           (h-power    1)
1008                           (m-inf     (u_inf))
1009                           (m-tau     (tau_u))
1010                           (h-inf     (v_inf))
1011                           (h-tau     (tau_v))
1012                           ))
1013                         
1014                         )
1015             
1016              (component (type pore)
1017                         (const  gbar  = 0.048 (unit S/cm2))
1018                         (output gbar ))
1019
1020             
1021              (component (type permeating-ion) (name na)
1022                         (e = (nernst (celsius nai nao 1)) (unit mV))
1023                         (output e ))
1024             
1025              ) ;; end Na current
1026
1027
1028   (component (type voltage-clamp) (name Na)
1029
1030              (const vchold   = -71)
1031              (const vcbase   = -60)
1032              (const vcinc    = 10)
1033              (const vcsteps  = 9)
1034              (const vchdur   = 30)
1035              (const vcbdur   = 100)
1036
1037              (component (type default-concentration) (name na)
1038                      (const cn    = 5)
1039                      (const cnout = 140)
1040                      (output cn cnout))
1041           
1042              (output vchold vcbase vcinc vcsteps vchdur vcbdur))
1043
1044   
1045   (component (type ionic-current) (name NaP )
1046             
1047              (input
1048               (nai from ion-pools)
1049               (nao from ion-pools)
1050               )
1051
1052              (defun linoid (x y)
1053                (if ((abs (x / y)) < 1e-6)
1054                    then (y * (1 - ((x / y) / 2)))
1055                    else (x / (exp (x / y) - 1))
1056                    ))
1057
1058
1059              (component (type gate)
1060
1061                         ;; rate constants
1062                         ( Q10 = (pow (3 ((celsius - 30) / 10))))
1063
1064                         (const Aalpha_m  = -0.91 (unit /mV-ms))
1065                         (const Kalpha_m  = -5    (unit mV))
1066                         (const V0alpha_m = -40   (unit mV))
1067
1068                         (const Abeta_m   = 0.62 (unit /mV-ms))
1069                         (const Kbeta_m   = 5    (unit mV))
1070                         (const V0beta_m  = -40  (unit mV))
1071
1072                         (const V0_minf   = -43  (unit mV))
1073                         (const B_minf    = 5    (unit mV))
1074
1075                         ;; rate functions
1076
1077                         (defun alpha_m (v Q10)
1078                           (Q10 * Aalpha_m * linoid( ( v - V0alpha_m ) Kalpha_m) ))
1079
1080                         (defun beta_m (v Q10)
1081                           (Q10 * Abeta_m * linoid( ( v - V0beta_m )  Kbeta_m ) ))
1082
1083                         (hh-ionic-gate
1084                          (NaP  ;; ion name: exported variables will be of the form {ion}_{id}
1085                           (initial-m  ((1 / (1 + exp ((neg (v - V0_minf)) / B_minf)))))
1086                           (m-power    1)
1087                           (h-power    0)
1088                           (m-inf     ((1 / (1 + exp ((neg (v - V0_minf)) / B_minf)))))
1089                           (m-tau     ((5 / (alpha_m (v Q10) + beta_m (v Q10)) )))
1090                           ))
1091                         )
1092
1093              (component (type pore)
1094                         (const  gbar  = 0.00019 (unit S/cm2))
1095                         (output gbar ))
1096
1097              (component (type permeating-ion) (name na)
1098                         (e = (nernst (celsius nai nao 1)) (unit mV))
1099                         (output e ))
1100             
1101              ) ;; end NaP current
1102
1103
1104   (component (type voltage-clamp) (name NaP)
1105
1106              (const vchold   = -71)
1107              (const vcbase   = -60)
1108              (const vcinc    = 10)
1109              (const vcsteps  = 9)
1110              (const vchdur   = 30)
1111              (const vcbdur   = 100)
1112
1113              (component (type default-concentration) (name na)
1114                      (const cn    = 5)
1115                      (const cnout = 140)
1116                      (output cn cnout))
1117           
1118              (output vchold vcbase vcinc vcsteps vchdur vcbdur))
1119
1120   
1121   (component (type ionic-current) (name NaR )
1122             
1123              (input
1124               (nai from ion-pools)
1125               (nao from ion-pools)
1126               )
1127
1128              (component (type gate)
1129
1130                         ;; rate constants
1131                         (Q10 = (pow (3 ((celsius - 20) / 10))))
1132
1133                         (const Aalpha_s     = -0.00493 (unit /ms))
1134                         (const V0alpha_s    = -4.48754 (unit mV))
1135                         (const Kalpha_s     = -6.81881 (unit mV))
1136                         (const Shiftalpha_s = 0.00008  (unit /ms))
1137
1138                         (const Abeta_s     = 0.01558   (unit /ms))
1139                         (const V0beta_s    = 43.97494  (unit mV))
1140                         (const Kbeta_s     = 0.10818   (unit mV))
1141                         (const Shiftbeta_s = 0.04752   (unit /ms))
1142
1143                         (const Aalpha_f  = 0.31836   (unit /ms))
1144                         (const V0alpha_f = -80       (unit mV))
1145                         (const Kalpha_f  = -62.52621 (unit mV))
1146
1147                         (const Abeta_f  = 0.01014    (unit /ms))
1148                         (const V0beta_f = -83.3332   (unit mV))
1149                         (const Kbeta_f  = 16.05379   (unit mV))
1150
1151                         ;; rate functions
1152                         (defun alpha_s (v Q10)
1153                           (Q10 * (Shiftalpha_s + ((Aalpha_s * ((v + V0alpha_s) )) / (exp ((v + V0alpha_s) / Kalpha_s) - 1)))))
1154
1155                         (defun beta_s (v Q10)
1156                           (let ((x1 ((v + V0beta_s) / Kbeta_s)))
1157                             (Q10 * (Shiftbeta_s + Abeta_s * ((v + V0beta_s) / (exp((if (x1 > 200) then 200 else x1)) - 1))))))
1158
1159                         (defun alpha_f (v Q10)
1160                           (Q10 * Aalpha_f * exp( ( v - V0alpha_f ) / Kalpha_f)))
1161
1162                         (defun beta_f (v Q10)
1163                           (Q10 * Abeta_f * exp( ( v - V0beta_f ) / Kbeta_f )  ))
1164
1165                         (s_inf = ((alpha_s (v Q10)) / (alpha_s (v Q10) + beta_s (v Q10)) ))
1166                         (tau_s = (1 / (alpha_s (v Q10) + beta_s (v Q10)) ))
1167                         (f_inf = ((alpha_f (v Q10)) / (alpha_f (v Q10) + beta_f (v Q10)) ))
1168                         (tau_f = (1 / (alpha_f (v Q10) + beta_f (v Q10)) ))
1169
1170                         (hh-ionic-gate
1171                          (NaR  ;; ion name: exported variables will be of the form {ion}_{id}
1172                           (initial-m  (s_inf))
1173                           (initial-h  (f_inf))
1174                           (m-power    1)
1175                           (h-power    1)
1176                           (m-inf     (s_inf))
1177                           (m-tau     (tau_s))
1178                           (h-inf     (f_inf))
1179                           (h-tau     (tau_f))
1180
1181                           ))
1182                         )
1183
1184              (component (type pore)
1185                         (const  gbar  = 0.0017 (unit S/cm2))
1186                         (output gbar ))
1187
1188              (component (type permeating-ion) (name na)
1189                         (e = (nernst (celsius nai nao 1)) (unit mV))
1190                         (output e ))
1191             
1192              ) ;; end Nar current
1193   
1194
1195   (component (type voltage-clamp) (name NaR)
1196
1197              (const vchold   = -71)
1198              (const vcbase   = -60)
1199              (const vcinc    = 10)
1200              (const vcsteps  = 9)
1201              (const vchdur   = 30)
1202              (const vcbdur   = 100)
1203
1204              (component (type default-concentration) (name na)
1205                      (const cn    = 5)
1206                      (const cnout = 140)
1207                      (output cn cnout))
1208           
1209              (output vchold vcbase vcinc vcsteps vchdur vcbdur))
1210
1211     
1212   (component (type ionic-current) (name Lkg)
1213             
1214              (component (type pore)
1215                         (const  gbar  = 21e-6 (unit S/cm2))
1216                         (output gbar))
1217             
1218              (component (type permeating-ion) (name non-specific)
1219                         (const e = -55 (unit mV))
1220                         (output e ))
1221             
1222              ) ;; end leak current
1223   )
1224
1225   ;; Following are templates for various driver scripts used to run this model
1226   (
1227
1228    ("syn_template.hoc" () 
1229#<<EOF
1230/*******Cerebellar Golgi Cell Model **********
1231
1232Developers:    Sergio Solinas & Egidio D'Angelo
1233Code contributors:  Thierry Neius, Shyam Diwakar, Lia Forti
1234Data Analysis: Sergio Solinas
1235
1236Work Progress: April 2004 - May 2007
1237
1238Developed At:  Università Degli Studi Di Pavia
1239               Dipartimento Di Scienze Fisiologiche
1240               Pavia - Italia
1241               
1242Model Published in:
1243             Sergio M. Solinas, Lia Forti, Elisabetta Cesana,
1244             Jonathan Mapelli, Erik De Schutter and Egidio D`Angelo (2008)
1245             Computational reconstruction of pacemaking and intrinsic
1246             electroresponsiveness in cerebellar golgi cells
1247             Frontiers in Cellular Neuroscience 2:2
1248
1249
1250********************************************/
1251
1252// In this configuration the ion channels
1253// were not corrected for the Liquid Junction potential.
1254// The ion reversal potential were corrected in agreement
1255// with the voltage shift.
1256// See Table 1 Solinas et al. 2008 Frontiers Neurosci 2:2
1257
1258begintemplate GoC
1259public soma,elec,seal,dend,nclist,syn
1260public SpikeTrain, RT_Vm, EL_Vm
1261
1262create soma
1263create dend[3]
1264
1265objref syn[3]
1266objref nclist, SpikeTrain, netcon, nil
1267objref RT_Vm, EL_Vm
1268
1269proc init() {
1270
1271    nclist = new List()
1272    RT_Vm = new Vector()
1273    EL_Vm = new Vector()
1274    create soma
1275    soma {
1276        nseg = 1
1277        diam = 27 // 22 pF Dieudonne98
1278        L = 27
1279        Ra = 100 // From Roth&Hausser2000
1280        celsius = 23
1281       
1282        insert Golgi_Lkg
1283       
1284        insert Golgi_Na
1285        insert Golgi_NaR
1286        insert Golgi_NaP
1287       
1288        insert Golgi_CaHVA
1289        insert Golgi_CaLVA
1290               
1291        insert Golgi_KV
1292        insert Golgi_KM
1293        insert Golgi_KA
1294
1295        insert Golgi_KCa
1296        insert Golgi_SK2
1297       
1298        insert Golgi_HCN1
1299        insert Golgi_HCN2
1300
1301        insert Golgi_ca
1302        insert Golgi_ca2
1303        insert Golgi_na
1304        insert Golgi_k
1305       
1306        cai0_ca_ion = 50e-6
1307        ca2i0_ca2_ion = cai0_ca_ion
1308       
1309        cai = cai0_ca_ion
1310        ca2i = cai
1311        ca2o = cao
1312
1313        nai0_na_ion = 5
1314        nai = nai0_na_ion
1315        ki0_k_ion = 140
1316        ki = ki0_k_ion
1317 
1318        SpikeTrain = new Vector()
1319        netcon = new NetCon(&v(0.5),nil)
1320        netcon.threshold=-20
1321        netcon.record(SpikeTrain)
1322       
1323        RT_Vm.record(&v(0.5))
1324    }
1325   
1326    create dend[3]
1327    for i=0,2 {
1328        dend[i] {
1329            nseg = 10
1330            diam = 3
1331            L = 113
1332            Ra = 100
1333            celsius = 23
1334           
1335            insert Golgi_Lkg
1336
1337            syn[i] = new Golgi_AMPA(0)
1338        }
1339        connect dend[i](0), soma(1)     
1340    }
1341   
1342
1343}
1344
1345                   
1346endtemplate GoC
1347EOF
1348)
1349
1350
1351    ("template.hoc" () 
1352#<<EOF
1353/*******Cerebellar Golgi Cell Model **********
1354
1355Developers:    Sergio Solinas & Egidio D'Angelo
1356Code contributors:  Thierry Neius, Shyam Diwakar, Lia Forti
1357Data Analysis: Sergio Solinas
1358
1359Work Progress: April 2004 - May 2007
1360
1361Developed At:  Università Degli Studi Di Pavia
1362               Dipartimento Di Scienze Fisiologiche
1363               Pavia - Italia
1364               
1365Model Published in:
1366             Sergio M. Solinas, Lia Forti, Elisabetta Cesana,
1367             Jonathan Mapelli, Erik De Schutter and Egidio D`Angelo (2008)
1368             Computational reconstruction of pacemaking and intrinsic
1369             electroresponsiveness in cerebellar golgi cells
1370             Frontiers in Cellular Neuroscience 2:2
1371
1372
1373********************************************/
1374
1375// In this configuration the ion channels
1376// were not corrected for the Liquid Junction potential.
1377// The ion reversal potential were corrected in agreement
1378// with the voltage shift.
1379// See Table 1 Solinas et al. 2008 Frontiers Neurosci 2:2
1380
1381begintemplate GoC
1382public soma,elec,seal,dend
1383public SpikeTrain, RT_Vm, EL_Vm
1384
1385create soma
1386create dend[3]
1387
1388objref SpikeTrain, netcon, nil
1389objref RT_Vm, EL_Vm
1390
1391proc init() {
1392
1393    RT_Vm = new Vector()
1394    EL_Vm = new Vector()
1395    create soma
1396    soma {
1397        nseg = 1
1398        diam = 27 // 22 pF Dieudonne98
1399        L = 27
1400        Ra = 100 // From Roth&Hausser2000
1401        celsius = 23
1402       
1403        insert Golgi_Lkg
1404       
1405        insert Golgi_Na
1406        insert Golgi_NaR
1407        insert Golgi_NaP
1408       
1409        insert Golgi_CaHVA
1410        insert Golgi_CaLVA
1411               
1412        insert Golgi_KV
1413        insert Golgi_KM
1414        insert Golgi_KA
1415
1416        insert Golgi_KCa
1417        insert Golgi_SK2
1418       
1419        insert Golgi_HCN1
1420        insert Golgi_HCN2
1421
1422        insert Golgi_ca
1423        insert Golgi_ca2
1424        insert Golgi_na
1425        insert Golgi_k
1426       
1427        cai0_ca_ion = 50e-6
1428        ca2i0_ca2_ion = cai0_ca_ion
1429       
1430        cai = cai0_ca_ion
1431        ca2i = cai
1432        ca2o = cao
1433
1434        nai0_na_ion = 5
1435        nai = nai0_na_ion
1436        ki0_k_ion = 140
1437        ki = ki0_k_ion
1438 
1439        SpikeTrain = new Vector()
1440        netcon = new NetCon(&v(0.5),nil)
1441        netcon.threshold=-20
1442        netcon.record(SpikeTrain)
1443       
1444        RT_Vm.record(&v(0.5))
1445    }
1446   
1447    create dend[3]
1448    for i=0,2 {
1449        dend[i] {
1450            nseg = 10
1451            diam = 3
1452            L = 113
1453            Ra = 100
1454            celsius = 23
1455           
1456            insert Golgi_Lkg
1457        }
1458        connect dend[i](0), soma(1)     
1459    }
1460   
1461
1462}
1463
1464                   
1465endtemplate GoC
1466EOF
1467)
1468
1469
1470    ("reduced_template.hoc" () 
1471#<<EOF
1472/*******Cerebellar Golgi Cell Model **********
1473
1474Reduced single-compartment cerebellar Golgi cell model
1475
1476Developers:    Sergio Solinas & Egidio D'Angelo
1477Code contributors:  Thierry Neius, Shyam Diwakar, Lia Forti
1478Data Analysis: Sergio Solinas
1479
1480Work Progress: April 2004 - May 2007
1481
1482Developed At:  Università Degli Studi Di Pavia
1483               Dipartimento Di Scienze Fisiologiche
1484               Pavia - Italia
1485               
1486Model Published in:
1487             Sergio M. Solinas, Lia Forti, Elisabetta Cesana,
1488             Jonathan Mapelli, Erik De Schutter and Egidio D`Angelo (2008)
1489             Computational reconstruction of pacemaking and intrinsic
1490             electroresponsiveness in cerebellar golgi cells
1491             Frontiers in Cellular Neuroscience 2:2
1492
1493
1494********************************************/
1495
1496// In this configuration the ion channels
1497// were not corrected for the Liquid Junction potential.
1498// The ion reversal potential were corrected in agreement
1499// with the voltage shift.
1500// See Table 1 Solinas et al. 2008 Frontiers Neurosci 2:2
1501
1502begintemplate GoC
1503public soma,axon,elec,seal,dend,syn
1504public nclist, SpikeTrain, RT_Vm, EL_Vm
1505
1506create soma
1507
1508objref nclist, SpikeTrain, netcon, nil
1509objref RT_Vm, EL_Vm
1510objref syn[3]
1511
1512proc init() {
1513
1514    nclist = new List()
1515    RT_Vm = new Vector()
1516    EL_Vm = new Vector()
1517    create soma
1518    soma {
1519        nseg = 1
1520        diam = 27 // 22 pF Dieudonne98
1521        L = 27
1522        Ra = 100 // From Roth&Hausser2000
1523        celsius = 23
1524       
1525        insert Golgi_Lkg
1526       
1527        insert Golgi_Na
1528        insert Golgi_NaR
1529        insert Golgi_NaP
1530       
1531        insert Golgi_CaHVA
1532        insert Golgi_CaLVA
1533               
1534        insert Golgi_KV
1535        insert Golgi_KM
1536        insert Golgi_KA
1537
1538        insert Golgi_KCa
1539        insert Golgi_SK2
1540       
1541        insert Golgi_HCN1
1542        insert Golgi_HCN2
1543
1544        insert Golgi_ca
1545        insert Golgi_ca2
1546        insert Golgi_na
1547        insert Golgi_k
1548
1549        cai0_ca_ion = 50e-6
1550        ca2i0_ca2_ion = cai0_ca_ion
1551       
1552        cai = cai0_ca_ion
1553        ca2i = cai
1554        ca2o = cao
1555       
1556        ena =  87.39
1557        ek  = -84.69
1558
1559        for i=0,2 {
1560            syn[i] = new Golgi_AMPA(0)
1561        }
1562
1563        SpikeTrain = new Vector()
1564        netcon = new NetCon(&v(0.5),nil)
1565        netcon.threshold=-20
1566        netcon.record(SpikeTrain)
1567       
1568        RT_Vm.record(&v(0.5))
1569    }
1570}
1571
1572                   
1573endtemplate GoC
1574EOF
1575)
1576
1577    ("netstim_protocol.hoc" () 
1578#<<EOF
1579load_file("nrngui.hoc")
1580load_file("Golgi_protocol.ses")
1581
1582// Load the Golgi cell template
1583xopen("Golgi_syn_template.hoc")
1584
1585objref goc
1586goc = new GoC()
1587
1588prelength = 1000
1589mainlength = 6000
1590
1591
1592//**********************************************************************
1593proc simulate() { local preDT, mainDT, logsize  localobj logfile, tlog, Vlog, iklog, inalog, icalog, ica2log
1594   
1595    mainDT = 0.025
1596    preDT = 0.025
1597   
1598    dt = preDT
1599    tstop = prelength
1600    run()
1601   
1602   
1603    if ( stoprun == 1) return
1604   
1605    dt = mainDT
1606    continuerun(prelength + mainlength)
1607    if ( stoprun == 1) return
1608   
1609    logfile=$o1
1610    tlog=$o2
1611    Vlog=$o3
1612    iklog=$o4
1613    inalog=$o5
1614    icalog=$o6
1615    ica2log=$o7
1616   
1617    logsize = tlog.size()
1618   
1619    for i=0,tlog.size()-1 {
1620        logfile.printf("%g %g %g %g %g\n", tlog.x[i], Vlog.x[i], iklog.x[i], inalog.x[i], icalog.x[i], ica2log.x[i])
1621    }
1622
1623}
1624
1625
1626//*************User-Interface*******************************************
1627
1628nrnsecmenu(0.5, 1)
1629
1630xpanel("Spontaneous firing")
1631xvalue("Time for Initialization", "prelength")
1632xvalue("Main duration", "mainlength")
1633
1634xvalue("dt", "dt")
1635xvalue("t", "t")
1636xlabel("")
1637xbutton("Start", "simulate()")
1638xbutton("Stop", "stoprun = 1")
1639xpanel()
1640
1641vec_sizes = (prelength+mainlength)/dt + 1       // recorded traces are all this size
1642
1643strdef logfn
1644objref ns, iklog, inalog, icalog, ica2log, Vlog, tlog, logfile
1645
1646objref stim_rates // input stimulus rates
1647
1648nrates = 5
1649stim_rates = new Vector (nrates)
1650stim_rates.x[0] = 1
1651stim_rates.x[1] = 20
1652stim_rates.x[2] = 50
1653stim_rates.x[3] = 100
1654stim_rates.x[4] = 200
1655
1656for (k=0; k < nrates; k = k + 1) {
1657
1658  ns = new NetStim()
1659  ns.start    = 1000
1660  ns.number   = 100000
1661  ns.noise    = 1
1662  ns.interval = (1 / stim_rates.x[k]) * 1000
1663
1664  goc.nclist.append (new NetCon(ns,goc.syn[0],-20,1,0.1))
1665  goc.nclist.append (new NetCon(ns,goc.syn[1],-20,12,0.1))
1666  goc.nclist.append (new NetCon(ns,goc.syn[2],-20,15,0.1))
1667
1668  iklog = new Vector(vec_sizes)
1669  iklog.record (&goc.soma.ik(0.5))
1670
1671  inalog = new Vector(vec_sizes)
1672  inalog.record (&goc.soma.ina(0.5))
1673
1674  icalog = new Vector(vec_sizes)
1675  icalog.record (&goc.soma.ica(0.5))
1676
1677  ica2log = new Vector(vec_sizes)
1678  ica2log.record (&goc.soma.ica2(0.5))
1679
1680  Vlog = new Vector(vec_sizes)
1681  Vlog.record (&goc.soma.v(0.5))
1682
1683  tlog = new Vector(vec_sizes,0)
1684  tlog.record (&t)
1685
1686  sprint (logfn, "Golgi_netstim_protocol_%d_Hz.dat", stim_rates.x[k] )
1687
1688  logfile = new File()
1689  logfile.wopen ( logfn  )
1690
1691  simulate(logfile,tlog,Vlog,iklog,inalog,icalog,ica2log)
1692  logfile.close()
1693}
1694
1695quit()
1696EOF
1697)
1698
1699
1700    ("reduced_netstim_protocol.hoc" () 
1701#<<EOF
1702load_file("nrngui.hoc")
1703load_file("Golgi_protocol.ses")
1704
1705// Load the Golgi cell template
1706xopen("Golgi_reduced_template.hoc")
1707
1708objref goc
1709goc = new GoC()
1710
1711prelength = 1000
1712mainlength = 5000
1713
1714
1715//**********************************************************************
1716proc simulate() { local preDT, mainDT, logsize  localobj logfile, tlog, Vlog, iklog, inalog, icalog, ica2log
1717   
1718    mainDT = 0.025
1719    preDT = 0.025
1720   
1721    dt = preDT
1722    tstop = prelength
1723    run()
1724   
1725   
1726    if ( stoprun == 1) return
1727   
1728    dt = mainDT
1729    continuerun(prelength + mainlength)
1730    if ( stoprun == 1) return
1731   
1732    logfile=$o1
1733    tlog=$o2
1734    Vlog=$o3
1735    iklog=$o4
1736    inalog=$o5
1737    icalog=$o6
1738    ica2log=$o7
1739   
1740    logsize = tlog.size()
1741   
1742    for i=0,tlog.size()-1 {
1743        logfile.printf("%g %g %g %g %g\n", tlog.x[i], Vlog.x[i], iklog.x[i], inalog.x[i], icalog.x[i], ica2log.x[i])
1744    }
1745
1746}
1747
1748
1749//*************User-Interface*******************************************
1750
1751nrnsecmenu(0.5, 1)
1752
1753xpanel("Spontaneous firing")
1754xvalue("Time for Initialization", "prelength")
1755xvalue("Main duration", "mainlength")
1756
1757xvalue("dt", "dt")
1758xvalue("t", "t")
1759xlabel("")
1760xbutton("Start", "simulate()")
1761xbutton("Stop", "stoprun = 1")
1762xpanel()
1763
1764vec_sizes = (prelength+mainlength)/dt + 1       // recorded traces are all this size
1765
1766strdef logfn
1767objref ns, iklog, inalog, icalog, ica2log, Vlog, tlog, logfile
1768
1769objref stim_rates // input stimulus rates
1770
1771nrates = 5
1772stim_rates = new Vector (nrates)
1773stim_rates.x[0] = 1
1774stim_rates.x[1] = 20
1775stim_rates.x[2] = 50
1776stim_rates.x[3] = 100
1777stim_rates.x[4] = 200
1778
1779for (k=0; k < nrates; k = k + 1) {
1780
1781   ns = new NetStim()
1782   ns.start    = 1000
1783   ns.number   = 1000000
1784   ns.noise    = 1
1785   ns.interval = (1 / stim_rates.x[k]) * 1000
1786   ns.seed (1)
1787
1788   goc.nclist.append (new NetCon(ns,goc.syn[0],-20,1,0.1))
1789   goc.nclist.append (new NetCon(ns,goc.syn[1],-20,12,0.1))
1790   goc.nclist.append (new NetCon(ns,goc.syn[2],-20,15,0.1))
1791
1792   iklog = new Vector(vec_sizes)
1793   iklog.record (&goc.soma.ik(0.5))
1794
1795   inalog = new Vector(vec_sizes)
1796   inalog.record (&goc.soma.ina(0.5))
1797
1798   icalog = new Vector(vec_sizes)
1799   icalog.record (&goc.soma.ica(0.5))
1800
1801   ica2log = new Vector(vec_sizes)
1802   ica2log.record (&goc.soma.ica2(0.5))
1803
1804   Vlog = new Vector(vec_sizes)
1805   Vlog.record (&goc.soma.v(0.5))
1806
1807   tlog = new Vector(vec_sizes,0)
1808   tlog.record (&t)
1809
1810   sprint (logfn, "Golgi_reduced_netstim_protocol_%d_Hz.dat", stim_rates.x[k])
1811
1812   logfile = new File()
1813   logfile.wopen ( logfn )
1814
1815   simulate(logfile,tlog,Vlog,iklog,inalog,icalog,ica2log)
1816   logfile.close()
1817}
1818
1819quit()
1820EOF
1821)
1822
1823    ("protocol.hoc" () 
1824#<<EOF
1825load_file("nrngui.hoc")
1826load_file("Golgi_protocol.ses")
1827
1828// Load the Golgi cell template
1829xopen("Golgi_template.hoc")
1830
1831objref goc
1832goc = new GoC()
1833
1834
1835prelength = 1000
1836mainlength = 5000
1837
1838
1839//**********************************************************************
1840proc simulate() { local preDT, mainDT, logsize  localobj logfile, tlog, Vlog, iklog, inalog, icalog, ica2log
1841   
1842    mainDT = 0.025
1843    preDT = 0.025
1844   
1845    dt = preDT
1846    tstop = prelength
1847    run()
1848   
1849   
1850    if ( stoprun == 1) return
1851   
1852    dt = mainDT
1853    continuerun(prelength + mainlength)
1854    if ( stoprun == 1) return
1855   
1856    logfile=$o1
1857    tlog=$o2
1858    Vlog=$o3
1859    iklog=$o4
1860    inalog=$o5
1861    icalog=$o6
1862    ica2log=$o7
1863   
1864    logsize = tlog.size()
1865   
1866    for i=0,tlog.size()-1 {
1867        logfile.printf("%g %g %g %g %g\n", tlog.x[i], Vlog.x[i], iklog.x[i], inalog.x[i], icalog.x[i], ica2log.x[i])
1868    }
1869
1870}
1871
1872
1873//*************User-Interface*******************************************
1874
1875nrnsecmenu(0.5, 1)
1876
1877xpanel("Spontaneous firing")
1878xvalue("Time for Initialization", "prelength")
1879xvalue("Main duration", "mainlength")
1880
1881xvalue("dt", "dt")
1882xvalue("t", "t")
1883xlabel("")
1884xbutton("Start", "simulate()")
1885xbutton("Stop", "stoprun = 1")
1886xpanel()
1887
1888vec_sizes = (prelength+mainlength)/dt + 1       // recorded traces are all this size
1889
1890stimdur = 1000
1891
1892objectvar stim1
1893goc.soma stim1 = new IClamp(0.5)
1894stim1.del = prelength
1895stim1.dur = stimdur
1896stim1.amp = 1.
1897
1898objectvar stim2
1899goc.soma stim2 = new IClamp(0.5)
1900stim2.del = prelength+1*stimdur
1901stim2.dur = stimdur
1902stim2.amp = 2.
1903
1904objectvar stim3
1905goc.soma stim3 = new IClamp(0.5)
1906stim3.del = prelength+2*stimdur
1907stim3.dur = stimdur
1908stim3.amp = 3.
1909
1910objectvar stim4
1911goc.soma stim4 = new IClamp(0.5)
1912stim4.del = prelength+3*stimdur
1913stim4.dur = stimdur
1914stim4.amp = 4.
1915
1916objref  iklog
1917iklog = new Vector(vec_sizes)
1918iklog.record (&goc.soma.ik(0.5))
1919
1920objref  inalog
1921inalog = new Vector(vec_sizes)
1922inalog.record (&goc.soma.ina(0.5))
1923
1924objref  icalog
1925icalog = new Vector(vec_sizes)
1926icalog.record (&goc.soma.ica(0.5))
1927
1928objref  ica2log
1929ica2log = new Vector(vec_sizes)
1930ica2log.record (&goc.soma.ica2(0.5))
1931
1932objref  Vlog
1933Vlog = new Vector(vec_sizes)
1934Vlog.record (&goc.soma.v(0.5))
1935
1936objref tlog
1937tlog = new Vector(vec_sizes,0)
1938tlog.record (&t)
1939
1940objref logfile
1941logfile = new File()
1942logfile.wopen ( "Golgi_protocol.dat" )
1943
1944simulate(logfile,tlog,Vlog,iklog,inalog,icalog,ica2log)
1945logfile.close()
1946
1947
1948quit()
1949EOF
1950)
1951
1952    ("reduced_protocol.hoc" () 
1953#<<EOF
1954load_file("nrngui.hoc")
1955load_file("Golgi_protocol.ses")
1956
1957// Load the Golgi cell template
1958xopen("Golgi_reduced_template.hoc")
1959objref goc
1960goc = new GoC()
1961
1962
1963prelength = 1000
1964mainlength = 5000
1965
1966
1967//**********************************************************************
1968proc simulate() { local preDT, mainDT, logsize  localobj logfile, tlog, Vlog, iklog, inalog, icalog, ica2log
1969   
1970    mainDT = 0.025
1971    preDT = 0.025
1972   
1973    dt = preDT
1974    tstop = prelength
1975    run()
1976   
1977   
1978    if ( stoprun == 1) return
1979   
1980    dt = mainDT
1981    continuerun(prelength + mainlength)
1982    if ( stoprun == 1) return
1983   
1984    logfile=$o1
1985    tlog=$o2
1986    Vlog=$o3
1987    iklog=$o4
1988    inalog=$o5
1989    icalog=$o6
1990    ica2log=$o7
1991   
1992    logsize = tlog.size()
1993   
1994    for i=0,tlog.size()-1 {
1995        logfile.printf("%g %g %g %g %g\n", tlog.x[i], Vlog.x[i], iklog.x[i], inalog.x[i], icalog.x[i], ica2log.x[i])
1996    }
1997
1998}
1999
2000
2001//*************User-Interface*******************************************
2002
2003nrnsecmenu(0.5, 1)
2004
2005xpanel("Spontaneous firing")
2006xvalue("Time for Initialization", "prelength")
2007xvalue("Main duration", "mainlength")
2008
2009xvalue("dt", "dt")
2010xvalue("t", "t")
2011xlabel("")
2012xbutton("Start", "simulate()")
2013xbutton("Stop", "stoprun = 1")
2014xpanel()
2015
2016vec_sizes = (prelength+mainlength)/dt + 1       // recorded traces are all this size
2017
2018stimdur = 1000
2019
2020objectvar stim1
2021goc.soma stim1 = new IClamp(0.5)
2022stim1.del = prelength
2023stim1.dur = stimdur
2024stim1.amp = 1.
2025
2026objectvar stim2
2027goc.soma stim2 = new IClamp(0.5)
2028stim2.del = prelength+1*stimdur
2029stim2.dur = stimdur
2030stim2.amp = 2.
2031
2032objectvar stim3
2033goc.soma stim3 = new IClamp(0.5)
2034stim3.del = prelength+2*stimdur
2035stim3.dur = stimdur
2036stim3.amp = 3.
2037
2038objectvar stim4
2039goc.soma stim4 = new IClamp(0.5)
2040stim4.del = prelength+3*stimdur
2041stim4.dur = stimdur
2042stim4.amp = 4.
2043
2044
2045objref  iklog
2046iklog = new Vector(vec_sizes)
2047iklog.record (&goc.soma.ik(0.5))
2048
2049objref  inalog
2050inalog = new Vector(vec_sizes)
2051inalog.record (&goc.soma.ina(0.5))
2052
2053objref  icalog
2054icalog = new Vector(vec_sizes)
2055icalog.record (&goc.soma.ica(0.5))
2056
2057objref  ica2log
2058ica2log = new Vector(vec_sizes)
2059ica2log.record (&goc.soma.ica2(0.5))
2060
2061objref  Vlog
2062Vlog = new Vector(vec_sizes)
2063Vlog.record (&goc.soma.v(0.5))
2064
2065objref tlog
2066tlog = new Vector(vec_sizes,0)
2067tlog.record (&t)
2068
2069objref logfile
2070logfile = new File()
2071logfile.wopen ( "Golgi_reduced_protocol.dat" )
2072
2073simulate(logfile,tlog,Vlog,iklog,inalog,icalog,ica2log)
2074logfile.close()
2075
2076
2077quit()
2078EOF
2079)
2080    ("protocol.ses" () 
2081#<<EOF
2082{load_file("nrngui.hoc")}
2083objectvar save_window_, rvp_
2084objectvar scene_vector_[4]
2085objectvar ocbox_, ocbox_list_, scene_, scene_list_
2086{ocbox_list_ = new List()  scene_list_ = new List()}
2087{pwman_place(0,0,0)}
2088{
2089save_window_ = new Graph(0)
2090save_window_.size(0,5500,-80,50)
2091scene_vector_[2] = save_window_
2092{save_window_.view(0, -80, 7000, 130, 477, 9, 774.9, 232.3)}
2093graphList[0].append(save_window_)
2094save_window_.save_name("graphList[0].")
2095save_window_.addexpr("v(.5)", 1, 1, 0.8, 0.9, 2)
2096}
2097{doNotify()}
2098EOF
2099)
2100
2101    (".hoc" () 
2102#<<EOF
2103
2104create soma
2105access soma
2106
2107insert {{model_name}}
2108
2109soma {
2110        nseg = 1
2111
2112{% with diam = default(diam, 5.65), L = default(L, 5.65), celsius = default(celsius, 23) %}
2113        diam = {{diam}}
2114        L = {{L}}
2115        celsius = {{celsius}}
2116{% endwith %}
2117
2118        cm = 1
2119
2120        Ra = 100
2121    }
2122EOF
2123)
2124
2125    ("plot_neuron_original_iclamp.m" () 
2126#<<EOF
2127
2128CaHVA0       = load ("original/CaHVA.dat");
2129CaHVA1       = load ("nemo_nmodl_generated/CaHVA.dat");
2130
2131CaLVA0       = load ("original/CaLVA.dat");
2132CaLVA1       = load ("nemo_nmodl_generated/CaLVA.dat");
2133
2134HCN10       = load ("original/HCN1.dat");
2135HCN11       = load ("nemo_nmodl_generated/HCN1.dat");
2136
2137HCN20       = load ("original/HCN2.dat");
2138HCN21       = load ("nemo_nmodl_generated/HCN2.dat");
2139
2140KA0       = load ("original/KA.dat");
2141KA1       = load ("nemo_nmodl_generated/KA.dat");
2142
2143KCa0       = load ("original/KCa.dat");
2144KCa1       = load ("nemo_nmodl_generated/KCa.dat");
2145
2146KM0      = load ("original/KM.dat");
2147KM1      = load ("nemo_nmodl_generated/KM.dat");
2148
2149KV0      = load ("original/KV.dat");
2150KV1      = load ("nemo_nmodl_generated/KV.dat");
2151
2152SK20      = load ("original/SK2.dat");
2153SK21      = load ("nemo_nmodl_generated/SK2.dat");
2154
2155Na0      = load ("original/Na.dat");
2156Na1      = load ("nemo_nmodl_generated/Na.dat");
2157
2158NaR0      = load ("original/NaR.dat");
2159NaR1      = load ("nemo_nmodl_generated/NaR.dat");
2160
2161NaP0      = load ("original/NaP.dat");
2162NaP1      = load ("nemo_nmodl_generated/NaP.dat");
2163
2164
2165subplot(3,4,1);
2166plot(CaHVA0(:,1),CaHVA0(:,2),CaHVA1(:,1),CaHVA1(:,2),'linewidth',2);
2167title ("CaHVA current");
2168
2169subplot(3,4,2);
2170plot(CaLVA0(:,1),CaLVA0(:,2),CaLVA1(:,1),CaLVA1(:,2),'linewidth',2);
2171title ("CaLVA current");
2172
2173subplot(3,4,3);
2174plot(HCN10(:,1),HCN10(:,2),HCN11(:,1),HCN11(:,2),'linewidth',2);
2175title ("HCN1 current");
2176
2177subplot(3,4,4);
2178plot(HCN20(:,1),HCN20(:,2),HCN21(:,1),HCN21(:,2),'linewidth',2);
2179title ("HCN2 current");
2180
2181subplot(3,4,5);
2182plot(KA0(:,1),KA0(:,2),KA1(:,1),KA1(:,2),'linewidth',2);
2183title ("KA current");
2184
2185subplot(3,4,6);
2186plot(KCa0(:,1),KCa0(:,2),KCa1(:,1),KCa1(:,2),'linewidth',2);
2187title ("KCa current");
2188
2189subplot(3,4,7);
2190plot(KM0(:,1),KM0(:,2),KM1(:,1),KM1(:,2),'linewidth',2);
2191title ("KM current");
2192
2193subplot(3,4,8);
2194plot(KV0(:,1),KV0(:,2),KV1(:,1),KV1(:,2),'linewidth',2);
2195title ("KV current");
2196
2197subplot(3,4,9);
2198plot(SK20(:,1),SK20(:,2),SK21(:,1),SK21(:,2),'linewidth',2);
2199title ("SK2 current");
2200
2201subplot(3,4,10);
2202plot(Na0(:,1),Na0(:,2),Na1(:,1),Na1(:,2),'linewidth',2);
2203title ("Na current");
2204
2205subplot(3,4,11);
2206plot(NaR0(:,1),NaR0(:,2),NaR1(:,1),NaR1(:,2),'linewidth',2);
2207title ("NaR current");
2208
2209subplot(3,4,12);
2210plot(NaP0(:,1),NaP0(:,2),NaP1(:,1),NaP1(:,2),'linewidth',2);
2211title ("NaP current");
2212
2213print  ("NEURON_Original_Iclamp.eps", "-depsc");
2214EOF
2215)
2216
2217
2218
2219    ("plot_neuron_original_vclamp.m" () 
2220#<<EOF
2221
2222CaHVA0       = load ("original/CaHVA.dat");
2223CaHVA1       = load ("nemo_nmodl_generated/CaHVA.dat");
2224
2225CaLVA0       = load ("original/CaLVA.dat");
2226CaLVA1       = load ("nemo_nmodl_generated/CaLVA.dat");
2227
2228HCN10       = load ("original/HCN1.dat");
2229HCN11       = load ("nemo_nmodl_generated/HCN1.dat");
2230
2231HCN20       = load ("original/HCN2.dat");
2232HCN21       = load ("nemo_nmodl_generated/HCN2.dat");
2233
2234KA0       = load ("original/KA.dat");
2235KA1       = load ("nemo_nmodl_generated/KA.dat");
2236
2237KCa0       = load ("original/KCa.dat");
2238KCa1       = load ("nemo_nmodl_generated/KCa.dat");
2239
2240KM0      = load ("original/KM.dat");
2241KM1      = load ("nemo_nmodl_generated/KM.dat");
2242
2243KV0      = load ("original/KV.dat");
2244KV1      = load ("nemo_nmodl_generated/KV.dat");
2245
2246SK20      = load ("original/SK2.dat");
2247SK21      = load ("nemo_nmodl_generated/SK2.dat");
2248
2249Na0      = load ("original/Na.dat");
2250Na1      = load ("nemo_nmodl_generated/Na.dat");
2251
2252NaR0      = load ("original/NaR.dat");
2253NaR1      = load ("nemo_nmodl_generated/NaR.dat");
2254
2255NaP0      = load ("original/NaP.dat");
2256NaP1      = load ("nemo_nmodl_generated/NaP.dat");
2257
2258
2259subplot(3,4,1);
2260plot(CaHVA0(:,1),CaHVA0(:,2),CaHVA1(:,1),CaHVA1(:,2),'linewidth',2);
2261title ("CaHVA current");
2262
2263subplot(3,4,2);
2264plot(CaLVA0(:,1),CaLVA0(:,2),CaLVA1(:,1),CaLVA1(:,2),'linewidth',2);
2265title ("CaLVA current");
2266
2267subplot(3,4,3);
2268plot(HCN10(:,1),HCN10(:,2),HCN11(:,1),HCN11(:,2),'linewidth',2);
2269title ("HCN1 current");
2270
2271subplot(3,4,4);
2272plot(HCN20(:,1),HCN20(:,2),HCN21(:,1),HCN21(:,2),'linewidth',2);
2273title ("HCN2 current");
2274
2275subplot(3,4,5);
2276plot(KA0(:,1),KA0(:,2),KA1(:,1),KA1(:,2),'linewidth',2);
2277title ("KA current");
2278
2279subplot(3,4,6);
2280plot(KCa0(:,1),KCa0(:,2),KCa1(:,1),KCa1(:,2),'linewidth',2);
2281title ("KCa current");
2282
2283subplot(3,4,7);
2284plot(KM0(:,1),KM0(:,2),KM1(:,1),KM1(:,2),'linewidth',2);
2285title ("KM current");
2286
2287subplot(3,4,8);
2288plot(KV0(:,1),KV0(:,2),KV1(:,1),KV1(:,2),'linewidth',2);
2289title ("KV current");
2290
2291subplot(3,4,9);
2292plot(SK20(:,1),SK20(:,2),SK21(:,1),SK21(:,2),'linewidth',2);
2293title ("SK2 current");
2294
2295subplot(3,4,10);
2296plot(Na0(:,1),Na0(:,2),Na1(:,1),Na1(:,2),'linewidth',2);
2297title ("Na current");
2298
2299subplot(3,4,11);
2300plot(NaR0(:,1),NaR0(:,2),NaR1(:,1),NaR1(:,2),'linewidth',2);
2301title ("NaR current");
2302
2303subplot(3,4,12);
2304plot(NaP0(:,1),NaP0(:,2),NaP1(:,1),NaP1(:,2),'linewidth',2);
2305title ("NaP current");
2306
2307print  ("NEURON_Original_Vclamp.eps", "-depsc");
2308
2309EOF
2310)
2311
2312    ("nestmodule_bootstrap.sh" () 
2313#<<EOF
2314#!/bin/sh
2315
2316echo "Bootstrapping Golgi module..."
2317
2318if test -d autom4te.cache ; then
2319# we must remove this cache, because it
2320# may screw up things if configure is run for
2321# different platforms.
2322  echo "  -> Removing old automake cache ..."
2323  rm -rf autom4te.cache
2324fi
2325
2326echo "  -> Running aclocal ..."
2327aclocal
2328
2329echo "  -> Running libtoolize ..."
2330if [ `uname -s` = Darwin ] ; then
2331# libtoolize is glibtoolize on OSX
2332  LIBTOOLIZE=glibtoolize
2333else 
2334  LIBTOOLIZE=libtoolize
2335fi
2336
2337libtool_major=`$LIBTOOLIZE --version | head -n1 | cut -d\) -f2 | cut -d\. -f1`
2338$LIBTOOLIZE --force --copy --ltdl
2339
2340echo "  -> Re-running aclocal ..."
2341if test $libtool_major -le 2; then
2342  aclocal --force
2343else
2344  aclocal --force -I $(pwd)/libltdl/m4
2345fi
2346
2347echo "  -> Running autoconf ..."
2348autoconf
2349
2350# autoheader must run before automake
2351echo "  -> Running autoheader ..."
2352autoheader
2353
2354echo "  -> Running automake ..."
2355automake --foreign --add-missing --force-missing --copy
2356
2357echo "Done."
2358
2359EOF
2360)
2361
2362    ("nestmodule_configure.ac" () 
2363#<<EOF
2364AC_PREREQ(2.52)
2365
2366AC_INIT(golgimodule, 1.0, raikov@oist.jp)
2367
2368# These variables are exported to include/config.h
2369GOLGIMODULE_MAJOR=1
2370GOLGIMODULE_MINOR=0
2371GOLGIMODULE_PATCHLEVEL=0
2372
2373# Exporting source and build directories requires full path names.
2374# Thus we have to expand.
2375# Here, we are in top build dir, since source dir must exist, we can just
2376# move there and call pwd
2377if test "x$srcdir" = x ; then
2378  PKGSRCDIR=`pwd`
2379else
2380  PKGSRCDIR=`cd $srcdir && pwd`
2381fi
2382PKGBUILDDIR=`pwd`
2383
2384# If this is not called, install-sh will be put into .. by bootstrap.sh
2385# moritz, 06-26-06
2386AC_CONFIG_AUX_DIR(.)
2387
2388AM_INIT_AUTOMAKE(nest, $GOLGIMODULE_VERSION)
2389
2390# obtain host system type; HEP 2004-12-20
2391AC_CANONICAL_HOST
2392
2393# ------------------------------------------------------------------------
2394# Handle options
2395#
2396# NOTE: No programs/compilations must be run in this section;
2397#       otherwise CFLAGS and CXXFLAGS may take on funny default
2398#       values.
2399#       HEP 2004-12-20
2400# ------------------------------------------------------------------------
2401
2402# nest-config
2403NEST_CONFIG=`which nest-config`
2404AC_ARG_WITH(nest,[  --with-nest=script  nest-config script including path],
2405[
2406  if test "$withval" != yes; then
2407    NEST_CONFIG=$withval
2408  else
2409    AC_MSG_ERROR([--with-nest-config expects the nest-config script as argument. See README for details.])
2410  fi
2411])
2412
2413# -------------------------------------------
2414# END Handle options
2415# -------------------------------------------
2416
2417# sundials-config
2418SUNDIALS_CONFIG=`which sundials-config`
2419AC_ARG_WITH(sundials,[  --with-sundials=script  sundials-config script including path],
2420[
2421  if test "$withval" != yes; then
2422    SUNDIALS_CONFIG=$withval
2423#  else
2424#    AC_MSG_ERROR([--with-sundials-config expects the sundials-config script as argument. See README for details.])
2425  fi
2426])
2427
2428
2429# does nest-config work
2430AC_MSG_CHECKING([for nest-config ])
2431AC_CHECK_FILE($NEST_CONFIG, HAVE_NEST=yes,
2432              AC_MSG_ERROR([No usable nest-config was found. You may want to use --with-nest-config.]))
2433AC_MSG_RESULT(found)
2434
2435AC_MSG_CHECKING([for sundials-config ])
2436AC_CHECK_FILE($SUNDIALS_CONFIG, HAVE_SUNDIALS=yes,
2437              AC_MSG_WARN([No usable sundials-config was found. You may want to use --with-sundials-config.]))
2438AC_MSG_RESULT(found)
2439
2440# the following will crash if nest-config does not run
2441# careful, lines below must not break
2442AC_MSG_CHECKING([for NEST directory information ])
2443NEST_PREFIX=`$NEST_CONFIG --prefix`
2444NEST_CPPFLAGS=`$NEST_CONFIG --cflags`
2445NEST_COMPILER=`$NEST_CONFIG --compiler`
2446if test $prefix = NONE; then prefix=`$NEST_CONFIG --prefix`; fi
2447AC_MSG_RESULT($NEST_CPPFLAGS)
2448
2449
2450AC_MSG_CHECKING([for SUNDIALS preprocessor flags ])
2451SUNDIALS_CPPFLAGS="`$SUNDIALS_CONFIG -m cvode -t s -l c -s cppflags`"
2452AC_MSG_RESULT($SUNDIALS_CPPFLAGS)
2453
2454AC_MSG_CHECKING([for SUNDIALS linker options ])
2455SUNDIALS_LDFLAGS="`$SUNDIALS_CONFIG -m cvode -t s -l c -s libs` -lblas -llapack"
2456AC_MSG_RESULT($SUNDIALS_LDFLAGS)
2457
2458# Set the platform-dependent compiler flags based on the canonical
2459# host string.  These flags are placed in AM_{C,CXX}FLAGS.  If
2460# {C,CXX}FLAGS are given as environment variables, then they are
2461# appended to the set of automatically chosen flags.  After
2462# {C,CXX}FLAGS have been read out, they must be cleared, since
2463# system-dependent defaults will otherwise be placed into the
2464# Makefiles.  HEP 2004-12-20.
2465
2466# Before we can determine the proper compiler flags, we must know
2467# which compiler we are using.  Since the pertaining AC macros run the
2468# compiler and set CFLAGS, CXXFLAGS to system-dependent values, we
2469# need to save command line/enviroment settings of these variables
2470# first. AC_AIX must run before the compiler is run, so we must run it
2471# here.
2472# HEP 2004-12-21
2473
2474GOLGIMODULE_SAVE_CXXFLAGS=$CXXFLAGS
2475
2476# Must first check if we are on AIX
2477AC_AIX
2478
2479# Check for C++ compiler, looking for the same compiler
2480# used with NEST
2481AC_PROG_CXX([ $NEST_COMPILER ])
2482
2483# the following is makeshift, should have the macro set proper
2484# GOLGIMODULE_SET_CXXFLAGS
2485AM_CXXFLAGS=$GOLGIMODULE_SAVE_CXXFLAGS
2486CXXFLAGS=
2487
2488## Configure C environment
2489
2490AC_PROG_LD
2491AC_PROG_INSTALL
2492
2493AC_LIBLTDL_CONVENIENCE     ## put libltdl into a convenience library
2494AC_PROG_LIBTOOL            ## use libtool
2495AC_CONFIG_SUBDIRS(libltdl) ## also configure subdir containing libltdl
2496
2497#-- Set the language to C++
2498AC_LANG_CPLUSPLUS
2499
2500#-- Look for programs needed in the Makefile
2501AC_PROG_CXXCPP
2502AM_PROG_LIBTOOL
2503AC_PATH_PROGS([MAKE],[gmake make],[make])
2504
2505# ---------------------------------------------------------------
2506# Configure directories to be built
2507# ---------------------------------------------------------------
2508
2509PKGDATADIR=$datadir/$PACKAGE
2510PKGDOCDIR=$datadir/doc/$PACKAGE
2511
2512# set up directories from which to build help
2513# second line replaces space with colon as separator
2514HELPDIRS="$PKGSRCDIR $PKGSRCDIR/sli"
2515HELPDIRS=`echo $HELPDIRS | tr " " ":"`
2516
2517#-- Replace these variables in *.in
2518AC_SUBST(HAVE_NEST)
2519AC_SUBST(NEST_CONFIG)
2520AC_SUBST(NEST_CPPFLAGS)
2521AC_SUBST(NEST_COMPILER)
2522AC_SUBST(NEST_PREFIX)
2523AC_SUBST(HELPDIRS)
2524AC_SUBST(PKGSRCDIR)
2525AC_SUBST(PKGBUILDDIR)
2526AC_SUBST(PKGDATADIR)
2527AC_SUBST(PKGDOCDIR)
2528AC_SUBST(KERNEL)
2529AC_SUBST(HOST)
2530AC_SUBST(SED)
2531AC_SUBST(LD)
2532AC_SUBST(host_os)
2533AC_SUBST(host_cpu)
2534AC_SUBST(host_vendor)
2535AC_SUBST(AS)
2536AC_SUBST(CXX)
2537AC_SUBST(AR)
2538AC_SUBST(ARFLAGS)
2539AC_SUBST(CXX_AR)
2540AC_SUBST(AM_CXXFLAGS)
2541AC_SUBST(AM_CFLAGS)
2542AC_SUBST(MAKE)
2543AC_SUBST(MAKE_FLAGS)
2544AC_SUBST(INCLTDL)
2545AC_SUBST(LIBLTDL)
2546AC_SUBST(SUNDIALS_CONFIG)
2547AC_SUBST(SUNDIALS_CPPFLAGS)
2548AC_SUBST(SUNDIALS_LDFLAGS)
2549
2550AM_CONFIG_HEADER(golgimodule_config.h:golgimodule_config.h.in)
2551AC_CONFIG_FILES(Makefile)
2552
2553# -----------------------------------------------
2554# Create output
2555# -----------------------------------------------
2556AC_OUTPUT
2557
2558
2559# -----------------------------------------------
2560# Report, after output at end of configure run
2561# Must come after AC_OUTPUT, so that it is
2562# displayed after libltdl has been configured
2563# -----------------------------------------------
2564
2565echo
2566echo "-------------------------------------------------------"
2567echo "Golgi module Configuration Summary"
2568echo "-------------------------------------------------------"
2569echo
2570echo "C++ compiler        : $CXX"
2571echo "C++ compiler flags  : $AM_CXXFLAGS"
2572echo "NEST compiler flags : $NEST_CPPFLAGS"
2573echo "SUNDIALS compiler flags : $SUNDIALS_CPPFLAGS"
2574echo "SUNDIALS linker flags : $SUNDIALS_LDFLAGS"
2575
2576# these variables will still contain '${prefix}'
2577# we want to have the versions where this is resolved, too:
2578eval eval eval  PKGDOCDIR_AS_CONFIGURED=$PKGDOCDIR
2579eval eval eval  PKGDATADIR_AS_CONFIGURED=$PKGDATADIR
2580
2581echo
2582echo "-------------------------------------------------------"
2583echo
2584echo "You can build and install Golgi module now, using"
2585echo "  make"
2586echo "  make install"
2587echo
2588echo "Golgi module will be installed to:"
2589echo -n "  "; eval eval echo "$libdir"
2590echo
2591
2592EOF
2593)
2594
2595    ("nestmodule_makefile.am" () 
2596#<<EOF
2597
2598# Automake file for external dynamic modules for NEST
2599#
2600# Hans Ekkehard Plesser, April 2008
2601# Automake file for the Developer Module
2602#
2603# libgolgimodule is built as a normal, installable library.
2604# It will be installed to $prefix/lib by make install.
2605#
2606# Headers from this directory are not to be installed upon
2607# make install. They are therefore included in _SOURCES.
2608
2609
2610libdir= @libdir@/nest
2611
2612lib_LTLIBRARIES=      golgimodule.la libgolgimodule.la
2613
2614golgimodule_la_CXXFLAGS= @AM_CXXFLAGS@
2615golgimodule_la_SOURCES=  golgimodule.cpp      golgimodule.h      \
2616                      Golgi.cpp Golgi.h
2617golgimodule_la_LDFLAGS=  -module
2618
2619libgolgimodule_la_CXXFLAGS= $(golgimodule_la_CXXFLAGS) -DLINKED_MODULE
2620libgolgimodule_la_SOURCES=  $(golgimodule_la_SOURCES)
2621
2622MAKEFLAGS= @MAKE_FLAGS@
2623
2624AM_CPPFLAGS= @NEST_CPPFLAGS@ \
2625             @SUNDIALS_CPPFLAGS@ \
2626             @INCLTDL@     
2627
2628AM_LDFLAGS = @SUNDIALS_LDFLAGS@
2629
2630.PHONY: install-slidoc
2631
2632nobase_pkgdata_DATA=\
2633        golgimodule.sli
2634
2635install-slidoc:
2636        NESTRCFILENAME=/dev/null $(DESTDIR)$(NEST_PREFIX)/bin/sli --userargs="@HELPDIRS@" $(NEST_PREFIX)/share/nest/sli/install-help.sli
2637
2638install-data-hook: install-exec install-slidoc
2639
2640EXTRA_DIST= sli
2641
2642EOF
2643)
2644
2645    ("nestmodule.cpp" () 
2646#<<EOF
2647/*
2648 *  golgimodule.cpp
2649 *  This file is part of NEST.
2650 *
2651 *  Copyright (C) 2008 by
2652 *  The NEST Initiative
2653 *
2654 *  See the file AUTHORS for details.
2655 *
2656 *  Permission is granted to compile and modify
2657 *  this file for non-commercial use.
2658 *  See the file LICENSE for details.
2659 *
2660 */
2661
2662// include necessary NEST headers
2663//#include "config.h"
2664#include "network.h"
2665#include "model.h"
2666#include "dynamicloader.h"
2667#include "genericmodel.h"
2668#include "generic_connector.h"
2669#include "booldatum.h"
2670#include "integerdatum.h"
2671#include "tokenarray.h"
2672#include "exceptions.h"
2673#include "sliexceptions.h"
2674#include "nestmodule.h"
2675
2676// include headers with your own stuff
2677#include "golgimodule.h"
2678#include "Golgi.h"
2679
2680// -- Interface to dynamic module loader ---------------------------------------
2681
2682/*
2683 * The dynamic module loader must be able to find your module.
2684 * You make the module known to the loader by defining an instance of your
2685 * module class in global scope. This instance must have the name
2686 *
2687 * <modulename>_LTX_mod
2688 *
2689 * The dynamicloader can then load modulename and search for symbol "mod" in it.
2690 */
2691 
2692golginest::GolgiModule golgimodule_LTX_mod;
2693
2694// -- DynModule functions ------------------------------------------------------
2695
2696golginest::GolgiModule::GolgiModule()
2697  {
2698#ifdef LINKED_MODULE
2699     // register this module at the dynamic loader
2700     // this is needed to allow for linking in this module at compile time
2701     // all registered modules will be initialized by the main app's dynamic loader
2702     nest::DynamicLoaderModule::registerLinkedModule(this);
2703#endif     
2704   }
2705
2706golginest::GolgiModule::~GolgiModule()
2707   {
2708   }
2709
2710   const std::string golginest::GolgiModule::name(void) const
2711   {
2712     return std::string("Golgi Module"); // Return name of the module
2713   }
2714
2715   const std::string golginest::GolgiModule::commandstring(void) const
2716   {
2717     /* 1. Tell interpreter that we provide the C++ part of GolgiModule with the
2718           current revision number.
2719        2. Instruct the interpreter to check that golgimodule-init.sli exists,
2720           provides at least version 1.0 of the SLI interface to GolgiModule, and
2721           to load it.
2722      */
2723     return std::string(
2724       "/golgimodule /C++ ($Revision: 8512 $) provide-component "
2725       "/golgimodule /SLI (7165) require-component"
2726       );
2727   }
2728
2729   /* BeginDocumentation
2730      Name: StepPatternConnect - Connect sources and targets with a stepping pattern
2731     
2732      Synopsis:
2733      [sources] source_step [targets] target_step synmod StepPatternConnect -> n_connections
2734     
2735      Parameters:
2736      [sources]     - Array containing GIDs of potential source neurons
2737      source_step   - Make connection from every source_step'th neuron
2738      [targets]     - Array containing GIDs of potential target neurons
2739      target_step   - Make connection to every target_step'th neuron
2740      synmod        - The synapse model to use (literal, must be key in synapsedict)
2741      n_connections - Number of connections made
2742     
2743      Description:
2744      This function subsamples the source and target arrays given with steps
2745      source_step and target_step, beginning with the first element in each array,
2746      and connects the selected nodes.
2747     
2748      Example:
2749      /first_src 0 /network_size get def
2750      /last_src /iaf_neuron 20 Create def  % nodes  1 .. 20
2751      /src [first_src last_src] Range def
2752      /last_tgt /iaf_neuron 10 Create def  % nodes 21 .. 30
2753      /tgt [last_src 1 add last_tgt] Range def
2754     
2755      src 6 tgt 4 /drop_odd_spike StepPatternConnect
2756 
2757      This connects nodes [1, 7, 13, 19] as sources to nodes [21, 25,
2758      29] as targets using synapses of type drop_odd_spike, and
2759      returning 12 as the number of connections.  The following
2760      command will print the connections (you must paste the SLI
2761      command as one long line):
2762
2763      src { /s Set << /source s /synapse_type /static_synapse >> FindConnections { GetStatus /target get } Map dup length 0 gt { cout s <- ( -> ) <- exch <-- endl } if ; } forall
2764      1 -> [21 25 29]
2765      7 -> [21 25 29]
2766      13 -> [21 25 29]
2767      19 -> [21 25 29]
2768     
2769      Remark:
2770      This function is only provided as an example for how to write your own
2771      interface function.
2772     
2773      Author:
2774      Hans Ekkehard Plesser
2775     
2776      SeeAlso:
2777      Connect, ConvergentConnect, DivergentConnect
2778   */
2779   void golginest::GolgiModule::StepPatternConnect_Vi_i_Vi_i_lFunction::execute(SLIInterpreter *i) const
2780   {
2781     // Check if we have (at least) five arguments on the stack.
2782     i->assert_stack_load(5);
2783
2784     // Retrieve source, source step, target, target step from the stack
2785     const TokenArray sources = getValue<TokenArray> (i->OStack.pick(4)); // bottom
2786     const long src_step      = getValue<long>       (i->OStack.pick(3));
2787     const TokenArray targets = getValue<TokenArray> (i->OStack.pick(2));
2788     const long tgt_step      = getValue<long>       (i->OStack.pick(1)); 
2789     const Name synmodel_name = getValue<std::string>(i->OStack.pick(0)); // top
2790     
2791     // Obtain synapse model index
2792     const Token synmodel
2793       = nest::NestModule::get_network().get_synapsedict().lookup(synmodel_name);
2794     if ( synmodel.empty() )
2795       throw nest::UnknownSynapseType(synmodel_name.toString());
2796     const nest::index synmodel_id = static_cast<nest::index>(synmodel);
2797
2798     // Build a list of targets with the given step
2799     TokenArray selected_targets;
2800     for ( size_t t = 0 ; t < targets.size() ; t += tgt_step )
2801       selected_targets.push_back(targets[t]);
2802     
2803     // Now connect all appropriate sources to this list of targets
2804     size_t Nconn = 0;  // counts connections
2805     for ( size_t s = 0 ; s < sources.size() ; s += src_step )
2806     {
2807       // We must first obtain the GID of the source as integer
2808       const nest::long_t sgid = getValue<nest::long_t>(sources[s]);
2809
2810       // nest::network::divergent_connect() requires weight and delay arrays. We want to use
2811       // default values from the synapse model, so we pass empty arrays.
2812       nest::NestModule::get_network().divergent_connect(sgid, selected_targets,
2813                                                         TokenArray(), TokenArray(),
2814                                                         synmodel_id);
2815       Nconn += selected_targets.size();
2816     }
2817
2818     // We get here only if none of the operations above throws and exception.
2819     // Now we can safely remove the arguments from the stack and push Nconn
2820     // as our result.
2821     i->OStack.pop(5);
2822     i->OStack.push(Nconn);
2823     
2824     // Finally, we pop the call to this functions from the execution stack.
2825     i->EStack.pop();
2826   }
2827
2828  //-------------------------------------------------------------------------------------
2829
2830  void golginest::GolgiModule::init(SLIInterpreter *i, nest::Network*)
2831  {
2832    /* Register a neuron or device model.
2833       Give node type as template argument and the name as second argument.
2834       The first argument is always a reference to the network.
2835       Return value is a handle for later unregistration.
2836    */
2837       nest::register_model<nest::Golgi>(nest::NestModule::get_network(),
2838                                            "Golgi");
2839
2840  }  // GolgiModule::init()
2841
2842 
2843
2844EOF
2845)
2846
2847    ("nestmodule.h" () 
2848#<<EOF
2849/*
2850 *  golgimodule.h
2851 *
2852 *  This file is part of NEST.
2853 *
2854 *  Copyright (C) 2008 by
2855 *  The NEST Initiative
2856 *
2857 *  See the file AUTHORS for details.
2858 *
2859 *  Permission is granted to compile and modify
2860 *  this file for non-commercial use.
2861 *  See the file LICENSE for details.
2862 *
2863 */
2864
2865#ifndef GOLGIMODULE_H
2866#define GolgiMODULE_H
2867
2868#include "dynmodule.h"
2869#include "slifunction.h"
2870
2871namespace nest
2872{
2873  class Network;
2874}
2875
2876// Put your stuff into your own namespace.
2877namespace golginest {
2878 
2879/**
2880 * Class defining your model.
2881 * @note For each model, you must define one such class, with a unique name.
2882 */
2883class GolgiModule : public DynModule
2884{
2885public:
2886
2887  // Interface functions ------------------------------------------
2888 
2889  /**
2890   * @note The constructor registers the module with the dynamic loader.
2891   *       Initialization proper is performed by the init() method.
2892   */
2893  GolgiModule();
2894 
2895  /**
2896   * @note The destructor does not do much in modules. Proper "downrigging"
2897   *       is the responsibility of the unregister() method.
2898   */
2899  ~GolgiModule();
2900
2901  /**
2902   * Initialize module by registering models with the network.
2903   * @param SLIInterpreter* SLI interpreter
2904   * @param nest::Network*  Network with which to register models
2905   * @note  Parameter Network is needed for historical compatibility
2906   *        only.
2907   */
2908  void init(SLIInterpreter*, nest::Network*);
2909
2910  /**
2911   * Return the name of your model.
2912   */
2913  const std::string name(void) const;
2914 
2915  /**
2916   * Return the name of a sli file to execute when golgimodule is loaded.
2917   * This mechanism can be used to define SLI commands associated with your
2918   * module, in particular, set up type tries for functions you have defined.
2919   */
2920  const std::string commandstring(void) const;
2921     
2922public:
2923 
2924  // Classes implementing your functions -----------------------------
2925 
2926  /**
2927   * Implement a function for a step-pattern-based connection.
2928   * @note What this function does is described in the SLI documentation
2929   *       in the cpp file.
2930   * @note The mangled name indicates this function expects the following
2931   *       arguments on the stack (bottom first): vector of int, int,
2932   *       vector of int, int.
2933   * @note You must define a member object in your module class
2934   *       of the function class. execute() is later invoked on this
2935   *       member.
2936   */
2937  class StepPatternConnect_Vi_i_Vi_i_lFunction: public SLIFunction
2938     {
2939     public:
2940       void execute(SLIInterpreter *) const;
2941     };
2942
2943     StepPatternConnect_Vi_i_Vi_i_lFunction stepPatternConnect_Vi_i_Vi_i_lFunction;
2944  };
2945} // namespace golginest
2946
2947#endif
2948
2949EOF
2950)
2951
2952    ("golgimodule.sli" () 
2953#<<EOF
2954/*
2955 * Initialization file for GolgiModule.
2956 * Run automatically when GolgiModule is loaded.
2957 */
2958
2959M_DEBUG (golgimodule.sli) (Initializing SLI support for GolgiModule.) message
2960
2961/golgimodule /SLI ($Revision: 7918 $) provide-component
2962/golgimodule /C++ (7165) require-component
2963
2964/StepPatternConnect [ /arraytype /integertype /arraytype /integertype /literaltype ]
2965{
2966  StepPatternConnect_Vi_i_Vi_i_l
2967} def
2968
2969EOF
2970)
2971
2972    ("testrun.sli" () 
2973#<<EOF
2974
2975/dt 0.025 def
2976
2977ResetKernel
29780
2979 <<
2980   /resolution  dt
2981 >> SetStatus
2982
2983(golgimodule) Install
2984/neuron /Golgi Create def
2985/neuron_params << /V_m -65.0 >> def
2986neuron neuron_params SetStatus
2987
2988/stepinput /step_current_generator Create def
2989/vlog /voltmeter Create def
2990/vlog_parameters << /interval dt /to_file true /to_memory false >> def
2991vlog vlog_parameters SetStatus
2992
2993stepinput neuron Connect
2994vlog neuron Connect
2995
2996/step_current_parameters << /amplitude_times [0.0 1000.0 2000.0 3000.0 4000.0 5000.0] /amplitude_values [0.0 1000.0 2000.0 3000.0 4000.0 0.0 ] >> def
2997stepinput step_current_parameters SetStatus
2998
29996000.0 Simulate
3000
3001EOF
3002)
3003
3004    ("netstim.sli" () 
3005#<<EOF
3006
3007(golgimodule) Install
3008
3009[ 1.0 20.0 50.0 100.0 200.0 ]
3010
3011{
3012
3013/i exch def
3014
3015(Golgi_netstim_) i cvs (_Hz) join join /label exch def
3016
3017/dt 0.025 def
3018
3019/seeds 0 /total_num_virtual_procs get array 1 add def
3020
3021ResetKernel
30220
3023 <<
3024   /resolution  dt
3025   /rng_seeds seeds
3026 >> SetStatus
3027
3028/neuron /Golgi Create def
3029/neuron_params << /V_m -65.0 >> def
3030neuron neuron_params SetStatus
3031
3032/input1 /poisson_generator Create def
3033/input2 /poisson_generator Create def
3034/input3 /poisson_generator Create def
3035
3036/vlog /voltmeter Create def
3037
3038/poisson_parameters1 << /rate i Hz /start 1000.0 >> def
3039input1 poisson_parameters1 SetStatus
3040
3041/poisson_parameters2 << /rate i Hz /start 1000.0 >> def
3042input2 poisson_parameters2 SetStatus
3043
3044/poisson_parameters3 << /rate i Hz /start 1000.0 >> def
3045input3 poisson_parameters3 SetStatus
3046
3047/vlog_parameters << /interval dt /label label /to_file true  /to_memory false >> def
3048vlog vlog_parameters SetStatus
3049
3050/conn_parameters1 << /receptor_type 1 /weight 0.1 /delay 1.0 >> def
3051input1 neuron conn_parameters1 Connect
3052
3053/conn_parameters2 << /receptor_type 1 /weight 0.1 /delay 12.0 >> def
3054input2 neuron conn_parameters2 Connect
3055
3056/conn_parameters3 << /receptor_type 1 /weight 0.1 /delay 15.0 >> def
3057input3 neuron conn_parameters3 Connect
3058
3059vlog neuron Connect
3060
30616000.0 Simulate
3062
3063i
3064}
3065Map
3066EOF
3067)
3068
3069))
Note: See TracBrowser for help on using the repository browser.