Changeset 366

Show
Ignore:
Timestamp:
06/16/08 17:32:56 (5 months ago)
Author:
apdavison
Message:

Started reimplementing the neuron2 module from scratch

Files:

Legend:

Unmodified
Added
Removed
Modified
Copied
Moved

  • trunk/src/common.py

    r359 r366  
    5353    # Note that there is a NotImplementedError built-in exception we could use 
    5454    raise Exception("Unimplemented abstract method: %s" % _function_id(obj, 1)) 
     55 
     56def is_listlike(obj): 
     57    return hasattr(obj, "__len__") and not isinstance(obj, basestring) 
    5558 
    5659def build_translations(*translation_list): 
    …… 
    9396     
    9497    non_parameter_attributes = ('parent', '_cellclass', 'cellclass', 
    95                                 '_position', 'position', 'hocname') 
     98                                '_position', 'position', 'hocname', '_cell') 
    9699     
    97100    def __init__(self): 

  • trunk/src/neuron/__init__.py

    r365 r366  
    242242    return syn_objref 
    243243 
    244 def checkParams(param, val=None): 
    245     """Check parameters are of valid types, normalise the different ways of 
    246        specifying parameters and values by putting everything in a dict. 
    247        Called by set() and Population.set().""" 
    248     if isinstance(param, str): 
    249         if isinstance(val, float) or isinstance(val, int): 
    250             param_dict = {param:float(val)} 
    251         elif isinstance(val,(str, list)): 
    252             param_dict = {param:val} 
    253         else: 
    254             raise common.InvalidParameterValueError 
    255     elif isinstance(param, dict): 
    256         param_dict = param 
    257     else: 
    258         raise common.InvalidParameterValueError 
    259     return param_dict 
     244#def checkParams(param, val=None): 
     245#    """Check parameters are of valid types, normalise the different ways of 
     246#       specifying parameters and values by putting everything in a dict. 
     247#       Called by set() and Population.set().""" 
     248#    if isinstance(param, str): 
     249#        if isinstance(val, float) or isinstance(val, int): 
     250#            param_dict = {param:float(val)} 
     251#        elif isinstance(val,(str, list)): 
     252#            param_dict = {param:val} 
     253#        else: 
     254#            raise common.InvalidParameterValueError 
     255#    elif isinstance(param, dict): 
     256#        param_dict = param 
     257#    else: 
     258#        raise common.InvalidParameterValueError 
     259#    return param_dict 
    260260 
    261261class HocToPy: 

  • trunk/src/neuron2/__init__.py

    r342 r366  
    22""" 
    33nrnpython implementation of the PyNN API. 
    4  
    5 This is an attempt at a parallel-enabled implementation.                                                                 
     4                                                       
    65$Id:__init__.py 188 2008-01-29 10:03:59Z apdavison $ 
    76""" 
    87__version__ = "$Rev: 191 $" 
    98 
    10 import neuron 
    11 from pyNN import __path__ as pyNN_path 
    12 neuron.h.nrn_load_dll("%s/hoc/i686/.libs/libnrnmech.so" % pyNN_path[0]) # put this in setup()? 
    13  
    149from pyNN.random import * 
    15 from math import * 
    16 from pyNN import common 
     10from pyNN.neuron2.utility import * 
     11from pyNN import common, utility 
    1712from pyNN.neuron2.cells import * 
    1813from pyNN.neuron2.connectors import * 
    1914from pyNN.neuron2.synapses import * 
    2015 
    21 import os.path 
    22 import types 
     16from math import * 
    2317import sys 
    2418import numpy 
    2519import logging 
    26  
    27 gid           = 0 
    28 ncid          = 0 
    29 gidlist       = [] 
    30 vfilelist     = {} 
    31 spikefilelist = {} 
    32 dt            = 0.1 
    33 running       = False 
    34 initialised   = False 
    35  
    36 # ============================================================================== 
    37 #   Utility classes and functions 
    38 # ============================================================================== 
    39  
    40 class ID(common.ID): 
    41     """ 
    42     Instead of storing ids as integers, we store them as ID objects, 
    43     which allows a syntax like: 
    44         p[3,4].tau_m = 20.0 
    45     where p is a Population object. The question is, how big a memory/performance 
    46     hit is it to replace integers with ID objects? 
    47     """ 
    48      
    49     def __init__(self, n): 
    50         common.ID.__init__(self, n) 
    51         self.hocname = None 
    52         self.model_obj = None 
    53      
    54     def __getattr__(self, name): 
    55         """Note that this currently does not translate units.""" 
    56         if type(self.cellclass) == type and issubclass(self.cellclass, common.StandardCellType): 
    57             translated_name = self.cellclass.translations[name][0] 
    58         else: 
    59             translated_name = name 
    60         if self.hocname: 
    61             return getattr(neuron.h, '%s.%s' % (self.hocname, translated_name)) 
    62         else: 
    63             return getattr(neuron.h, 'cell%d.%s' % (int(self), translated_name)) 
    64      
    65     def setParameters(self,**parameters): 
    66         # We perform a call to the low-level function set() of the API. 
    67         # If the cellclass is not defined in the ID object, we have an error: 
    68         if (self.cellclass == None): 
    69             raise Exception("Unknown cellclass") 
    70         else: 
    71             #Otherwise we use the ID one. Nevertheless, here we have a small problem in the 
    72             #parallel framework. Suppose a population is created, distributed among 
    73             #several nodes. Then a call like cell[i, j].set() should be performed only on the 
    74             #node who owns the cell. To do that, if the node doesn't have the cell, a call to set() 
    75             #does nothing... 
    76             set(self, self.cellclass, parameters) 
    77  
    78     def getParameters(self): 
    79         params = {} 
    80         for k in self.cellclass.translations.keys(): 
    81             params[k] = self.__getattr__(k) 
    82         return params 
    83  
    84 def list_standard_models(): 
    85     return [obj for obj in globals().values() if isinstance(obj, type) and issubclass(obj, common.StandardCellType)] 
    86  
    87 # ============================================================================== 
    88 #   Module-specific functions and classes (not part of the common API) 
    89 # ============================================================================== 
    90  
    91 def _translate_synapse_type(synapse_type, weight=None): 
    92     """ 
    93     If synapse_type is given (not None), it is used to determine whether the 
    94     synapse is excitatory or inhibitory. 
    95     Otherwise, the synapse type is inferred from the sign of the weight. 
    96     Much testing needed to check if this behaviour matches nest and pcsim. 
    97     """ 
    98      
    99     if synapse_type: 
    100         if synapse_type == 'excitatory': 
    101             syn_objref = "esyn" 
    102         elif synapse_type == 'inhibitory': 
    103             syn_objref = "isyn" 
    104         else: 
    105             # More sophisticated treatment needed once we have more sophisticated synapse 
    106             # models, e.g. NMDA... 
    107             #raise common.InvalidParameterValueError, synapse_type, "valid types are 'excitatory' or 'inhibitory'" 
    108             syn_objref = synapse_type 
    109     else: 
    110         if weight is None or weight >= 0.0: 
    111             syn_objref = "esyn" 
    112         else: 
    113             syn_objref = "isyn" 
    114     return syn_objref 
    115  
    116 def checkParams(param, val=None): 
    117     """Check parameters are of valid types, normalise the different ways of 
    118        specifying parameters and values by putting everything in a dict. 
    119        Called by set() and Population.set().""" 
    120     if isinstance(param, str): 
    121         if isinstance(val, float) or isinstance(val, int): 
    122             param_dict = {param:float(val)} 
    123         elif isinstance(val,(str, list)): 
    124             param_dict = {param:val} 
    125         else: 
    126             raise common.InvalidParameterValueError 
    127     elif isinstance(param, dict): 
    128         param_dict = param 
    129     else: 
    130         raise common.InvalidParameterValueError 
    131     return param_dict 
    132                  
     20import operator 
     21 
     22# Global variables 
     23gid_counter = 0 
     24initialised = False 
     25running = False 
     26quit_on_end = True 
     27recorder_list = [] 
     28 
    13329# ============================================================================== 
    13430#   Functions for simulation set-up and control 
    13531# ============================================================================== 
    13632 
    137 def setup(timestep=0.1, min_delay=0.1, max_delay=0.1, debug=False,**extra_params): 
     33def setup(timestep=0.1, min_delay=0.1, max_delay=10.0, debug=False,**extra_params): 
    13834    """ 
    13935    Should be called at the very beginning of a script. 
    …… 
    14137    simulator but not by others. 
    14238    """ 
    143     global dt, nhost, myid, _min_delay, logger, initialised, pc 
    144     dt = timestep 
    145     _min_delay = min_delay 
    146      
    147     # Initialisation of the log module. To write in the logfile, simply enter 
    148     # logging.critical(), logging.debug(), logging.info(), logging.warning()  
    149     if debug: 
    150         logging.basicConfig(level=logging.DEBUG, 
    151                     format='%(asctime)s %(levelname)s %(message)s', 
    152                     filename='neuron2.log', 
    153                     filemode='w') 
    154     else: 
    155         logging.basicConfig(level=logging.INFO, 
    156                     format='%(asctime)s %(levelname)s %(message)s', 
    157                     filename='neuron2.log', 
    158                     filemode='w') 
    159          
    160     logging.info("Initialization of NEURON (use setup(.., debug=True) to see a full logfile)") 
    161      
    162     # All the objects that will be used frequently in the hoc code are declared in the setup 
    163     neuron.h.dt = dt 
    164     if initialised: 
    165         pass 
    166     else: 
    167         #hoc_commands = [ 
    168         #    'tmp = xopen("%s")' % os.path.join(pyNN_path[0],'hoc','standardCells.hoc'), 
    169         #    'tmp = xopen("%s")' % os.path.join(pyNN_path[0],'hoc','odict.hoc'), 
    170         #    'objref pc', 
    171         #    'pc = new ParallelContext()', 
    172         #    'dt = %f' % dt, 
    173         #    'create dummy_section', 
    174         #    'access dummy_section', 
    175         #    'objref netconlist, nil', 
    176         #    'netconlist = new List()',  
    177         #    'strdef cmd', 
    178         #    'strdef fmt',  
    179         #    'objref nc',  
    180         #    'objref rng', 
    181         #    'objref cell'] 
    182         neuron.xopen("%s" % os.path.join(pyNN_path[0],'hoc','standardCells.hoc') ) 
    183         neuron.xopen("%s" % os.path.join(pyNN_path[0],'hoc','odict.hoc') ) # to suppress - odict no longer needed? 
     39    global initialised, quit_on_end, running, pc 
     40    if not initialised: 
     41        h('min_delay = 0') 
     42        h('tstop = 0') 
    18443        pc = neuron.ParallelContext() 
    185          
    186         #---Experimental--- Optimize the simulation time ? / Reduce inter-processors exchanges ? 
    187         #hoc_commands += [ 
    188         #    'tmp   = pc.spike_compress(1,0)'] 
    18944        pc.spike_compress(1,0) 
    190         if extra_params.has_key('use_cvode') and extra_params['use_cvode'] == True: 
    191             #hoc_commands += [ 
    192             #    'objref cvode', 
    193             #    'cvode = new CVode()', 
    194             #    'cvode.active(1)'] 
    195             cvode = neuron.CVode() 
    196             cvode.active(1) 
    197          
    198     #hoc_execute(hoc_commands,"--- setup() ---") 
    199     nhost = int(pc.nhost()) 
    200     if nhost < 2: 
    201         nhost = 1; myid = 0 
    202     else: 
    203         myid = int(pc.id()) 
    204     print "\nHost #%d of %d" % (myid+1, nhost) 
    205      
    206     initialised = True 
    207     return myid 
     45        cvode = neuron.CVode() 
     46        utility.init_logging("neuron2.log.%d" % rank(), debug) 
     47        logging.info("Initialization of NEURON (use setup(.., debug=True) to see a full logfile)") 
     48    h.dt = timestep 
     49    h.tstop = 0 
     50    h.min_delay = min_delay 
     51    running = False 
     52    if 'quit_on_end' in extra_params: 
     53        quit_on_end = extra_params['quit_on_end'] 
     54    if extra_params.has_key('use_cvode'): 
     55        cvode.active(int(extra_params['use_cvode'])) 
     56    return rank() 
    20857 
    20958def end(compatible_output=True): 
    21059    """Do any necessary cleaning up before exiting.""" 
    211     global logfile, myid, vfilelist, spikefilelist 
    212     #hoc_commands = [] 
    213     print vfilelist 
    214     if len(vfilelist) > 0: 
    215         #hoc_commands = ['objref fileobj', 
    216         #                'fileobj = new File()'] 
    217         tstop = neuron.h.tstop 
    218         header = "# dt = %g\n# n = %d\n" % (dt, int(tstop/dt)) 
    219         for filename, cell_list in vfilelist.items(): 
    220             fileobj = neuron.File(filename, 'w') 
    221             fileobj.write(header) 
    222             for cell in cell_list: 
    223                 fmt = "%s\t%d\n" % ("%.6g", cell.gid) 
    224                 cell.vtrace.printf(fileobj.hoc_obj, fmt) 
    225             fileobj.close() 
    226     if len(spikefilelist) > 0: 
    227         header = "# dt = %g\n# "% dt 
    228         for filename, cell_list in spikefilelist.items(): 
    229             fileobj = neuron.File(filename, 'w') 
    230             for cell in cell_list: 
    231                 fmt = "%s\t%d\n" % ("%.2f", cell.gid) 
    232                 cell.spiketimes.printf(fileobj.hoc_obj, fmt) 
    233             fileobj.close() 
     60    for recorder in recorder_list: 
     61        recorder.write(gather=False, compatible_output=compatible_output) 
    23462    pc.runworker() 
    23563    pc.done() 
    236     neuron.quit() 
    237     logging.info("Finishing up with NEURON.") 
    238     sys.exit(0) 
    239  
     64    if quit_on_end: 
     65        logging.info("Finishing up with NEURON.") 
     66        h.quit() 
     67         
    24068def run(simtime): 
    24169    """Run the simulation for simtime ms.""" 
    24270    global running 
     71    logging.info("Running the simulation for %d ms" % simtime) 
    24372    if not running: 
    24473        running = True 
    245         #neuron.h.tstop = simtime 
    246         print "dt        = %f" % dt 
    247         print "min delay = ", pc.set_maxstep(100) 
    248         #neuron.finitialize() 
    249         neuron.init() 
    250     print "tstop     = ", simtime 
    251     neuron.h('tstop = %g' % simtime) 
    252     pc.psolve(simtime) 
     74        local_minimum_delay = pc.set_maxstep(10) 
     75        h.finitialize() 
     76        h.tstop = 0 
     77        logging.debug("local_minimum_delay on host #%d = %g" % (rank(), local_minimum_delay)) 
     78        if num_processes() > 1: 
     79            assert local_minimum_delay >= get_min_delay(),\ 
     80                   "There are connections with delays (%g) shorter than the minimum delay (%g)" % (local_minimum_delay, get_min_delay()) 
     81    h.tstop = simtime 
     82    pc.psolve(h.tstop) 
     83    return get_current_time() 
     84 
     85# ============================================================================== 
     86#   Functions returning information about the simulation state 
     87# ============================================================================== 
     88 
     89def get_current_time(): 
     90    """Return the current time in the simulation.""" 
     91    return h.t 
     92 
     93def get_time_step(): 
     94    return h.dt 
     95common.get_time_step = get_time_step 
     96 
     97def get_min_delay(): 
     98    return h.min_delay 
     99common.get_min_delay = get_min_delay 
     100 
     101def num_processes(): 
     102    return int(pc.nhost()) 
     103 
     104def rank(): 
     105    """Return the MPI rank.""" 
     106    return int(pc.id()) 
     107 
     108def list_standard_models(): 
     109    return [obj for obj in globals().values() if isinstance(obj, type) and issubclass(obj, common.StandardCellType)] 
     110 
     111class ID(int, common.IDMixin): 
     112     
     113    def __init__(self, n): 
     114        int.__init__(n) 
     115        common.IDMixin.__init__(self) 
     116     
     117    def _build_cell(self, celltype, parent=None): 
     118        gid = int(self) 
     119        self._cell = celltype.model(**celltype.parameters)  # create the cell object 
     120        pc.set_gid2node(gid, rank())                        # assign the gid to this node 
     121        nc = neuron.NetCon(self._cell.source, None)         # } associate the cell spike source 
     122        pc.cell(gid, nc.hoc_obj)                            # } with the gid (using a temporary NetCon) 
     123        self.parent = parent 
     124     
     125    def get_native_parameters(self): 
     126        D = {} 
     127        for name in ('tau_m', 'cm', 'v_rest', 'v_thresh', 't_refrac', 
     128                     'i_offset', 'v_reset', 'v_init', 'tau_e', 'tau_i'): 
     129            D[name] = getattr(self._cell, name) 
     130        return D 
     131     
     132    def set_native_parameters(self, parameters): 
     133        for name, val in parameters.items(): 
     134            setattr(self._cell, name, val) 
     135             
    253136 
    254137# ============================================================================== 
    255138#   Low-level API for creating, connecting and recording from individual neurons 
    256139# ============================================================================== 
    257  
    258 def create(cellclass, param_dict=None, n=1): 
    259     """ 
    260     Create n cells all of the same type. 
    261     If n > 1, return a list of cell ids/references. 
    262     If n==1, return just the single id. 
    263     """ 
    264     global gid, gidlist, nhost, myid, pc, nc 
    265      
    266     assert n > 0, 'n must be a positive integer' 
    267     #if isinstance(cellclass, type): 
    268     #    celltype = cellclass(param_dict) 
    269     #    hoc_name = celltype.hoc_name 
    270     #    hoc_commands, argstr = _hoc_arglist([celltype.parameters]) 
    271     #elif isinstance(cellclass, str): 
    272     #    hoc_name = cellclass 
    273     #    hoc_commands, argstr = _hoc_arglist([param_dict]) 
    274     #argstr = argstr.strip().strip(',') 
    275   
    276     # round-robin partitioning 
    277     newgidlist = [i+myid for i in range(gid, gid+n, nhost) if i < gid+n-myid] 
    278     cell_list = [] 
    279     for cell_id in newgidlist: 
    280         #hoc_commands += ['tmp = pc.set_gid2node(%d,%d)' % (cell_id, myid), 
    281         #                 'objref cell%d' % cell_id, 
    282         #                 'cell%d = new %s(%s)' % (cell_id, hoc_name, argstr), 
    283         #                 'tmp = cell%d.connect2target(nil, nc)' % cell_id, 
    284         #                 #'nc = new NetCon(cell%d.source, nil)' % cell_id, 
    285         #                 'tmp = pc.cell(%d, nc)' % cell_id] 
    286         cell = cellclass(param_dict)           # create the cell object 
    287         cell.gid = cell_id                    # and assign its gid 
    288         pc.set_gid2node(cell.gid, myid)       # assign this gid to this node 
    289         nc = neuron.NetCon(cell.source, None) # } associate the cell spike source 
    290         pc.cell(cell.gid, nc.hoc_obj)         # } with the gid (using a temporary NetCon) 
    291         cell_list.append(cell) 
    292  
    293     gidlist.extend(newgidlist) 
    294     #cell_list = [ID(i) for i in range(gid, gid+n)] 
    295     #for id in cell_list: 
    296     #    id.cellclass = cellclass 
    297     gid = gid+n 
    298     if n == 1: 
    299         cell_list = cell_list[0] 
    300     return cell_list 
    301  
    302 def connect(source, target, weight=None, delay=None, synapse_type=None, p=1, rng=None): 
    303     """Connect a source of spikes to a synaptic target. source and target can 
    304     both be individual cells or lists of cells, in which case all possible 
    305     connections are made with probability p, using either the random number 
    306     generator supplied, or the default rng otherwise. 
    307     Weights should be in nA or µS.""" 
    308     global ncid, gid, gidlist, _min_delay, pc 
    309     if type(source) != types.ListType: 
    310         source = [source] 
    311     if type(target) != types.ListType: 
    312         target = [target] 
    313     if weight is None:  weight = 0.0 
    314     if delay  is None:  delay = _min_delay 
    315     syn_objref = _translate_synapse_type(synapse_type, weight) 
    316     nc_start = ncid 
    317     hoc_commands = [] 
    318     conn_list = [] 
    319     for tgt in target: 
    320         if tgt.gid > gid or tgt.gid < 0: # or not isinstance(tgt, int): 
    321             raise common.ConnectionError, "Postsynaptic cell id %s does not exist." % str(tgt.gid) 
    322         else: 
    323             if tgt.gid in gidlist: # only create connections to cells that exist on this machine 
    324                 if p < 1: 
    325                     if rng: # use the supplied RNG 
    326                         rarr = self.rng.uniform(0,1, len(source)) 
    327                     else:   # use the default RNG 
    328                         rarr = numpy.random.uniform(0,1, len(source)) 
    329                 for j, src in enumerate(source): 
    330                     if src.gid > gid or src.gid < 0: # or not isinstance(src, int): 
    331                         raise common.ConnectionError, "Presynaptic cell id %s does not exist." % str(src.gid) 
    332                     else: 
    333                         if p >= 1.0 or rarr[j] < p: # might be more efficient to vectorise the latter comparison 
    334                              
    335                             nc = pc.gid_connect(src.gid, 
    336                                                 getattr(tgt, syn_objref).hoc_obj) 
    337                             nc.delay = delay 
    338                             nc.weight[0] = weight 
    339                             conn_list.append(nc) 
    340                             ncid += 1 
    341                             print nc, weight, nc.weight[0] 
    342     #hoc_execute(hoc_commands, "--- connect(%s,%s) ---" % (str(source), str(target))) 
    343     #return range(nc_start, ncid) # why not return a list of NetCon objects, instead of ids? 
    344     if len(conn_list) == 1: 
    345         conn_list = conn_list[0] 
    346     return conn_list 
    347  
    348 def set(cells, cellclass, param, val=None): 
    349     """Set one or more parameters of an individual cell or list of cells. 
    350     param can be a dict, in which case val should not be supplied, or a string 
    351     giving the parameter name, in which case val is the parameter value. 
    352     cellclass must be supplied for doing translation of parameter names.""" 
    353     global gidlist 
    354      
    355     param_dict = checkParams(param, val) 
    356  
    357     if type(cellclass) == type and issubclass(cellclass, common.StandardCellType): 
    358         param_dict = cellclass({}).translate(param_dict) 
    359     if not isinstance(cells, list): 
    360         cells = [cells]     
    361     hoc_commands = [] 
    362     for param, val in param_dict.items(): 
    363         if isinstance(val, str): 
    364             fmt = 'pc.gid2cell(%d).%s = "%s"' 
    365         elif isinstance(val, list): 
    366             cmds, argstr = _hoc_arglist([val]) 
    367             hoc_commands += cmds 
    368             fmt = 'pc.gid2cell(%d).%s = %s' 
    369             val = argstr 
    370         else: 
    371             fmt = 'pc.gid2cell(%d).%s = %g' 
    372         for cell in cells: 
    373             if cell in gidlist: 
    374                 hoc_commands += [fmt % (cell, param, val), 
    375                                  'tmp = pc.gid2cell(%d).param_update()' % cell] 
    376     hoc_execute(hoc_commands, "--- set() ---") 
    377  
    378 def record(source, filename): 
    379     """Record spikes to a file. source can be an individual cell or a list of 
    380     cells.""" 
    381     # would actually like to be able to record to an array and choose later 
    382     # whether to write to a file. 
    383     global spikefilelist, gidlist 
    384     if type(source) != types.ListType: 
    385         source = [source] 
    386     hoc_commands = [] 
    387     if not spikefilelist.has_key(filename): 
    388         spikefilelist[filename] = [] 
    389     for src in source: 
    390         if src in gidlist: 
    391             #hoc_commands += ['tmp = cell%d.record(1)' % src] 
    392             src.record(1) 
    393             spikefilelist[filename] += [src] # writing to file is done in end() 
    394     #hoc_execute(hoc_commands, "---record() ---") 
    395  
    396 def record_v(source, filename): 
    397     """ 
    398     Record membrane potential to a file. source can be an individual cell or 
    399     a list of cells.""" 
    400     # would actually like to be able to record to an array and 
    401     # choose later whether to write to a file. 
    402     global vfilelist, gidlist 
    403     if type(source) != types.ListType: 
    404         source = [source] 
    405     hoc_commands = [] 
    406     if not vfilelist.has_key(filename): 
    407         vfilelist[filename] = [] 
    408     for src in source: 
    409         if src.gid in gidlist: 
    410             #hoc_commands += ['tmp = cell%d.record_v(1,%g)' % (src, dt)] 
    411             src.record_v(1) 
    412             vfilelist[filename] += [src] # writing to file is done in end() 
    413     #hoc_execute(hoc_commands, "---record_v() ---") 
    414140 
    415141# ============================================================================== 
    …… 
    441167        label is an optional name for the population. 
    442168        """ 
    443         global gid, myid, nhost, gidlist, fullgidlist 
     169        global gid_counter 
     170        common.Population.__init__(self, dims, cellclass, cellparams, label) 
     171        self.recorders = {'spikes': Recorder('spikes', population=self), 
     172                          'v': Recorder('v', population=self)} 
     173        self.first_id = gid_counter 
     174        self.last_id = gid_counter + self.size 
     175        self.label = self.label or 'population%d' % Population.nPop 
    444176         
    445         common.Population.__init__(self, dims, cellclass, cellparams, label) 
    446         #if self.ndim > 1: 
    447         #    for i in range(1, self.ndim): 
    448         #        if self.dim[i] != self.dim[0]: 
    449         #            raise common.InvalidDimensionsError, "All dimensions must be the same size (temporary restriction)." 
    450  
    451         # set the steps list, used by the __getitem__() method. 
    452         self.steps = [1]*self.ndim 
    453         for i in xrange(self.ndim-1): 
    454             for j in range(i+1, self.ndim): 
    455                 self.steps[i] *= self.dim[j] 
    456  
    457         if isinstance(cellclass, type): 
    458             self.celltype = cellclass(cellparams) 
    459             self.cellparams = self.celltype.parameters 
    460             hoc_name = self.celltype.hoc_name 
    461         elif isinstance(cellclass, str): # not a standard model 
    462             hoc_name = cellclass 
     177        # Build the arrays of cell ids 
     178        # Cells on the local node are represented as ID objects, other cells by integers 
     179        # All are stored in a single numpy array for easy lookup by address 
     180        # The local cells are also stored in a list, for easy iteration 
     181        self._all_ids = numpy.array([id for id in range(self.first_id, self.last_id)], ID) #.reshape(self.dim) 
     182        # _mask_local is used to extract those elements from arrays that apply to the cells on the current node, e.g. for tset() 
     183        self._mask_local = self._all_ids%num_processes()==0 # round-robin distribution of cells between nodes 
     184        ## _map_local is used to map addresses to the index within the list of local cells 
     185        #self._map_local = numpy.where(self._mask_local, self._all_ids/int(num_processes()), None) 
     186        for i,(id,is_local) in enumerate(zip(self._all_ids, self._mask_local)): 
     187            if is_local: 
     188                self._all_ids[i] = ID(id) 
     189        # _local_ids is a list containing only those cells that exist on the current node 
     190        self._local_ids = self._all_ids[self._mask_local] 
     191        self._all_ids = self._all_ids.reshape(self.dim) 
     192        self._mask_local = self._mask_local.reshape(self.dim) 
    463193         
    464         if self.cellparams is not None: 
    465             hoc_commands, argstr = _hoc_arglist([self.cellparams]) 
    466             argstr = argstr.strip().strip(',') 
    467         else: 
    468             hoc_commands = [] 
    469             argstr = '' 
    470      
    471         if not self.label: 
    472             self.label = 'population%d' % Population.nPop 
    473         self.hoc_label = self.label.replace(" ","_") 
    474          
    475         self.record_from = { 'spiketimes': [], 'vtrace': [] } 
    476          
    477          
    478         # Now the gid and cellclass are stored as instance of the ID class, which will allow a syntax like 
    479         # p[i, j].set(param, val). But we have also to deal with positions : a population needs to know ALL the positions 
    480         # of its cells, and not only those of the cells located on a particular node (i.e in self.gidlist). So 
    481         # each population should store what we call a "fullgidlist" with the ID of all the cells in the populations  
    482         # (and therefore their positions) 
    483         self.fullgidlist = numpy.array([ID(i) for i in range(gid, gid+self.size) if i < gid+self.size], ID) 
    484         self.cell = self.fullgidlist 
    485          
    486         # self.gidlist is now derived from self.fullgidlist since it contains only the cells of the population located on 
    487         # the node 
    488         self.gidlist     = [self.fullgidlist[i+myid] for i in range(0, len(self.fullgidlist), nhost) if i < len(self.fullgidlist)-myid] 
    489         self.gid_start   = gid 
    490  
    491         # Write hoc commands 
    492         hoc_commands += ['objref %s' % self.hoc_label, 
    493                          '%s = new List()' % self.hoc_label] 
    494  
    495         for cell_id in self.gidlist: 
    496             hoc_commands += ['tmp = pc.set_gid2node(%d,%d)' % (cell_id, myid), 
    497                              'cell = new %s(%s)' % (hoc_name, argstr), 
    498                              #'nc = new NetCon(cell.source, nil)', 
    499                              'tmp = cell.connect2target(nil, nc)', 
    500                              'tmp = pc.cell(%d, nc)' % cell_id, 
    501                              'tmp = %s.append(cell)' %(self.hoc_label)]        
    502         hoc_execute(hoc_commands, "--- Population[%s].__init__() ---" %self.label) 
     194        self.celltype = cellclass(cellparams) 
     195        for id in self._local_ids: 
     196            id._build_cell(self.celltype, parent=self) 
     197        gid_counter += self.size 
    503198        Population.nPop += 1 
    504         gid = gid+self.size 
    505  
    506         # We add the gidlist of the population to the global gidlist 
    507         gidlist += self.gidlist 
    508          
    509         # By default, the positions of the cells are their coordinates, given by the locate() 
    510         # method. Note that each node needs to know all the positions of all the cells  
    511         # in the population 
    512         for cell_id in self.fullgidlist: 
    513             cell_id.parent = self 
    514             #cell_id.setPosition(self.locate(cell_id)) 
    515                      
    516         # On the opposite, each node has to know only the precise hocname of its cells, if we 
    517         # want to be able to use the low level set() function 
    518         for cell_id in self.gidlist: 
    519             cell_id.hocname = "%s.o(%d)" % (self.hoc_label, self.gidlist.index(cell_id)) 
    520  
     199        logging.info(self.describe('Creating Population "%(label)s" of shape %(dim)s, '+ 
     200                                   'containing `%(celltype)s`s with indices between %(first_id)s and %(last_id)s')) 
     201        logging.debug(self.describe()) 
     202     
    521203    def __getitem__(self, addr): 
    522204        """Return a representation of the cell with coordinates given by addr, 
    …… 
    526208             p[2,3] is equivalent to p.__getitem__((2,3)). 
    527209        """ 
    528  
    529         global gidlist 
    530  
    531         # What we actually pass around are gids. 
    532210        if isinstance(addr, int): 
    533211            addr = (addr,) 
    534         if len(addr) != len(self.dim): 
     212        if len(addr) == self.ndim: 
     213            id = self._all_ids[addr] 
     214        else: 
    535215            raise common.InvalidDimensionsError, "Population has %d dimensions. Address was %s" % (self.ndim, str(addr)) 
    536         index = 0 
    537         for i, s in zip(addr, self.steps): 
    538             index += i*s 
    539         id = index + self.gid_start 
    540         assert addr == self.locate(id), 'index=%s addr=%s id=%s locate(id)=%s' % (index, addr, id, self.locate(id)) 
    541         # We return the gid as an ID object. Note that each instance of Populations 
    542         # distributed on several node can give the ID object, because fullgidlist is duplicated 
    543         # and common to all the node (not the case of global gidlist, or self.gidlist) 
    544         return self.fullgidlist[index] 
    545  
     216        if addr != self.locate(id): 
     217            raise IndexError, 'Invalid cell address %s' % str(addr) 
     218        return id 
     219     
    546220    def __iter__(self): 
    547         return self.__gid_gen() 
     221        """Iterator over cell ids.""" 
     222        return iter(self._local_ids) 
    548223 
    549224    def __address_gen(self): 
    …… 
    552227        returning addresses. 
    553228        """ 
    554         for i in self.gidlist: 
     229        for i in self.__iter__(): 
    555230            yield self.locate(i) 
    556      
    557     def __gid_gen(self): 
    558         """ 
    559         Generator to produce an iterator over all cells on this node, 
    560         returning gids. 
    561         """ 
    562         for i in self.gidlist: 
    563             yield i 
    564          
     231 
    565232    def addresses(self): 
     233        """Iterator over cell addresses.""" 
    566234        return self.__address_gen() 
    567      
     235 
    568236    def ids(self): 
    569         return self.__gid_gen() 
    570          
     237        """Iterator over cell ids.""" 
     238        return self.__iter__() 
     239     
    571240    def locate(self, id): 
    572241        """Given an element id in a Population, return the coordinates. 
    …… 
    574243                         7 9 
    575244        """ 
    576         # id should be a gid 
    577         assert isinstance(id, int), "id is %s, not int" % type(id) 
    578         id -= self.gid_start 
     245        id = id - self.first_id 
    579246        if self.ndim == 3: 
    580247            rows = self.dim[1]; cols = self.dim[2] 
    581248            i = id/(rows*cols); remainder = id%(rows*cols) 
    582249            j = remainder/cols; k = remainder%cols 
    583             coords = (i, j, k) 
     250            coords = (i,j,k) 
    584251        elif self.ndim == 2: 
    585252            cols = self.dim[1] 
    586253            i = id/cols; j = id%cols 
    587             coords = (i, j) 
     254            coords = (i,j) 
    588255        elif self.ndim == 1: 
    589256            coords = (id,) 
    …… 
    591258            raise common.InvalidDimensionsError 
    592259        return coords 
     260 
     261    def get(self, parameter_name, as_array=False): 
     262        """ 
     263        Get the values of a parameter for every cell in the population. 
     264        """ 
     265        values = [getattr(cell, parameter_name) for cell in self] 
     266        if as_array: 
     267            values = numpy.array(values) 
     268        return values 
    593269 
    594270    def set(self, param, val=None): 
    …… 
    601277             p.set({'tau_m':20,'v_rest':-65}) 
    602278        """ 
    603         param_dict = checkParams(param, val) 
    604         if isinstance(self.celltype, common.StandardCellType): 
    605             param_dict = self.celltype.translate(param_dict) 
    606  
    607         strfmt  = '%s.object(tmp).%s = "%s"' % (self.hoc_label,"%s","%s") 
    608         numfmt  = '%s.object(tmp).%s = %s' % (self.hoc_label,"%s","%g") 
    609         listfmt = '%s.object(tmp).%s = %s' % (self.hoc_label,"%s","%s") 
    610         for param, val in param_dict.items(): 
    611             if isinstance(val, str): 
    612                 fmt = strfmt 
    613             elif isinstance(val, list): 
    614                 cmds, argstr = _hoc_arglist([val]) 
    615                 hoc_commands += cmds 
    616                 fmt = listfmt 
    617                 val = argstr 
     279        if isinstance(param, str): 
     280            if isinstance(val, (str, float, int)): 
     281                param_dict = {param: float(val)} 
     282            elif isinstance(val, (list, numpy.ndarray)): 
     283                param_dict = {param: val} 
    618284            else: 
    619                 fmt = numfmt 
    620             # We do the loop in hoc, to speed up the code 
    621             loop = "for tmp = 0, %d" %(len(self.gidlist)-1) 
    622             cmd  = fmt % (param, val) 
    623             hoc_commands = ['cmd="%s { %s success = %s.object(tmp).param_update()}"' %(loop, cmd, self.hoc_label), 
    624                             'success = execute1(cmd)'] 
    625         hoc_execute(hoc_commands, "--- Population[%s].__set()__ ---" %self.label) 
    626  
    627     def tset(self, parametername, valueArray): 
     285                raise common.InvalidParameterValueError 
     286        elif isinstance(param, dict): 
     287            param_dict = param 
     288        else: 
     289            raise common.InvalidParameterValueError 
     290        logging.info("%s.set(%s)", self.label, param_dict) 
     291        for cell in self: 
     292            cell.set_parameters(**param_dict) 
     293     
     294    def tset(self, parametername, value_array): 
    628295        """ 
    629296        'Topographic' set. Set the value of parametername to the values in 
    630297        value_array, which must have the same dimensions as the Population. 
    631298        """ 
    632         if self.dim == valueArray.shape: 
    633             values = numpy.reshape(valueArray, valueArray.size) 
    634             values = values.take(numpy.array(self.gidlist)-self.gid_start) # take just the values for cells on this machine 
    635             assert len(values) == len(self.gidlist) 
    636             if isinstance(self.celltype, common.StandardCellType): 
    637                 parametername = self.celltype.translate({parametername: values[0]}).keys()[0] 
    638             hoc_commands = [] 
    639             fmt = '%s.object(%s).%s = %s' % (self.hoc_label, "%d", parametername, "%g") 
    640             for i, val in enumerate(values): 
    641                 try: 
    642                     hoc_commands += [fmt % (i, val), 
    643                                      'success = %s.object(%d).param_update()' % (self.hoc_label, i)] 
    644                 except TypeError: 
    645                     raise common.InvalidParameterValueError, "%s is not a numeric value" % str(val) 
    646             hoc_execute(hoc_commands, "--- Population[%s].__tset()__ ---" %self.label) 
    647         else: 
    648             raise common.InvalidDimensionsError 
     299        if self.dim == value_array.shape: # the values are numbers or non-array objects 
     300            local_values = value_array[self._mask_local] 
     301        elif len(value_array.shape) == len(self.dim)+1: # the values are themselves 1D arrays 
     302            #values = numpy.reshape(value_array, (self.dim, value_array.size/self.size)) 
     303            local_values = value_array[self._mask_local] # not sure this works 
     304        else: 
     305            raise common.InvalidDimensionsError, "Population: %s, value_array: %s" % (str(self.dim), 
     306                                                                                      str(value_array.shape)) 
     307        ##local_indices = numpy.array([cell.gid for cell in self]) - self.first_id 
     308        ##values = values.take(local_indices) # take just the values for cells on this machine 
     309        assert local_values.size == self._local_ids.size, "%d != %d" % (local_values.size, self._local_ids.size) 
     310         
     311        logging.info("%s.tset('%s', array(shape=%s, min=%s, max=%s)", 
     312                     self.label, parametername, value_array.shape, 
     313                     value_array.min(), value_array.max()) 
     314        # Set the values for each cell 
     315        for cell, val in zip(self, local_values.flat): 
     316            setattr(cell, parametername, val) 
    649317 
    650318    def rset(self, parametername, rand_distr): 
    …… 
    653321        rand_distr, which should be a RandomDistribution object. 
    654322        """ 
     323        # Note that we generate enough random numbers for all cells on all nodes 
     324        # but use only those relevant to this node. This ensures that the 
     325        # sequence of random numbers does not depend on the number of nodes, 
     326        # provided that the same rng with the same seed is used on each node. 
    655327        if isinstance(rand_distr.rng, NativeRNG): 
    656             if isinstance(self.celltype, common.StandardCellType): 
    657                 parametername = self.celltype.translate({parametername: 0}).keys()[0] 
    658             paramfmt = "%g,"*len(rand_distr.parameters); paramfmt = paramfmt.strip(',') 
    659             distr_params = paramfmt % tuple(rand_distr.parameters) 
    660             hoc_commands = ['rng = new Random(%d)' % 0 or distribution.rng.seed, 
    661                             'tmp = rng.%s(%s)' % (rand_distr.name, distr_params)] 
    662             # We do the loop in hoc, to speed up the code 
    663             loop = "for tmp = 0, %d" %(len(self.gidlist)-1) 
    664             cmd = '%s.object(tmp).%s = rng.repick()' % (self.hoc_label, parametername) 
    665             hoc_commands += ['cmd="%s { %s success = %s.object(tmp).param_update()}"' %(loop, cmd, self.hoc_label), 
    666                              'success = execute1(cmd)'] 
    667             hoc_execute(hoc_commands, "--- Population[%s].__rset()__ ---" %self.label)    
    668         else: 
    669             rarr = rand_distr.next(n=self.size) 
    670             rarr = rarr.reshape(self.dim) 
    671             hoc_comment("--- Population[%s].__rset()__ --- " %self.label) 
    672             self.tset(parametername, rarr) 
    673  
    674     def _call(self, methodname, arguments): 
    675         """ 
    676         Calls the method methodname(arguments) for every cell in the population. 
    677         e.g. p.call("set_background","0.1") if the cell class has a method 
    678         set_background(). 
    679         """ 
    680         raise Exception("Method not yet implemented") 
    681         ## Not sure this belongs in the API, because cell classes only have 
    682         ## parameters/attributes, not methods. 
    683  
    684     def _tcall(self, methodname, objarr): 
    685         """ 
    686         `Topographic' call. Calls the method methodname() for every cell in the  
    687         population. The argument to the method depends on the coordinates of the 
    688         cell. objarr is an array with the same dimensions as the Population. 
    689         e.g. p.tcall("memb_init", vinitArray) calls 
    690         p.cell[i][j].memb_init(vInitArray[i][j]) for all i, j. 
    691         """ 
    692         raise Exception("Method not yet implemented")     
     328            rng = h.Random(rand_distr.rng.seed or 0) 
     329            native_rand_distr = getattr(rng, rand_distr.name) 
     330            rarr = [native_rand_distr(*rand_distr.parameters)] + [rng.repick() for i in range(self._all_ids.size-1)] 
     331        else: 
     332            rarr = rand_distr.next(n=self._all_ids.size) 
     333        rarr = numpy.array(rarr) 
     334        logging.info("%s.rset('%s', %s)", self.label, parametername, rand_distr) 
     335        for cell,val in zip(self, rarr[self._mask_local.flatten()]): 
     336            setattr(cell, parametername, val) 
     337 
     338    def randomInit(self, rand_distr): 
     339        """ 
     340        Set initial membrane potentials for all the cells in the population to 
     341        random values. 
     342        """ 
     343        self.rset('v_init', rand_distr) 
    693344 
    694345    def __record(self, record_what, record_from=None, rng=None): 
    …… 
    697348        """ 
    698349        global myid 
    699         hoc_commands = [] 
    700350        fixed_list=False 
    701  
    702351        if isinstance(record_from, list): #record from the fixed list specified by user 
    703352            fixed_list=True 
    704353        elif record_from is None: # record from all cells: