""" This module handles virtually all the necessary processes for setting up stimulation and recording from electrodes. One more basic function (analysis.set_recording) that is used in this module is defined in the module analysis. This is so because this function is not necessarily used by an electrode """ import os import numpy as np from collections import OrderedDict from neuron import h import anatomy import geometry as geo import analysis as als import workspace as ws class CuffRing(): """ Ring of a cuff electrode Input parameters: z: Position of the ring along the z-axis d: Diameter npads: Number of pads thts: Angular positions of the pads Note about the stimulation set-up: this class has no method set_stimulation because some rings won't have any stimulation at all if it is not defined in the configuration files """ def __init__(self, ring_number, z, center, d, npads, thts): self.ring_number = ring_number self.z = z self.center = center self.d = d self.npads = npads self.thts = thts # Dictionary with the information about all the injections self.injections = {} # Dictionary with the information about all the recordings self.recordings = {} # Locate the contact points with the nerve self.locate_pads_on_nerve() def locate_pads_on_nerve(self): """ Find the locations in the nerve corresponding to the pads """ # Iterate over the ring's pads self.contact_points = [] for i_pad in range(self.npads): # The pad's angular position determines which cable cell(s) # is (are) stimulated # Select which of the targets is in contact # with a pad 'i_pad' according to its angular position angdiff = np.abs(ws.contour_angles - self.thts[i_pad]) which_cables = np.where(angdiff == angdiff.min())[0] # Iterate over possible cables ncables_here = len(which_cables) for i_cable in which_cables: # For cable_cell i, see where this is try: i_sec, zpos = anatomy.locate(ws.seclens[i_cable], self.z) except TypeError: msg = 'ERROR: A ring of an electrode could not be placed on the nerve' ws.log(msg) ws.terminate() # Add the point point = { 'cable': i_cable, 'section': i_sec, 'z': zpos } try: self.contact_points[i_pad].append(point) except IndexError: self.contact_points.append([point]) def assign_injections(self, pad, protocol): """ Assign current injections to the pads """ # Dictionary with the information about all the injections self.injections[pad] = [] # Iterate over the contact points of this pad # NO! Just one pad is allowed. Otherwise, the current with # intensity 'amp', which is meant for only one injection point, # would be added to all the self.contact_points, so we would # have (len(self.contact_points) - 1) more current injections # than we intended point = self.contact_points[pad][0] # IClamp object(s) # Only square pulses are supported for now if protocol['type'] == 'square pulses': square_pulses(self.injections[pad], point, protocol) if protocol['type'] == 'ground': ground_segment(point) def set_recordings(self): """ Assign the recording vectors to each pad """ # Iterate over all pads on this ring for pad in range(self.npads): # Dictionary with the information about all the recordings self.recordings[pad] = {} # Only one is allowed for each pad in this case point = self.contact_points[pad][0] # Set the recording for this pad set_recording(self.recordings[pad], point, 'vext[0]') class CuffElectrode(): """ Cuff electrode with as many rings as one wants Input parameters: z: Position (center) along the z-axis nrings: Number of rings rng_sp: Separation between rings nppr: Number of pads per ring l: Length of the cuff d_in: Inner diameter thickness: Insulator thickness rho: Insulator resistivity (Ohm*cm) """ def __init__(self, name): self.name = name settings = ws.electrodes_settings[name] self.type = settings['type'] self.role = settings['role'] self.z = settings['z-position'] self.nrings = settings['number of rings'] self.rng_sp = settings['inter-ring distance'] self.nppr = settings['pads per ring'] self.length = settings['length'] self.thickness = settings['thickness'] self.rho = settings['insulator resistivity'] self.center = ws.nvt.circumcenter self.d_in = ws.nvt.circumdiameter self.d_out = self.d_in + 2. * self.thickness if self.role == 'stimulation': self.stimulation_protocol = \ settings['stimulation protocol'] # Necessary by-products # Ends of the cuff z = self.z l = self.length nrings = self.nrings left_end = z - 0.5 * l rght_end = z + 0.5 * l # Positions of the rings rings_halfspan = 0.5 * (nrings - 1) * self.rng_sp z_ring_left = z - rings_halfspan z_ring_rght = z + rings_halfspan self.rings_poss = np.linspace(z_ring_left, z_ring_rght, nrings) # Angular positions of the pads on the rings self.thts = settings['angular positions'] # Consistency check if len(self.thts) != self.nppr: msg = 'ERROR: Different number of pads and angles for them' ws.log(msg) ws.terminate() # Warning if size mismatch if z_ring_left < left_end or z_ring_rght > rght_end: msg = 'ERROR: Cuff electrode rings outside the cuff.' ws.log(msg) ws.terminate() # Store variables into attributes self.left_end = left_end self.rght_end = rght_end self.z_ring_left = z_ring_left self.z_ring_rght = z_ring_rght # Create rings self.create_rings() def create_rings(self): """ Create the rings """ self.rings = [] for i in range(self.nrings): z = self.rings_poss[i] self.rings.append(CuffRing(i, z, self.center, self.d_in, self.nppr, self.thts)) def set_stimulation(self): """ Set up the current injections according to the assigned protocol """ sp = self.stimulation_protocol for pad, protocol in sp.items(): ring_number, pad_no = [int(x) for x in pad.split(',')] self.rings[ring_number].assign_injections(pad_no, protocol) def set_recordings(self): """ Set up the recording vectors on the pads """ for ring in self.rings: ring.set_recordings() class PointElectrodeSet(): """ Point electrode set for intracellular current injections, extracellular current sources, recordings or connections to ground """ def __init__(self, name, index): self.name = name self.index = index self.settings = ws.electrodes_settings[name] self.type = self.settings['type'] self.role = self.settings['role'] def set_stim_info(self): """ Set up just a basic thing about stimulation info """ # Stimulation protocol if self.role == 'stimulation': try: self.stimulation_protocol = \ self.settings['stimulation protocol'][self.index] except KeyError: self.stimulation_protocol = \ self.settings['stimulation protocol'] if self.role == 'ground': self.stimulation_protocol = { 'type': 'ground' } def set_stimulation(self): """ Set up the current injections according to the assigned protocol """ # Set up the injection if we are using the RN for all types of electrodes, or if we are using a "simple" model, only for intracellular injections if (ws.settings["nerve model"] == "resistor network") | (ws.settings["nerve model"] == "simple" and self.type == "intracellular"): # Dictionary containing the information about all # the current injections self.injections = {} # Injections list # There's only one 'pad' for an intracellular electrode self.injections[0] = [] # Fetch protocol protocol = self.stimulation_protocol # Iterate over all the points in this set for point in self.points: # Create IClamp object(s) # Only square pulses are supported for now if protocol['type'] == 'square pulses': ws.log('Setting up injection on %s'%str(point)) square_pulses(self.injections[0], point, protocol) # Connections to ground if protocol['type'] == 'ground': ws.log('Setting up ground on %s'%str(point)) ground_segment(point) def set_recordings(self): """ Set up the recording vectors in the points of this electrode set """ self.recordings = {} self.recordings[0] = {} for point in self.points: set_recording(self.recordings[0], point, 'v') class IntracellularElectrodeSet(PointElectrodeSet): """ Intracellular electrode for intracellular current injections or recordings """ def __init__(self, name, index=0): PointElectrodeSet.__init__(self, name, 0) settings = self.settings try: self.z = settings['z-position'] except KeyError: self.axons = settings['axons'] try: self.node = settings['node of Ranvier'] except KeyError: self.node = settings['internode'] self.position = settings['position'] else: # For a node of Ranvier, inject current in the middle self.position = 0.5 else: # Warning if position mismatch if self.z < 0. or self.z > ws.length: ws.log('ERROR: Electrode is outside the limits of the nerve.') ws.terminate() # Note: the warnings or errors for targetting a non-existing # node will be automatically raised by NEURON # Set points self.points = [] # Setting up these points is necessary to set up stimulation and recording # Besides, it gives the process uniformity across different types of electrodes if self.axons[0] == "all": self.axons = range(ws.naxons_total) for axon in self.axons: axon = int(axon) if ws.settings["nerve model"] == "resistor network": cable_i = ws.nNAELC + axon elif ws.settings["nerve model"] == "simple": cable_i = axon k = 0 if "node of Ranvier" in settings.keys(): for j, sec in enumerate(ws.sections[cable_i]): secname = sec.name() if "NODE" in secname or "node" in secname: if k == self.node: break k += 1 self.points.append( { 'cable': cable_i, 'section': j, 'z': self.position } ) self.set_stim_info() class ExtracellularPointElectrode(PointElectrodeSet): """ Extracellular point electrode for current sources, recordings or connections to ground """ def __init__(self, name, index=0): PointElectrodeSet.__init__(self, name, index) settings = self.settings try: self.xyz = settings['point'] except KeyError: self.xyz = settings['points'][self.index] self.x, self.y, self.z = self.xyz # Find the cable this corresponds to dists = geo.dist((self.x, self.y), (ws.nvt.x, ws.nvt.y)) these = np.where(dists == dists.min())[0] # print('Points:', these) # Choose only one cable cable_i = these[0] # print('Cable:', ws.nvt.cables[cable_i]) # Locate the section and segment (z value) of the cable j_sec, z_on_j = anatomy.locate(ws.seclens[cable_i], self.z) # print('Section and segment:', j_sec, z_on_j) # Cable and point of the cable to work on self.points = [ { 'cable': cable_i, 'section': j_sec, 'z': z_on_j } ] self.set_stim_info() class ExtracellularPointElectrodeSet(): """ Set of ExtracellularPointElectrode instances """ def __init__(self, name): self.name = name self.role = ws.electrodes_settings[name]["role"] self.type = "extracellular points set" self.electrodes = [] def set_stimulation(self): """ Create the list of extracellular point electrodes """ for i, (x, y, z) in ws.electrodes_settings[self.name]["points"].items(): # Create the ExtracellularPointElectrode instance and store # it pe = ExtracellularPointElectrode(self.name, i) pe.set_stimulation() self.electrodes.append(pe) class GroundPointSet(): """ Set of points which are connected to ground """ def __init__(self, name): self.name = name self.role = "ground" self.type = "ground points set" def set_stimulation(self): """ Connect grounds """ if ws.settings["nerve model"] == "resistor network": # Go straigth to connect the wanted points to ground for (x, y, z) in ws.electrodes_settings[self.name]["points"]: # Find the cable this corresponds to dists = geo.dist((x, y), (ws.nvt.x, ws.nvt.y)) these = np.where(dists == dists.min())[0] # print('Points:', these) # Choose only one cable cable_i = these[0] # print('Cable:', ws.nvt.cables[cable_i]) # Locate the section and segment (z value) of the cable j_sec, z_on_j = anatomy.locate(ws.seclens[cable_i], z) # Cable and point of the cable to work on point = { 'cable': cable_i, 'section': j_sec, 'z': z_on_j } # Ground it ground_segment(point) ######################################################################## # Once the classes are defined, define the types of electrodes electrode_types = { "cuff": CuffElectrode, "intracellular": IntracellularElectrodeSet, "extracellular point": ExtracellularPointElectrode, "extracellular points set": ExtracellularPointElectrodeSet, "ground points set": GroundPointSet } ######################################################################## # Functions def create_electrodes(): """ Create the electrodes which are specified in the electrodes.json file """ ws.electrodes = {} for name, el in ws.electrodes_settings.items(): print('Electrode: %s'%name) ws.electrodes[name] = electrode_types[el['type']](name) def set_stimulation(): """ Set up the stimulation protocols described in the json files. This information is already in the electrodes attributes """ for electrode in ws.electrodes.values(): if electrode.role == "stimulation": electrode.set_stimulation() def set_recordings(): """ Set up the stimulation protocols described in the json files. This information is already in the electrodes attributes """ # # Time # set_time_recording() # Electrode recordings for electrode in ws.electrodes.values(): if electrode.role == 'recording': # Tell ws that there are electrode recordings ws.settings["electrode recordings"] = True # Set the recordings electrode.set_recordings() def check_cuff_cover(zi): """ Check if a value over the z-axis is under the region of a cuff electrode """ for name, electrode in ws.electrodes.items(): if electrode.type == 'cuff': if electrode.left_end <= zi <= electrode.rght_end: return name return None def prepare_ground_paths(): """ Calculate the variables that need to be used to connect the outer axons to the distant ground, through the cuff insulators and the container """ # Shortennings nvt = ws.nvt container_rho = ws.container_settings['resistivity'] container_diam = ws.container_settings['diameter'] # Convert the resistivity to Ohm*um container_rho *= 1.e4 # Bases of the cylinder: conductivity to ground # Conductances p.u.area in both cases (mho/cm2) # Just direct grounding ws.xg_distal = 1.e9 ws.xg_proxim = 1.e9 # Actually, I want to think that there's medium at the bases # This base will have the cross-sectional area of the container, # not the nerve, so I need to use that to compute rg # ASSUMPTION *1 # Let's assume its thickness is the same as in the radial direction ws.xg_distal = [] ws.xg_proxim = [] # Iterate over cables for i in range(ws.counters['cables']): area = nvt.free_areas[i] # Area that corresponds to this cable for the base of the # container # NOTE: This is just an approximation that would be ideal # should the nerve be circular in its cross-section # If a nerve is not, it's still OK, because I use its # circumcircle area_augmented = area * (container_diam / nvt.circumdiameter) ** 2 # Conductance in mho # Note that the 'depth' (distance from segment to ground along # the z-axis, which means through the saline bath) is the # container's radius, according to our assumption above (ASSUMPTION *1) ws.xg_distal.append(area_augmented / (container_rho * 0.5 * container_diam)) ws.xg_proxim.append(ws.xg_distal[i]) # Thickness of the containing medium. There's one for each electrode ws.container_settings['thickness'] = {} # Round wall of the cylinder # Area that the current from one segment to ground has to cross (um2) # Length of a segment times its portion of the perimeter area = (ws.length / ws.NAELC_nseg) \ * np.pi * 0.5 * (nvt.circumdiameter + container_diam) / ws.npc # Resistance to ground p.u.area (Ohm*um2) # The distances must be in cm (they have to be converted from um) # No, it now must be in um again, becuase I've put container_rho in Ohm*um rg = {} rg['wout_cuff'] = container_rho * 0.5 * (container_diam - nvt.circumdiameter) rg['with_cuff'] = {} # Conductance to ground in mho (p.u.area (mho/um2) * area (um2)) ws.xg_rwall = {} # Without cuff ws.xg_rwall['wout_cuff'] = area / rg['wout_cuff'] # With cuff ws.xg_rwall['with_cuff'] = {} # Iterate over cuff electrodes for name, el in ws.electrodes.items(): if el.type == 'cuff': # Add the thickness from the electrodes thickness ws.container_settings['thickness'][name] = \ 0.5 * container_diam - (nvt.circumradius \ + ws.electrodes_settings[name]['thickness']) # Round wall of the cylinder # Resistance to ground p.u.area (Ohm*um2) # The resistivity has to be converted from Ohm*cm to Ohm*um rg['with_cuff'][name] = \ 1.e4 * ws.electrodes_settings[name]['insulator resistivity'] \ * ws.electrodes_settings[name]['thickness'] \ + container_rho * ws.container_settings['thickness'][name] # Conductance to ground in mho (p.u.area (mho/um2) * area (um2)) # Conductance ws.xg_rwall['with_cuff'][name] = area / rg['with_cuff'][name] def square_pulses(injections_, point, protocol): """ Assign current injections to the pads """ # Injection segment in hoc seg = ws.sections[point['cable']][point['section']](point['z']) # IClamp object(s) for i, amp in enumerate(protocol['currents']): # Gather the data dur = protocol['pulse durations'][i] delay = protocol['pulse onset times'][i] # Create IClamp object and add its properties injection = h.IClamp(seg) # Intensity from uA to nA injection.amp = amp * 1.e3 injection.dur = dur injection.delay = delay # Add this information as an attribute injections_.append({ 'injection id': ws.counters['current injections'], 'point': point, 'amp': amp, 'dur': dur, 'delay': delay }) # Update the counter ws.counters['current injections'] += 1 # Add the duration and delay to ws.totaldur ws.totaldur = max(ws.totaldur, delay + dur) # Append it to ws ws.injections.append(injection) ws.injections_info = ws.injections_info + injections_ ws.log('Injection added to: %s, %s, %s, %s'%(str(point), str(amp), str(dur), str(delay))) def ground_segment(point): """ This connects a segment directly to ground. It is obviously not a current injection protocol, but still is a necessary kind of protocol in case it's needed, I think. And this is a good place to define it. Because the paths or connections to ground have already been defined in this module previously and declared in h in biophysics.electrics, this replaces whatever path to ground this segment had, if any. If it didn't have it, it creates it. """ # Corresponding segment in hoc seg = ws.sections[point['cable']][point['section']](point['z']) # Level or layer where to point to level = ws.cables[point['cable']].properties['extlayers'] - 1 # Ground it directly seg.xg[level] = 1.e9 def set_recording(recordings_, point, var): """ Assign the recording vectors to each pad """ # This is the contact segment in hoc seg = ws.sections[point['cable']][point['section']](point['z']) # Set up the recording information recordings_.update({ 'recording id': ws.counters['recordings'], 'point': point, 'variable': var }) # Set up the recording object als.set_recording(seg, var)