Feature/devices

.gitignore

@@ -2,6 +2,7 @@
 __pycache__/
 *.py[cod]
 *$py.class
+*.h5
 
 # C extensions
 *.so

bsb_nest/device.py

@@ -14,10 +14,13 @@ class NestRule:
 @config.dynamic(attr_name="device", auto_classmap=True, default="external")
 class NestDevice(DeviceModel):
     weight = config.attr(type=float, required=True)
+    """weight of the connection between the device and its target"""
     delay = config.attr(type=float, required=True)
+    """delay of the transmission between the device and its target"""
     targetting = config.attr(
         type=types.or_(Targetting, NestRule), default=dict, call_default=True
     )
+    """Targets of the device, which should be either a population or a nest rule"""
 
     def get_target_nodes(self, adapter, simulation, simdata):
         if isinstance(self.targetting, Targetting):

bsb_nest/devices/__init__.py

@@ -1,2 +1,4 @@
+from .dc_generator import DCGenerator
+from .multimeter import Multimeter
 from .poisson_generator import PoissonGenerator
 from .spike_recorder import SpikeRecorder

bsb_nest/devices/dc_generator.py

@@ -0,0 +1,23 @@
+import nest
+from bsb import config
+
+from ..device import NestDevice
+
+
+@config.node
+class DCGenerator(NestDevice, classmap_entry="dc_generator"):
+    amplitude = config.attr(type=float, required=True)
+    """Current amplitude of the dc generator"""
+    start = config.attr(type=float, required=False, default=0.0)
+    """Activation time in ms"""
+    stop = config.attr(type=float, required=False, default=None)
+    """Deactivation time in ms. 
+        If not specified, generator will last until the end of the simulation."""
+
+    def implement(self, adapter, simulation, simdata):
+        nodes = self.get_target_nodes(adapter, simulation, simdata)
+        params = {"amplitude": self.amplitude, "start": self.start}
+        if self.stop is not None and self.stop > self.start:
+            params["stop"] = self.stop
+        device = self.register_device(simdata, nest.Create("dc_generator", params=params))
+        self.connect_to_nodes(device, nodes)

bsb_nest/devices/multimeter.py

@@ -0,0 +1,52 @@
+import nest
+import quantities as pq
+from bsb import ConfigurationError, _util, config, types
+from neo import AnalogSignal
+
+from ..device import NestDevice
+
+
+@config.node
+class Multimeter(NestDevice, classmap_entry="multimeter"):
+    weight = config.provide(1)
+    properties: list[str] = config.attr(type=types.list(str))
+    """List of properties to record in the Nest model."""
+    units: list[str] = config.attr(type=types.list(str))
+    """List of properties' units."""
+
+    def boot(self):
+        _util.assert_samelen(self.properties, self.units)
+        for i in range(len(self.units)):
+            if not self.units[i] in pq.units.__dict__.keys():
+                raise ConfigurationError(
+                    f"Unit {self.units[i]} not in the list of known units of quantities"
+                )
+
+    def implement(self, adapter, simulation, simdata):
+
+        nodes = self.get_target_nodes(adapter, simulation, simdata)
+        device = self.register_device(
+            simdata,
+            nest.Create(
+                "multimeter",
+                params={
+                    "interval": self.simulation.resolution,
+                    "record_from": self.properties,
+                },
+            ),
+        )
+        self.connect_to_nodes(device, nodes)
+
+        def recorder(segment):
+            for prop, unit in zip(self.properties, self.units):
+                segment.analogsignals.append(
+                    AnalogSignal(
+                        device.events[prop],
+                        units=pq.units.__dict__[unit],
+                        sampling_period=self.simulation.resolution * pq.ms,
+                        name=self.name,
+                        senders=device.events["senders"],
+                    )
+                )
+
+        simdata.result.create_recorder(recorder)

bsb_nest/devices/poisson_generator.py

@@ -8,11 +8,20 @@
 @config.node
 class PoissonGenerator(NestDevice, classmap_entry="poisson_generator"):
     rate = config.attr(type=float, required=True)
+    """Frequency of the poisson generator"""
+    start = config.attr(type=float, required=False, default=0.0)
+    """Activation time in ms"""
+    stop = config.attr(type=float, required=False, default=None)
+    """Deactivation time in ms. 
+        If not specified, generator will last until the end of the simulation."""
 
     def implement(self, adapter, simulation, simdata):
         nodes = self.get_target_nodes(adapter, simulation, simdata)
+        params = {"rate": self.rate, "start": self.start}
+        if self.stop is not None and self.stop > self.start:
+            params["stop"] = self.stop
         device = self.register_device(
-            simdata, nest.Create("poisson_generator", params={"rate": self.rate})
+            simdata, nest.Create("poisson_generator", params=params)
         )
         sr = nest.Create("spike_recorder")
         nest.Connect(device, sr)

tests/test_nest.py

@@ -2,7 +2,7 @@
 
 import nest
 import numpy as np
-from bsb import CastError
+from bsb import BootError, CastError, ConfigurationError
 from bsb.config import Configuration
 from bsb.core import Scaffold
 from bsb.services import MPI
@@ -162,6 +162,8 @@ def test_iaf_cond_alpha(self):
 
         spike_times_nest = spikeA.get("events")["times"]
 
+        duration = 1000
+        resolution = 0.1
         cfg = Configuration(
             {
                 "name": "test",
@@ -182,8 +184,8 @@ def test_iaf_cond_alpha(self):
                 "simulations": {
                     "test": {
                         "simulator": "nest",
-                        "duration": 1000,
-                        "resolution": 0.1,
+                        "duration": duration,
+                        "resolution": resolution,
                         "cell_models": {
                             "A": {
                                 "model": "iaf_cond_alpha",
@@ -199,7 +201,17 @@ def test_iaf_cond_alpha(self):
                                     "strategy": "cell_model",
                                     "cell_models": ["A"],
                                 },
-                            }
+                            },
+                            "voltmeter_A": {
+                                "device": "multimeter",
+                                "delay": resolution,
+                                "properties": ["V_m"],
+                                "units": ["mV"],
+                                "targetting": {
+                                    "strategy": "cell_model",
+                                    "cell_models": ["A"],
+                                },
+                            },
                         },
                     }
                 },
@@ -210,8 +222,130 @@ def test_iaf_cond_alpha(self):
         netw.compile()
         results = netw.run_simulation("test")
         spike_times_bsb = results.spiketrains[0]
+        self.assertTrue(np.unique(spike_times_bsb.annotations["senders"]) == 1)
+        membrane_potentials = results.analogsignals[0]
+        # last time point is not recorded because of recorder delay.
+        self.assertTrue(len(membrane_potentials) == duration / resolution - 1)
+        self.assertTrue(np.unique(membrane_potentials.annotations["senders"]) == 1)
+        defaults = nest.GetDefaults("iaf_cond_alpha")
+        # since current injected is positive, the V_m should be clamped between default
+        # initial V_m = -70mV and spike threshold V_th = -55 mV
+        self.assertAll(
+            (membrane_potentials <= defaults["V_th"])
+            * (membrane_potentials >= defaults["V_m"])
+        )
         self.assertClose(np.array(spike_times_nest), np.array(spike_times_bsb))
 
+    def test_multimeter_errors(self):
+        cfg = get_test_config("gif_pop_psc_exp")
+        sim_cfg = cfg.simulations.test_nest
+        sim_cfg.devices.update(
+            {
+                "voltmeter": {
+                    "device": "multimeter",
+                    "delay": 0.1,
+                    "properties": ["V_m", "I_syn"],
+                    "units": ["mV"],
+                    "targetting": {
+                        "strategy": "cell_model",
+                        "cell_models": ["gif_pop_psc_exp"],
+                    },
+                },
+            }
+        )
+        with self.assertRaises(BootError):
+            Scaffold(cfg, self.storage)
+
+        sim_cfg.devices.update(
+            {
+                "voltmeter": {
+                    "device": "multimeter",
+                    "delay": 0.1,
+                    "properties": ["V_m"],
+                    "units": ["bla"],
+                    "targetting": {
+                        "strategy": "cell_model",
+                        "cell_models": ["gif_pop_psc_exp"],
+                    },
+                },
+            }
+        )
+        with self.assertRaises(ConfigurationError):
+            Scaffold(cfg, self.storage)
+
+    def test_dc_generator(self):
+        duration = 100
+        resolution = 0.1
+        cfg = Configuration(
+            {
+                "name": "test",
+                "storage": {"engine": "hdf5"},
+                "network": {"x": 1, "y": 1, "z": 1},
+                "partitions": {"B": {"type": "layer", "thickness": 1}},
+                "cell_types": {"A": {"spatial": {"radius": 1, "count": 1}}},
+                "placement": {
+                    "placement_A": {
+                        "strategy": "bsb.placement.strategy.FixedPositions",
+                        "cell_types": ["A"],
+                        "partitions": ["B"],
+                        "positions": [[1, 1, 1]],
+                    }
+                },
+                "connectivity": {},
+                "after_connectivity": {},
+                "simulations": {
+                    "test": {
+                        "simulator": "nest",
+                        "duration": duration,
+                        "resolution": resolution,
+                        "cell_models": {
+                            "A": {
+                                "model": "iaf_cond_alpha",
+                                "constants": {
+                                    "V_reset": -70,  # V_m, E_L and V_reset are the same
+                                },
+                            }
+                        },
+                        "connection_models": {},
+                        "devices": {
+                            "dc_generator": {
+                                "device": "dc_generator",
+                                "delay": resolution,
+                                "weight": 1.0,
+                                "amplitude": 200,  # Low enough so the neuron does not spike
+                                "start": 50,
+                                "stop": 60,
+                                "targetting": {
+                                    "strategy": "cell_model",
+                                    "cell_models": ["A"],
+                                },
+                            },
+                            "voltmeter_A": {
+                                "device": "multimeter",
+                                "delay": resolution,
+                                "properties": ["V_m"],
+                                "units": ["mV"],
+                                "targetting": {
+                                    "strategy": "cell_model",
+                                    "cell_models": ["A"],
+                                },
+                            },
+                        },
+                    }
+                },
+            }
+        )
+
+        netw = Scaffold(cfg, self.storage)
+        netw.compile()
+        results = netw.run_simulation("test")
+        v_ms = np.array(results.analogsignals[0])[:, 0]
+        self.assertAll(v_ms[: int(50 / resolution) + 1] == -70)
+        self.assertAll(
+            v_ms[int(50 / resolution) + 1 : int(60 / resolution) + 1] > -70,
+            "Current injected should raise membrane potential",
+        )
+
     def test_nest_randomness(self):
         nest.ResetKernel()
         nest.resolution = 0.1
@@ -227,7 +361,6 @@ def test_nest_randomness(self):
         nest.Connect(A, spikeA)
         nest.Simulate(1000.0)
         spike_times_nest = spikeA.get("events")["times"]
-        print(spike_times_nest)
 
         conf = {
             "name": "test",
@@ -292,9 +425,46 @@ def test_nest_randomness(self):
                 "std": 20.0,
             },
         )
-        # Test with an unknown distribution
-        conf["simulations"]["test"]["cell_models"]["A"]["constants"]["V_m"][
-            "distribution"
-        ] = "bean"
+
+    def test_unknown_distribution(self):
+        conf = {
+            "name": "test",
+            "storage": {"engine": "hdf5"},
+            "network": {"x": 1, "y": 1, "z": 1},
+            "partitions": {"B": {"type": "layer", "thickness": 1}},
+            "cell_types": {"A": {"spatial": {"radius": 1, "count": 1}}},
+            "placement": {
+                "placement_A": {
+                    "strategy": "bsb.placement.strategy.FixedPositions",
+                    "cell_types": ["A"],
+                    "partitions": ["B"],
+                    "positions": [[1, 1, 1]],
+                }
+            },
+            "connectivity": {},
+            "after_connectivity": {},
+            "simulations": {
+                "test": {
+                    "simulator": "nest",
+                    "duration": 1000,
+                    "resolution": 0.1,
+                    "cell_models": {
+                        "A": {
+                            "model": "gif_cond_exp",
+                            "constants": {
+                                "I_e": 200.0,
+                                "V_m": {
+                                    "distribution": "bean",
+                                    "mean": -70,
+                                    "std": 20.0,
+                                },
+                            },
+                        }
+                    },
+                    "connection_models": {},
+                    "devices": {},
+                }
+            },
+        }
         with self.assertRaises(CastError):
             Configuration(conf)