rectified model for stimulus generating cell for NetPyNE

tvb_multiscale/tvb_netpyne/netpyne/mod/dynvecstim.mod

@@ -0,0 +1,121 @@
+:  Vector stream of events
+
+NEURON {
+	THREADSAFE
+	ARTIFICIAL_CELL DynamicVecStim
+	POINTER ptr
+}
+
+ASSIGNED {
+	index
+	intervalEnd
+	etime (ms)
+	ptr
+}
+
+
+INITIAL {
+	index = 0
+	inferNextEvent()
+	if (index > 0) {
+		net_send(etime - t, 1) : schedule spike at time etime
+	} else {
+		: no spikes. wait for inetrval end to be able to initialize next interval
+		net_send(intervalEnd - t, 3)
+	}
+}
+
+NET_RECEIVE (w) {
+	: flag 1 - emit a spike and prepare for the next one
+	: flag 2 - interval just started, prepare for the first spike
+	: flag 3 - interval ended, waiting for initialization of the next one
+	if (flag == 1 || flag == 2) { 
+		if (flag == 1) {
+			net_event(t) : emit spike
+		}
+		inferNextEvent()
+		if (index > 0) {
+			if (etime > t) {
+				net_send(etime - t, 1) : schedule next spike
+			} else {
+				printf("WARNING in DynVecStim. Spike at %f skipped. Make sure it doesn't go before previous interval ends.\n", etime)
+				net_send(0, 2)
+			}
+
+		} else {
+			: no more spikes. wait for inetrval end to be able to initialize next interval
+			LOCAL endTime
+			endTime = intervalEnd - t
+			intervalEnd = -1
+			net_send(endTime, 3)
+		}
+	} else if (flag == 3) {
+		if (intervalEnd == -1) {
+			: still waiting for next interval initialization
+			net_send(dt, 3)
+		} else {
+			: next interval is already initialized ( from play() )! Proceed to first spike in it
+			net_send(0, 2)
+		}
+	}
+}
+
+DESTRUCTOR {
+VERBATIM
+	void* vv = (void*)(_p_ptr);  
+        if (vv) {
+		hoc_obj_unref(*vector_pobj(vv));
+	}
+ENDVERBATIM
+}
+
+PROCEDURE inferNextEvent() {
+VERBATIM	
+  { void* vv; int i, size; double* px;
+	i = (int)index;
+	if (i >= 0) {
+		vv = (void*)(_p_ptr);
+		if (vv) {
+			size = vector_capacity(vv);
+			px = vector_vec(vv);
+			if (i < size) {
+				etime = px[i];
+				if (etime < intervalEnd) {
+					index += 1.;
+				} else {
+					printf("WARNING in DynVecStim. Spike at %f will be skipped, as well as all later spikes in this interval, since they go after the interval end\n", etime);
+					index = -1.;
+				}
+			}else{
+				index = -1.;
+			}
+		}else{
+			index = -1.;
+		}
+	}
+  }
+ENDVERBATIM
+}
+
+PROCEDURE play() {
+	index = 0
+VERBATIM
+	void** pv;
+	void* ptmp = NULL;
+	if (ifarg(2)) {
+		intervalEnd = *getarg(2);
+	} else {
+		printf("ERROR in DynVecStim! End of interval should be specified as second arg in play() !!");
+		exit(1);
+	}
+	if (ifarg(1)) {
+		ptmp = vector_arg(1);
+		hoc_obj_ref(*vector_pobj(ptmp));
+	}
+	pv = (void**)(&_p_ptr);
+	if (*pv) {
+		hoc_obj_unref(*vector_pobj(*pv));
+	}
+	*pv = ptmp;
+ENDVERBATIM
+}

tvb_multiscale/tvb_netpyne/netpyne/module.py

@@ -14,14 +14,13 @@ def __init__(self):
 
         # Make sure that all required mod-files are compiled (is there a better way to check?)
         try:
-            h.DynamicNetStim()
+            h.DynamicVecStim()
         except:
             print("NetPyNE couldn't find necessary MOD-files. Trying to compile..")
-            import sys
-            import os
+            import sys, os
             python_path = sys.executable.split("python")[0]
             tvb_multiscale_path = os.path.abspath(__file__).split("tvb_multiscale")[0]
-            os.system('%snrnivmodl %s/tvb_multiscale/tvb_netpyne/netpyne/mod' % (python_path, tvb_multiscale_path))
+            os.system(f'{python_path}nrnivmodl {tvb_multiscale_path}/tvb_multiscale/tvb_netpyne/netpyne/mod')
             h.nrn_load_dll('./x86_64/libnrnmech.so')
 
     def importModel(self, netParams, simConfig, dt, config):
@@ -35,8 +34,8 @@ def importModel(self, netParams, simConfig, dt, config):
         self.netParams = netParams
         self.simConfig = simConfig
 
-        # using DynamicNetStim model for artificial cells serving as stimuli
-        self.netParams.cellParams['art_NetStim'] = {'cellModel': 'DynamicNetStim'}
+        # using DynamicVecStim model for artificial cells serving as stimuli
+        self.netParams.cellParams['art_NetStim'] = {'cellModel': 'DynamicVecStim'}
 
     @property
     def dt(self):
@@ -202,13 +201,9 @@ def stimulate(self, length):
 
         intervalEnd = h.t + length
         for gid in allNeurons:
-            # currently used .mod implementation allows no more then 3 spikes during the interval.
-            # if for given cell they are less than 3, use -1 for the rest. If they are more, the rest will be lost. But for the reasonable spiking rates, this latter case is highly unlikely. 
             spikes = allNeuronsSpikes.get(gid, [])
-            spks = [-1] * 3
-            for i, spike in enumerate(spikes[:3]):
-                spks[i] = spike
-            sim.net.cells[gid].hPointp.set_next_spikes(intervalEnd, spks[0], spks[1], spks[2])
+            spikes = h.Vector(spikes)
+            sim.net.cells[gid].hPointp.play(spikes, intervalEnd)
 
     def finalize(self):
         if self.time < sim.cfg.duration: