Changeset 413 for branches

Show
Ignore:
Timestamp:
07/31/09 09:27:57 (3 years ago)
Author:
pierre
Message:

Add some implementations performed by Luc....

Location:
branches/interval/src
Files:
2 modified

Legend:

Unmodified
Added
Removed

  • branches/interval/src/nex/nex_wrapper.py

    r409 r413  
    1212  ## WARNING : the timestamps and fragmentStarts vectors for the continuous 
    1313  ## data are not loaded. (they seem excessively NeuroExplorer specific). 
     14 
     15  ## WARNING 2 : the times are provided as seconds in NeuroExplorer format. They are 
     16  ## converted in ms for NeuroTools compatibility 
    1417 
    1518  # init the dictionnary of functions to decode the data structure  
    … …  
    3437    exit() 
    3538 
    36   data_dict = dict() 
     39  output_dict = dict() 
    3740  # read the general data at the begining of the nex file 
    38   nvar = read_basic_data(bin_file, data_dict) 
     41  nvar = read_header(bin_file, output_dict) 
    3942  # read iteratively all the variables in the file 
    4043  for i in range(nvar): 
    4144    _type = struct.unpack('<i', bin_file.read(4))[0] 
    42     header   = read_temporary_variables(bin_file) 
     45    data   = read_current_object_data(bin_file) 
    4346    filePosition = bin_file.tell() 
    4447    # in function of the data type of the variable, append a different data structure 
    4548    try : 
    46       append_data_structure[_type](bin_file, data_dict, header, filePosition) 
     49      append_data_structure[_type](bin_file, output_dict, data, filePosition) 
    4750    except KeyError : 
    4851      print 'could not append data of type : ' + str(_type)    
    … …  
    5053 
    5154  bin_file.close() 
    52   return data_dict 
     55  return output_dict 
    5356 
    54 def read_basic_data(bin_file, data_dict) : 
     57def read_header(bin_file, output_dict) : 
    5558  # read the 'global' data of the nex file 
    56   data_dict['version'] = struct.unpack('<i',bin_file.read(4))[0] 
    57   data_dict['comment'] = bin_file.read(256).replace('\x00','') 
    58   data_dict['freq']    = float(struct.unpack('<d', bin_file.read(8))[0]) 
    59   data_dict['tbeg']    = struct.unpack('<i', bin_file.read(4))[0]/data_dict['freq'] 
    60   data_dict['tend']    = struct.unpack('<i', bin_file.read(4))[0]/data_dict['freq'] 
     59  output_dict['version'] = struct.unpack('<i',bin_file.read(4))[0] 
     60  output_dict['comment'] = bin_file.read(256).replace('\x00','') 
     61  output_dict['freq']    = float(struct.unpack('<d', bin_file.read(8))[0]) 
     62  output_dict['tbeg']    = struct.unpack('<i', bin_file.read(4))[0]/output_dict['freq']*1000 
     63  output_dict['tend']    = struct.unpack('<i', bin_file.read(4))[0]/output_dict['freq']*1000 
    6164   
    6265  # number of variables in file 
    6366  nvar = struct.unpack('<i', bin_file.read(4))[0] 
    6467   
    65   # skip location of next header and padding 
     68  # skip location of next data and padding 
    6669  bin_file.seek(260, 1) 
    6770 
    6871  spikes  = [] 
    6972  id_list = [] 
    70   data_dict["neurons"]    = SpikeList(spikes, id_list, interval=(data_dict['tbeg'], data_dict['tend'])) 
    71   data_dict["events"]     = SpikeList(spikes, id_list, interval=(data_dict['tbeg'], data_dict['tend'])) 
    72   data_dict["intervals"]  = None 
    73   data_dict["waves"]      = {} 
    74   data_dict["popvectors"] = {} 
    75   data_dict["contvars"]   = {} 
    76   data_dict["markers"]    = {} 
     73  output_dict["neurons"]    = {} 
     74  output_dict["events"]     = {} 
     75  output_dict["intervals"]  = None 
     76  output_dict["waves"]      = {} 
     77  output_dict["popvectors"] = {} 
     78  output_dict["contvars"]   = {} 
     79  output_dict["markers"]    = {} 
    7780  return nvar 
    7881 
    79 def read_temporary_variables(bin_file) : 
    80   header = dict() 
     82def read_current_object_data(bin_file) : 
     83  data = dict() 
    8184   
    82   header['varVersion'] = struct.unpack('<i', bin_file.read(4))[0] 
    83   header['name'] = bin_file.read(64).replace('\x00','') 
    84   header['offset'] = struct.unpack('<i', bin_file.read(4))[0] 
    85   header['n'] = struct.unpack('<i', bin_file.read(4))[0] 
    86   header['WireNumber'] = struct.unpack('<i', bin_file.read(4))[0] 
    87   header['UnitNumber'] = struct.unpack('<i', bin_file.read(4))[0] 
    88   header['Gain'] = struct.unpack('<i', bin_file.read(4))[0] 
    89   header['Filter'] = struct.unpack('<i', bin_file.read(4))[0] 
    90   header['XPos'] = struct.unpack('<d', bin_file.read(8))[0] 
    91   header['YPos'] = struct.unpack('<d', bin_file.read(8))[0] 
    92   header['WFrequency'] = struct.unpack('<d', bin_file.read(8))[0] # wf sampling fr. 
    93   header['ADtoMV'] = struct.unpack('<d', bin_file.read(8))[0] # coeff to convert from AD values to Millivolts. 
    94   header['NPointsWave'] = struct.unpack('<i', bin_file.read(4))[0] # number of points in each wave 
    95   header['NMarkers'] = struct.unpack('<i', bin_file.read(4))[0] # how many values are associated with each marker 
    96   header['MarkerLength'] = struct.unpack('<i', bin_file.read(4))[0] # how many characters are in each marker value 
    97   header['MVOffset'] = struct.unpack('<d', bin_file.read(8))[0] # coeff to shift AD values in Millivolts: mv = raw*ADtoMV+MVOffset 
    98   return header 
     85  data['varVersion'] = struct.unpack('<i', bin_file.read(4))[0] 
     86  data['name'] = bin_file.read(64).replace('\x00','') 
     87  data['offset'] = struct.unpack('<i', bin_file.read(4))[0] 
     88  data['n'] = struct.unpack('<i', bin_file.read(4))[0] 
     89  data['WireNumber'] = struct.unpack('<i', bin_file.read(4))[0] 
     90  data['UnitNumber'] = struct.unpack('<i', bin_file.read(4))[0] 
     91  data['Gain'] = struct.unpack('<i', bin_file.read(4))[0] 
     92  data['Filter'] = struct.unpack('<i', bin_file.read(4))[0] 
     93  data['XPos'] = struct.unpack('<d', bin_file.read(8))[0] 
     94  data['YPos'] = struct.unpack('<d', bin_file.read(8))[0] 
     95  data['WFrequency'] = struct.unpack('<d', bin_file.read(8))[0] # wf sampling fr. 
     96  data['ADtoMV'] = struct.unpack('<d', bin_file.read(8))[0] # coeff to convert from AD values to Millivolts. 
     97  data['NPointsWave'] = struct.unpack('<i', bin_file.read(4))[0] # number of points in each wave 
     98  data['NMarkers'] = struct.unpack('<i', bin_file.read(4))[0] # how many values are associated with each marker 
     99  data['MarkerLength'] = struct.unpack('<i', bin_file.read(4))[0] # how many characters are in each marker value 
     100  data['MVOffset'] = struct.unpack('<d', bin_file.read(8))[0] # coeff to shift AD values in Millivolts: mv = raw*ADtoMV+MVOffset 
     101  return data 
    99102 
    100 def build_neuron(bin_file, data_dict, header, filePosition) : 
    101   bin_file.seek(header['offset'],0) 
    102   spikes = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
    103   current_ST = SpikeTrain(spikes, interval=(data_dict['tbeg'], data_dict['tend'])) 
    104   data_dict["neurons"].append(header['name'], current_ST) 
     103def build_neuron(bin_file, output_dict, data, filePosition) : 
     104  bin_file.seek(data['offset'],0) 
     105  spikes = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     106  current_ST = SpikeTrain(spikes, interval = Interval([[output_dict['tbeg'],output_dict['tend']]])) 
     107  output_dict["neurons"][data['name']] = current_ST 
    105108  bin_file.seek(filePosition,0) 
    106109 
    107 def build_event(bin_file, data_dict, header, filePosition) : 
    108   bin_file.seek(header['offset'],0) 
    109   spikes = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
    110   current_ET = SpikeTrain(spikes, interval=(data_dict['tbeg'], data_dict['tend'])) 
    111   data_dict["events"].append(header['name'],current_ET) 
     110def build_event(bin_file, output_dict, data, filePosition) : 
     111  bin_file.seek(data['offset'],0) 
     112  spikes = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     113  current_ET = SpikeTrain(spikes, interval = Interval([[output_dict['tbeg'],output_dict['tend']]])) 
     114  output_dict["events"][data['name']] = current_ET 
    112115  bin_file.seek(filePosition,0)   
    113116 
    114 def build_interval(bin_file, data_dict, header, filePosition) : 
    115   bin_file.seek(header['offset'],0) 
    116   start_times = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
    117   stop_times  = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
     117def build_interval(bin_file, output_dict, data, filePosition) : 
     118  bin_file.seek(data['offset'],0) 
     119  start_times = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     120  stop_times  = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
    118121  sub_intervals = [] 
    119122  for t_start, t_stop in zip(start_times, stop_times): 
    120123      res += [(t_start, t_stop)] 
    121   print "TOOOOOOOOO" 
    122   data_dict["intervals"].append(header['name']+'_Interval', Interval(res)) 
     124  output_dict["intervals"].append(data['name'], Interval(res)) 
    123125  bin_file.seek(filePosition,0) 
    124126 
    125 def build_waves(bin_file, data_dict, header, filePosition) : 
    126   bin_file.seek(header['offset'],0)  
    127   data_dict['waves'][header['name']]['NPointsWave'] = data_dict['NPointsWave'] 
    128   data_dict['waves'][header['name']]['WFrequency'] = WFrequency 
    129   data = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
    130   data_dict['waves'][header['name']]['timestamps'] = SpikeTrain(data) 
    131   data = np.fromstring(bin_file.read(header['n']*header['NPointsWave']*2), np.int16) 
    132   data_dict['waves'][header['name']]['waveforms'] = np.reshape(data*header['ADtoMV'] + header['MVOffset'], (header['n'],header['NPointsWave']))  
     127def build_waves(bin_file, output_dict, data, filePosition) : 
     128  bin_file.seek(data['offset'],0)  
     129  output_dict['waves'][data['name']]['NPointsWave'] = output_dict['NPointsWave'] 
     130  output_dict['waves'][data['name']]['WFrequency'] = WFrequency 
     131  data = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     132  output_dict['waves'][data['name']]['timestamps'] = SpikeTrain(data) 
     133  data = np.fromstring(bin_file.read(data['n']*data['NPointsWave']*2), np.int16) 
     134  output_dict['waves'][data['name']]['waveforms'] = np.reshape(data*data['ADtoMV'] + data['MVOffset'], (data['n'],data['NPointsWave']))  
    133135  bin_file.seek(filePosition,0) 
    134136 
    135 def build_popvectors(bin_file, data_dict, header, filePosition) : 
    136   bin_file.seek(header['offset'],0)       
    137   data_dict['popvectors'][name]['weights'] = np.fromstring(bin_file.read(header['n']*8), np.double) 
     137def build_popvectors(bin_file, output_dict, data, filePosition) : 
     138  bin_file.seek(data['offset'],0)       
     139  output_dict['popvectors'][name]['weights'] = np.fromstring(bin_file.read(data['n']*8), np.double) 
    138140  bin_file.seek(filePosition,0)       
    139141 
    140 def build_analog(bin_file, data_dict, header, filePosition) : 
    141   bin_file.seek(header['offset'],0)       
    142   timestamps = np.fromstring(bin_file.read(header['n']*4), np.int32)/data_dict['freq'] 
    143   fragmentStarts = np.fromstring(bin_file.read(header['n']*4), np.int32) 
    144   t0 = timestamps[0] - fragmentStarts[0]/float(header['WFrequency'])  
    145   data_dict['contvars'][header['name']] = AnalogSignal(np.fromstring(bin_file.read(header['NPointsWave']*2), np.int16)*header['ADtoMV'] + header['MVOffset'], 1/float(header['WFrequency']), t_start=t0); 
     142def build_analog(bin_file, output_dict, data, filePosition) : 
     143  bin_file.seek(data['offset'],0)       
     144  timestamps = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     145  fragmentStarts = np.fromstring(bin_file.read(data['n']*4), np.int32) 
     146  t0 = timestamps[0] - fragmentStarts[0]/float(data['WFrequency'])*1000. 
     147  output_dict['contvars'][data['name']] = AnalogSignal(np.fromstring(bin_file.read(data['NPointsWave']*2), np.int16)*data['ADtoMV'] + data['MVOffset'], 1/float(data['WFrequency']), t_start=t0); 
    146148  bin_file.seek(filePosition,0)       
    147149 
    148 def build_markers(bin_file, data_dict, header, filePosition) : 
    149   bin_file.seek(header['offset'],0) 
    150   data_dict['markers'][header['name']] = dict() 
    151   data_dict['markers'][header['name']]['timestamps'] = np.fromstring(bin_file.read(header['n']*4), np.int32) 
    152   data_dict['markers'][header['name']]['values'] = [dict()]*header['NMarkers'] 
    153   for j in range(header['NMarkers']) : 
    154       data_dict['markers'][header['name']]['values'][j]['name'] = bin_file.read(64).replace('\x00','') 
    155       data_dict['markers'][header['name']]['values'][j]['strings'] = ['']*header['n'] 
    156       for p in range(header['n']) : 
    157           data_dict['markers'][header['name']]['values'][j]['strings'][p] = bin_file.read(header['MarkerLength']).replace('\x00','') 
     150def build_markers(bin_file, output_dict, data, filePosition) : 
     151  bin_file.seek(data['offset'],0) 
     152  output_dict['markers'][data['name']] = dict() 
     153  output_dict['markers'][data['name']]['timestamps'] = np.fromstring(bin_file.read(data['n']*4), np.int32)/output_dict['freq']*1000 
     154  output_dict['markers'][data['name']]['values'] = [dict()]*data['NMarkers'] 
     155  for j in range(data['NMarkers']) : 
     156      output_dict['markers'][data['name']]['values'][j]['name'] = bin_file.read(64).replace('\x00','') 
     157      output_dict['markers'][data['name']]['values'][j]['strings'] = ['']*data['n'] 
     158      for p in range(data['n']) : 
     159          output_dict['markers'][data['name']]['values'][j]['strings'][p] = bin_file.read(data['MarkerLength']).replace('\x00','') 
    158160  bin_file.seek(filePosition,0) 

  • branches/interval/src/signals/intervals.py

    r409 r413  
    11from NeuroTools import check_dependency 
    2  
    32HAVE_INTERVAL = check_dependency('interval') 
    43 
    … …  
    109class Interval(object): 
    1110    """ 
    12     Interval(start_times, end_times). 
    13  
    14     Inputs: 
    15         start_times - A list of the start times for all the sub intervals considered, in ms 
    16         stop_times  - A list of the stop times for all the sub intervals considered, in ms 
     11    Interval(sub_intervals_list). 
     12 
     13    Input: 
     14        sub_intervals_list - A list of [start,stop] subintervals 
    1715     
    1816    Examples: 
    19         >> itv = Interval([0,100,200,300],[50,150,250,350]) 
     17        >> itv = Interval([[0,50],[100,150],[200,250]]) 
    2018        >> itv.time_parameters() 
    21             0, 350 
     19            0, 250 
    2220    """ 
    2321     
    … …  
    2523        """ 
    2624        Constructor of the Interval object. 
    27  
    2825        """ 
    2926        if not type(sub_intervals) is list: 
    … …  
    4441            assert (len(item) == 2), "Intervals must be a list of tuple (t_start, t_stop) !" 
    4542            assert item[1] > item[0], "Intervals must have tuple with t_start < t_stop !" 
    46  
    47     def intersect(self, itv) : 
    48         result = Interval(self.sub_intervals) 
    49         result.interval_data = result.interval_data & itv.interval_data 
    50         return result 
    51  
    52     def union(self, itv) : 
    53         if HAVE_INTERVAL: 
    54             result = Interval(self.sub_intervals) 
    55             result.interval_data = result.interval_data & itv.interval_data 
    56             result.interval_data = result.interval_data | itv.interval_data 
    57         else: 
    58             result = Interval(self.sub_intervals) 
    59             result.start_times   = min(result.start_times, spiketrain.start_times) 
    60             result.end_times     = max(result.end_times, spiketrain.end_times) 
    61         return result 
    6243     
    6344    def __str__(self): 
    … …  
    7657        return self.interval_data & interval([i,j]) 
    7758 
     59    def __getstate__(self): 
     60        """ 
     61        Replaces during pickle the interval object by pickable data 
     62        """ 
     63        #import copy 
     64        #copy_of_self = copy.deepcopy(self) 
     65        #copy_of_self.interval_data = str(self.interval_data) 
     66        #print "setstate : copy_of_self - "+str(copy_of_self.interval_data) 
     67        #print "setstate : self - "+str(self.interval_data) 
     68        #return copy_of_self 
     69        result = self.__dict__.copy() 
     70        result['interval_data'] = str(list(self.interval_data)) 
     71        return result 
     72 
     73    def __setstate__(self, dict): 
     74        """ 
     75        Is called during unpickling. Regenerates the interval 
     76        from the pickled string 
     77        """ 
     78        self.__dict__ = dict 
     79        self.interval_data = eval("interval(*%s)" %self.interval_data) 
     80 
    7881    def time_parameters(self): 
    7982        """ 
    … …  
    101104        return Interval(self.start_times, self.end_times, self.t_start, self.t_stop) 
    102105 
    103  
    104     def offset_start(self, shift, from_stop=False) : 
     106    def offset_start(self, shift, from_stop=False, colapse=False) : 
    105107        """ 
    106108        Shift all stop times of the interval by shift (ms). If from_start is 
    … …  
    110112 
    111113        if from_stop == False : 
    112             if (interval_array[:,0] + shift > interval_array[:,1]).any() : 
    113                 raise Exception("the shift is too big to preserve start/stop temporal order") 
    114             interval_array[:,0] += shift 
     114            if not colapse : 
     115                if (interval_array[:,0] + shift > interval_array[:,1]).any() : 
     116                    raise Exception("the shift is too big to preserve start/stop temporal order") 
     117                interval_array[:,0] += shift 
     118            else : 
     119                interval_array[:,0] = numpy.minimum(interval_array[:,0]+shift,interval_array[:,1]) 
    115120        else : 
    116121            if shift > 0 : 
    … …  
    121126 
    122127 
    123     def offset_stop(self, shift, from_start=False) : 
     128    def offset_stop(self, shift, from_start=False, colapse=False) : 
    124129        """ 
    125130        Shift all stop times of the interval by shift (ms). If from_start is 
    … …  
    129134 
    130135        if from_start == False : 
    131             if (interval_array[:,1] + shift < interval_array[:,0]).any() : 
    132                 raise Exception("the shift is too big to preserve start/stop temporal order") 
    133             interval_array[:,1] += shift 
     136            if not colapse : 
     137                if (interval_array[:,1] + shift < interval_array[:,0]).any() : 
     138                    raise Exception("the shift is too big to preserve start/stop temporal order") 
     139                interval_array[:,1] += shift 
     140            else : 
     141                interval_array[:,1] = numpy.maximum(interval_array[:,1]+shift,interval_array[:,0]) 
    134142        else : 
    135143            if shift < 0 : 
    … …  
    158166        if HAVE_INTERVAL: 
    159167            for itv in self.interval_data : 
    160                 spikes_selector = spikes_selector + (times > itv[0])*(times <= itv[1]) 
     168                spikes_selector = spikes_selector + (times >= itv[0])*(times <= itv[1]) 
    161169        else: 
    162170            spikes_selector = (times >= self.t_start()) & (times <= self.t_stop()) 
    … …  
    176184 
    177185 
    178 #def build_psth(spiketrain, eventtrain, before_Dt, after_Dt, intervals=None): 
    179     #"""          
    180     #build a psth of the spikes around the events  
    181      
    182     #If intervals != None, the eventtrain is restricted to the intervals provided. 
    183     #""" 
    184     ## tested : generates a correct PSTH when loaded with following data : 
    185     ##spikes = [0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,100.,101.,102.,103.,104.,105.,106.,200.,201.,202.,203.,204.] 
    186     ##  
    187     #if intervals != None: 
    188         ## keep only the events that are included in IntervalTrain object 
    189         #eventtrain = intervals.of_spikes(eventtrain) 
    190  
    191     #nRepeats = len(eventtrain.spike_times) 
    192     ## accumuate spikes around the events.  
    193     #spikes_around_event = [] 
    194     #for event_time in eventtrain.spike_times : 
    195         #spikes_around_event.extend(np.extract((spiketrain.spike_times > event_time - before_Dt)*(spiketrain.spike_times <= event_time + after_Dt), spiketrain.spike_times) - event_time) 
    196  
    197     #return np.sort(spikes_around_event, kind="quicksort"), nRepeats 
     186    def intersect(self, itv) : 
     187        result = Interval(self.sub_intervals) 
     188        result.interval_data = result.interval_data & itv.interval_data 
     189        return result 
     190 
     191    def union(self, itv) : 
     192        if HAVE_INTERVAL: 
     193            result = Interval(self.sub_intervals) 
     194            result.interval_data = result.interval_data & itv.interval_data 
     195            result.interval_data = result.interval_data | itv.interval_data 
     196        else: 
     197            result = Interval(self.sub_intervals) 
     198            result.start_times   = min(result.start_times, spiketrain.start_times) 
     199            result.end_times     = max(result.end_times, spiketrain.end_times) 
     200        return result