spike2rawio.py 28 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688
  1. """
  2. Classe for reading data in CED spike2 files (.smr).
  3. This code is based on:
  4. - sonpy, written by Antonio Gonzalez <Antonio.Gonzalez@cantab.net>
  5. Disponible here ::
  6. http://www.neuro.ki.se/broberger/
  7. and sonpy come from :
  8. - SON Library 2.0 for MATLAB, written by Malcolm Lidierth at
  9. King's College London.
  10. See http://www.kcl.ac.uk/depsta/biomedical/cfnr/lidierth.html
  11. This IO support old (<v6) and new files (>v7) of spike2
  12. Author: Samuel Garcia
  13. """
  14. # from __future__ import unicode_literals is not compatible with numpy.dtype both py2 py3
  15. from .baserawio import (BaseRawIO, _signal_channel_dtype, _unit_channel_dtype,
  16. _event_channel_dtype)
  17. import numpy as np
  18. from collections import OrderedDict
  19. class Spike2RawIO(BaseRawIO):
  20. """
  21. """
  22. extensions = ['smr']
  23. rawmode = 'one-file'
  24. def __init__(self, filename='', take_ideal_sampling_rate=False, ced_units=True,
  25. try_signal_grouping=True):
  26. BaseRawIO.__init__(self)
  27. self.filename = filename
  28. self.take_ideal_sampling_rate = take_ideal_sampling_rate
  29. self.ced_units = ced_units
  30. self.try_signal_grouping = try_signal_grouping
  31. def _parse_header(self):
  32. # get header info and channel_info
  33. with open(self.filename, 'rb') as fid:
  34. self._global_info = read_as_dict(fid, headerDescription)
  35. info = self._global_info
  36. if info['system_id'] < 6:
  37. info['dtime_base'] = 1e-6
  38. info['datetime_detail'] = 0
  39. info['datetime_year'] = 0
  40. if info['system_id'] == 9:
  41. self._diskblock = DISKBLOCK
  42. else:
  43. self._diskblock = 1
  44. self._time_factor = info['us_per_time'] * info['dtime_base']
  45. self._channel_infos = []
  46. for chan_id in range(info['channels']):
  47. fid.seek(512 + 140 * chan_id)
  48. chan_info = read_as_dict(fid, channelHeaderDesciption1)
  49. if chan_info['kind'] in [1, 6]:
  50. dt = [('scale', 'f4'), ('offset', 'f4'), ('unit', 'S6'), ]
  51. chan_info.update(read_as_dict(fid, dt))
  52. elif chan_info['kind'] in [7, 9]:
  53. dt = [('min', 'f4'), ('max', 'f4'), ('unit', 'S6'), ]
  54. chan_info.update(read_as_dict(fid, dt))
  55. elif chan_info['kind'] in [4]:
  56. dt = [('init_low', 'u1'), ('next_low', 'u1'), ]
  57. chan_info.update(read_as_dict(fid, dt))
  58. if chan_info['kind'] in [1, 6, 7, 9]:
  59. if info['system_id'] < 6:
  60. chan_info.update(read_as_dict(fid, [('divide', 'i2')]))
  61. else:
  62. chan_info.update(read_as_dict(fid, [('interleave', 'i2')]))
  63. chan_info['type'] = dict_kind[chan_info['kind']]
  64. if chan_info['blocks'] == 0:
  65. chan_info['t_start'] = 0. # this means empty signals
  66. else:
  67. fid.seek(chan_info['firstblock'] * self._diskblock)
  68. block_info = read_as_dict(fid, blockHeaderDesciption)
  69. chan_info['t_start'] = float(block_info['start_time']) * \
  70. float(info['us_per_time']) * float(info['dtime_base'])
  71. self._channel_infos.append(chan_info)
  72. # get data blocks index for all channel
  73. # run through all data block of of channel to prepare chan to block maps
  74. self._memmap = np.memmap(self.filename, dtype='u1', offset=0, mode='r')
  75. self._all_data_blocks = {}
  76. self._by_seg_data_blocks = {}
  77. for chan_id, chan_info in enumerate(self._channel_infos):
  78. data_blocks = []
  79. ind = chan_info['firstblock'] * self._diskblock
  80. for b in range(chan_info['blocks']):
  81. block_info = self._memmap[ind:ind + 20].view(blockHeaderDesciption)[0]
  82. data_blocks.append((ind, block_info['items'], 0,
  83. block_info['start_time'], block_info['end_time']))
  84. ind = block_info['succ_block'] * self._diskblock
  85. data_blocks = np.array(data_blocks, dtype=[(
  86. 'pos', 'int32'), ('size', 'int32'), ('cumsum', 'int32'),
  87. ('start_time', 'int32'), ('end_time', 'int32')])
  88. data_blocks['pos'] += 20 # 20 is ths header size
  89. self._all_data_blocks[chan_id] = data_blocks
  90. self._by_seg_data_blocks[chan_id] = []
  91. # For all signal channel detect gaps between data block (pause in rec) so new Segment.
  92. # then check that all channel have the same gaps.
  93. # this part is tricky because we need to check that all channel have same pause.
  94. all_gaps_block_ind = {}
  95. for chan_id, chan_info in enumerate(self._channel_infos):
  96. if chan_info['kind'] in [1, 9]:
  97. data_blocks = self._all_data_blocks[chan_id]
  98. sig_size = np.sum(self._all_data_blocks[chan_id]['size'])
  99. if sig_size > 0:
  100. interval = int(np.round(get_sample_interval(info,
  101. chan_info) / self._time_factor))
  102. # detect gaps
  103. inter_block_sizes = data_blocks['start_time'][1:] - \
  104. data_blocks['end_time'][:-1]
  105. gaps_block_ind, = np.nonzero(inter_block_sizes > interval)
  106. all_gaps_block_ind[chan_id] = gaps_block_ind
  107. # find t_start/t_stop for each seg based on gaps indexe
  108. self._sig_t_starts = {}
  109. self._sig_t_stops = {}
  110. if len(all_gaps_block_ind) == 0:
  111. # this means no signal channels
  112. nb_segment = 1
  113. # loop over event/spike channel to get the min/max time
  114. t_start, t_stop = None, None
  115. for chan_id, chan_info in enumerate(self._channel_infos):
  116. data_blocks = self._all_data_blocks[chan_id]
  117. if data_blocks.size > 0:
  118. # if t_start is None or data_blocks[0]['start_time']<t_start:
  119. # t_start = data_blocks[0]['start_time']
  120. if t_stop is None or data_blocks[-1]['end_time'] > t_stop:
  121. t_stop = data_blocks[-1]['end_time']
  122. # self._seg_t_starts = [t_start]
  123. self._seg_t_starts = [0]
  124. self._seg_t_stops = [t_stop]
  125. else:
  126. all_nb_seg = np.array([v.size + 1 for v in all_gaps_block_ind.values()])
  127. assert np.all(all_nb_seg[0] == all_nb_seg), \
  128. 'Signal channel have differents pause so diffrents nb_segment'
  129. nb_segment = int(all_nb_seg[0])
  130. for chan_id, gaps_block_ind in all_gaps_block_ind.items():
  131. data_blocks = self._all_data_blocks[chan_id]
  132. self._sig_t_starts[chan_id] = []
  133. self._sig_t_stops[chan_id] = []
  134. for seg_ind in range(nb_segment):
  135. if seg_ind == 0:
  136. fisrt_bl = 0
  137. else:
  138. fisrt_bl = gaps_block_ind[seg_ind - 1] + 1
  139. self._sig_t_starts[chan_id].append(data_blocks[fisrt_bl]['start_time'])
  140. if seg_ind < nb_segment - 1:
  141. last_bl = gaps_block_ind[seg_ind]
  142. else:
  143. last_bl = data_blocks.size - 1
  144. self._sig_t_stops[chan_id].append(data_blocks[last_bl]['end_time'])
  145. in_seg_data_block = data_blocks[fisrt_bl:last_bl + 1]
  146. in_seg_data_block['cumsum'][1:] = np.cumsum(in_seg_data_block['size'][:-1])
  147. self._by_seg_data_blocks[chan_id].append(in_seg_data_block)
  148. self._seg_t_starts = []
  149. self._seg_t_stops = []
  150. for seg_ind in range(nb_segment):
  151. # there is a small delay between all channel so take the max/min for t_start/t_stop
  152. t_start = min(
  153. self._sig_t_starts[chan_id][seg_ind] for chan_id in self._sig_t_starts)
  154. t_stop = max(self._sig_t_stops[chan_id][seg_ind] for chan_id in self._sig_t_stops)
  155. self._seg_t_starts.append(t_start)
  156. self._seg_t_stops.append(t_stop)
  157. # create typed channels
  158. sig_channels = []
  159. unit_channels = []
  160. event_channels = []
  161. self.internal_unit_ids = {}
  162. for chan_id, chan_info in enumerate(self._channel_infos):
  163. if chan_info['kind'] in [1, 6, 7, 9]:
  164. if self.take_ideal_sampling_rate:
  165. sampling_rate = info['ideal_rate']
  166. else:
  167. sample_interval = get_sample_interval(info, chan_info)
  168. sampling_rate = (1. / sample_interval)
  169. name = chan_info['title']
  170. if chan_info['kind'] in [1, 9]:
  171. # AnalogSignal
  172. if chan_id not in self._sig_t_starts:
  173. continue
  174. units = chan_info['unit']
  175. if chan_info['kind'] == 1: # int16
  176. gain = chan_info['scale'] / 6553.6
  177. offset = chan_info['offset']
  178. sig_dtype = 'int16'
  179. elif chan_info['kind'] == 9: # float32
  180. gain = 1.
  181. offset = 0.
  182. sig_dtype = 'float32'
  183. group_id = 0
  184. sig_channels.append((name, chan_id, sampling_rate, sig_dtype,
  185. units, gain, offset, group_id))
  186. elif chan_info['kind'] in [2, 3, 4, 5, 8]:
  187. # Event
  188. event_channels.append((name, chan_id, 'event'))
  189. elif chan_info['kind'] in [6, 7]: # SpikeTrain with waveforms
  190. wf_units = chan_info['unit']
  191. if chan_info['kind'] == 6:
  192. wf_gain = chan_info['scale'] / 6553.6
  193. wf_offset = chan_info['offset']
  194. wf_left_sweep = chan_info['n_extra'] // 4
  195. elif chan_info['kind'] == 7:
  196. wf_gain = 1.
  197. wf_offset = 0.
  198. wf_left_sweep = chan_info['n_extra'] // 8
  199. wf_sampling_rate = sampling_rate
  200. if self.ced_units:
  201. # this is a hudge pain because need
  202. # to jump over all blocks
  203. data_blocks = self._all_data_blocks[chan_id]
  204. dt = get_channel_dtype(chan_info)
  205. unit_ids = set()
  206. for bl in range(data_blocks.size):
  207. ind0 = data_blocks[bl]['pos']
  208. ind1 = data_blocks[bl]['size'] * dt.itemsize + ind0
  209. raw_data = self._memmap[ind0:ind1].view(dt)
  210. marker = raw_data['marker'] & 255
  211. unit_ids.update(np.unique(marker))
  212. unit_ids = sorted(list(unit_ids))
  213. else:
  214. # All spike from one channel are group in one SpikeTrain
  215. unit_ids = ['all']
  216. for unit_id in unit_ids:
  217. unit_index = len(unit_channels)
  218. self.internal_unit_ids[unit_index] = (chan_id, unit_id)
  219. _id = "ch{}#{}".format(chan_id, unit_id)
  220. unit_channels.append((name, _id, wf_units, wf_gain, wf_offset,
  221. wf_left_sweep, wf_sampling_rate))
  222. sig_channels = np.array(sig_channels, dtype=_signal_channel_dtype)
  223. unit_channels = np.array(unit_channels, dtype=_unit_channel_dtype)
  224. event_channels = np.array(event_channels, dtype=_event_channel_dtype)
  225. if len(sig_channels) > 0:
  226. if self.try_signal_grouping:
  227. # try to group signals channel if same sampling_rate/dtype/...
  228. # it can raise error for some files (when they do not have signal length)
  229. common_keys = ['sampling_rate', 'dtype', 'units', 'gain', 'offset']
  230. characteristics = sig_channels[common_keys]
  231. unique_characteristics = np.unique(characteristics)
  232. self._sig_dtypes = {}
  233. for group_id, charact in enumerate(unique_characteristics):
  234. chan_grp_indexes, = np.nonzero(characteristics == charact)
  235. sig_channels['group_id'][chan_grp_indexes] = group_id
  236. # check same size for channel in groups
  237. for seg_index in range(nb_segment):
  238. sig_sizes = []
  239. for ind in chan_grp_indexes:
  240. chan_id = sig_channels[ind]['id']
  241. sig_size = np.sum(self._by_seg_data_blocks[chan_id][seg_index]['size'])
  242. sig_sizes.append(sig_size)
  243. sig_sizes = np.array(sig_sizes)
  244. assert np.all(sig_sizes == sig_sizes[0]),\
  245. 'Signal channel in groups do not have same size'\
  246. ', use try_signal_grouping=False'
  247. self._sig_dtypes[group_id] = np.dtype(charact['dtype'])
  248. else:
  249. # if try_signal_grouping fail the user can try to split each channel in
  250. # separate group
  251. sig_channels['group_id'] = np.arange(sig_channels.size)
  252. self._sig_dtypes = {s['group_id']: np.dtype(s['dtype']) for s in sig_channels}
  253. # fille into header dict
  254. self.header = {}
  255. self.header['nb_block'] = 1
  256. self.header['nb_segment'] = [nb_segment]
  257. self.header['signal_channels'] = sig_channels
  258. self.header['unit_channels'] = unit_channels
  259. self.header['event_channels'] = event_channels
  260. # Annotations
  261. self._generate_minimal_annotations()
  262. bl_ann = self.raw_annotations['blocks'][0]
  263. bl_ann['system_id'] = info['system_id']
  264. seg_ann = bl_ann['segments'][0]
  265. seg_ann['system_id'] = info['system_id']
  266. for c, sig_channel in enumerate(sig_channels):
  267. chan_id = sig_channel['id']
  268. anasig_an = seg_ann['signals'][c]
  269. anasig_an['physical_channel_index'] = self._channel_infos[chan_id]['phy_chan']
  270. anasig_an['comment'] = self._channel_infos[chan_id]['comment']
  271. for c, unit_channel in enumerate(unit_channels):
  272. chan_id, unit_id = self.internal_unit_ids[c]
  273. unit_an = seg_ann['units'][c]
  274. unit_an['physical_channel_index'] = self._channel_infos[chan_id]['phy_chan']
  275. unit_an['comment'] = self._channel_infos[chan_id]['comment']
  276. for c, event_channel in enumerate(event_channels):
  277. chan_id = int(event_channel['id'])
  278. ev_an = seg_ann['events'][c]
  279. ev_an['physical_channel_index'] = self._channel_infos[chan_id]['phy_chan']
  280. ev_an['comment'] = self._channel_infos[chan_id]['comment']
  281. def _source_name(self):
  282. return self.filename
  283. def _segment_t_start(self, block_index, seg_index):
  284. return self._seg_t_starts[seg_index] * self._time_factor
  285. def _segment_t_stop(self, block_index, seg_index):
  286. return self._seg_t_stops[seg_index] * self._time_factor
  287. def _check_channel_indexes(self, channel_indexes):
  288. if channel_indexes is None:
  289. channel_indexes = slice(None)
  290. channel_indexes = np.arange(self.header['signal_channels'].size)[channel_indexes]
  291. return channel_indexes
  292. def _get_signal_size(self, block_index, seg_index, channel_indexes):
  293. channel_indexes = self._check_channel_indexes(channel_indexes)
  294. chan_id = self.header['signal_channels'][channel_indexes[0]]['id']
  295. sig_size = np.sum(self._by_seg_data_blocks[chan_id][seg_index]['size'])
  296. return sig_size
  297. def _get_signal_t_start(self, block_index, seg_index, channel_indexes):
  298. channel_indexes = self._check_channel_indexes(channel_indexes)
  299. chan_id = self.header['signal_channels'][channel_indexes[0]]['id']
  300. return self._sig_t_starts[chan_id][seg_index] * self._time_factor
  301. def _get_analogsignal_chunk(self, block_index, seg_index, i_start, i_stop, channel_indexes):
  302. if i_start is None:
  303. i_start = 0
  304. if i_stop is None:
  305. i_stop = self._get_signal_size(block_index, seg_index, channel_indexes)
  306. channel_indexes = self._check_channel_indexes(channel_indexes)
  307. chan_index = channel_indexes[0]
  308. chan_id = self.header['signal_channels'][chan_index]['id']
  309. group_id = self.header['signal_channels'][channel_indexes[0]]['group_id']
  310. dt = self._sig_dtypes[group_id]
  311. raw_signals = np.zeros((i_stop - i_start, len(channel_indexes)), dtype=dt)
  312. for c, channel_index in enumerate(channel_indexes):
  313. # NOTE: this actual way is slow because we run throught
  314. # the file for each channel. The loop should be reversed.
  315. # But there is no garanty that channels shared the same data block
  316. # indexes. So this make the job too difficult.
  317. chan_header = self.header['signal_channels'][channel_index]
  318. chan_id = chan_header['id']
  319. data_blocks = self._by_seg_data_blocks[chan_id][seg_index]
  320. # loop over data blocks and get chunks
  321. bl0 = np.searchsorted(data_blocks['cumsum'], i_start, side='right') - 1
  322. bl1 = np.searchsorted(data_blocks['cumsum'], i_stop, side='right')
  323. ind = 0
  324. for bl in range(bl0, bl1):
  325. ind0 = data_blocks[bl]['pos']
  326. ind1 = data_blocks[bl]['size'] * dt.itemsize + ind0
  327. data = self._memmap[ind0:ind1].view(dt)
  328. if bl == bl1 - 1:
  329. # right border
  330. # be carfull that bl could be both bl0 and bl1!!
  331. border = data.size - (i_stop - data_blocks[bl]['cumsum'])
  332. if border > 0:
  333. data = data[:-border]
  334. if bl == bl0:
  335. # left border
  336. border = i_start - data_blocks[bl]['cumsum']
  337. data = data[border:]
  338. raw_signals[ind:data.size + ind, c] = data
  339. ind += data.size
  340. return raw_signals
  341. def _count_in_time_slice(self, seg_index, chan_id, lim0, lim1, marker_filter=None):
  342. # count event or spike in time slice
  343. data_blocks = self._all_data_blocks[chan_id]
  344. chan_info = self._channel_infos[chan_id]
  345. dt = get_channel_dtype(chan_info)
  346. nb = 0
  347. for bl in range(data_blocks.size):
  348. ind0 = data_blocks[bl]['pos']
  349. ind1 = data_blocks[bl]['size'] * dt.itemsize + ind0
  350. raw_data = self._memmap[ind0:ind1].view(dt)
  351. ts = raw_data['tick']
  352. keep = (ts >= lim0) & (ts <= lim1)
  353. if marker_filter is not None:
  354. keep2 = (raw_data['marker'] & 255) == marker_filter
  355. keep = keep & keep2
  356. nb += np.sum(keep)
  357. if ts[-1] > lim1:
  358. break
  359. return nb
  360. def _get_internal_timestamp_(self, seg_index, chan_id,
  361. t_start, t_stop, other_field=None, marker_filter=None):
  362. chan_info = self._channel_infos[chan_id]
  363. # data_blocks = self._by_seg_data_blocks[chan_id][seg_index]
  364. data_blocks = self._all_data_blocks[chan_id]
  365. dt = get_channel_dtype(chan_info)
  366. if t_start is None:
  367. # lim0 = 0
  368. lim0 = self._seg_t_starts[seg_index]
  369. else:
  370. lim0 = int(t_start / self._time_factor)
  371. if t_stop is None:
  372. # lim1 = 2**32
  373. lim1 = self._seg_t_stops[seg_index]
  374. else:
  375. lim1 = int(t_stop / self._time_factor)
  376. timestamps = []
  377. othervalues = []
  378. for bl in range(data_blocks.size):
  379. ind0 = data_blocks[bl]['pos']
  380. ind1 = data_blocks[bl]['size'] * dt.itemsize + ind0
  381. raw_data = self._memmap[ind0:ind1].view(dt)
  382. ts = raw_data['tick']
  383. keep = (ts >= lim0) & (ts <= lim1)
  384. if marker_filter is not None:
  385. keep2 = (raw_data['marker'] & 255) == marker_filter
  386. keep = keep & keep2
  387. timestamps.append(ts[keep])
  388. if other_field is not None:
  389. othervalues.append(raw_data[other_field][keep])
  390. if ts[-1] > lim1:
  391. break
  392. if len(timestamps) > 0:
  393. timestamps = np.concatenate(timestamps)
  394. else:
  395. timestamps = np.zeros(0, dtype='int16')
  396. if other_field is None:
  397. return timestamps
  398. else:
  399. if len(timestamps) > 0:
  400. othervalues = np.concatenate(othervalues)
  401. else:
  402. othervalues = np.zeros(0, dtype=dt.fields[other_field][0])
  403. return timestamps, othervalues
  404. def _spike_count(self, block_index, seg_index, unit_index):
  405. chan_id, unit_id = self.internal_unit_ids[unit_index]
  406. if self.ced_units:
  407. marker_filter = unit_id
  408. else:
  409. marker_filter = None
  410. lim0 = self._seg_t_starts[seg_index]
  411. lim1 = self._seg_t_stops[seg_index]
  412. return self._count_in_time_slice(seg_index, chan_id,
  413. lim0, lim1, marker_filter=marker_filter)
  414. def _get_spike_timestamps(self, block_index, seg_index, unit_index, t_start, t_stop):
  415. unit_header = self.header['unit_channels'][unit_index]
  416. chan_id, unit_id = self.internal_unit_ids[unit_index]
  417. if self.ced_units:
  418. marker_filter = unit_id
  419. else:
  420. marker_filter = None
  421. spike_timestamps = self._get_internal_timestamp_(seg_index,
  422. chan_id, t_start, t_stop,
  423. marker_filter=marker_filter)
  424. return spike_timestamps
  425. def _rescale_spike_timestamp(self, spike_timestamps, dtype):
  426. spike_times = spike_timestamps.astype(dtype)
  427. spike_times *= self._time_factor
  428. return spike_times
  429. def _get_spike_raw_waveforms(self, block_index, seg_index, unit_index, t_start, t_stop):
  430. unit_header = self.header['unit_channels'][unit_index]
  431. chan_id, unit_id = self.internal_unit_ids[unit_index]
  432. if self.ced_units:
  433. marker_filter = unit_id
  434. else:
  435. marker_filter = None
  436. timestamps, waveforms = self._get_internal_timestamp_(seg_index, chan_id,
  437. t_start, t_stop,
  438. other_field='waveform',
  439. marker_filter=marker_filter)
  440. waveforms = waveforms.reshape(timestamps.size, 1, -1)
  441. return waveforms
  442. def _event_count(self, block_index, seg_index, event_channel_index):
  443. event_header = self.header['event_channels'][event_channel_index]
  444. chan_id = int(event_header['id']) # because set to string in header
  445. lim0 = self._seg_t_starts[seg_index]
  446. lim1 = self._seg_t_stops[seg_index]
  447. return self._count_in_time_slice(seg_index, chan_id, lim0, lim1, marker_filter=None)
  448. def _get_event_timestamps(self, block_index, seg_index, event_channel_index, t_start, t_stop):
  449. event_header = self.header['event_channels'][event_channel_index]
  450. chan_id = int(event_header['id']) # because set to string in header
  451. chan_info = self._channel_infos[chan_id]
  452. if chan_info['kind'] == 5:
  453. timestamps, labels = self._get_internal_timestamp_(seg_index,
  454. chan_id, t_start, t_stop,
  455. other_field='marker')
  456. elif chan_info['kind'] == 8:
  457. timestamps, labels = self._get_internal_timestamp_(seg_index,
  458. chan_id, t_start, t_stop,
  459. other_field='label')
  460. else:
  461. timestamps = self._get_internal_timestamp_(seg_index,
  462. chan_id, t_start, t_stop, other_field=None)
  463. labels = np.zeros(timestamps.size, dtype='U')
  464. labels = labels.astype('U')
  465. durations = None
  466. return timestamps, durations, labels
  467. def _rescale_event_timestamp(self, event_timestamps, dtype):
  468. event_times = event_timestamps.astype(dtype)
  469. event_times *= self._time_factor
  470. return event_times
  471. def read_as_dict(fid, dtype):
  472. """
  473. Given a file descriptor (seek at the good place externally)
  474. and a numpy.dtype of the binary struct return a dict.
  475. Make conversion for strings.
  476. """
  477. dt = np.dtype(dtype)
  478. h = np.frombuffer(fid.read(dt.itemsize), dt)[0]
  479. info = OrderedDict()
  480. for k in dt.names:
  481. v = h[k]
  482. if dt[k].kind == 'S':
  483. v = v.decode('iso-8859-1')
  484. if len(v) > 0:
  485. l = ord(v[0])
  486. v = v[1:l + 1]
  487. info[k] = v
  488. return info
  489. def get_channel_dtype(chan_info):
  490. """
  491. Get dtype by kind.
  492. """
  493. if chan_info['kind'] == 1: # Raw signal
  494. dt = 'int16'
  495. elif chan_info['kind'] in [2, 3, 4]: # Event data
  496. dt = [('tick', 'i4')]
  497. elif chan_info['kind'] in [5]: # Marker data
  498. dt = [('tick', 'i4'), ('marker', 'i4')]
  499. elif chan_info['kind'] in [6]: # AdcMark data (waveform)
  500. dt = [('tick', 'i4'), ('marker', 'i4'),
  501. # ('adc', 'S%d' % chan_info['n_extra'])]
  502. ('waveform', 'int16', chan_info['n_extra'] // 2)]
  503. elif chan_info['kind'] in [7]: # RealMark data (waveform)
  504. dt = [('tick', 'i4'), ('marker', 'i4'),
  505. # ('real', 'S%d' % chan_info['n_extra'])]
  506. ('waveform', 'float32', chan_info['n_extra'] // 4)]
  507. elif chan_info['kind'] in [8]: # TextMark data
  508. dt = [('tick', 'i4'), ('marker', 'i4'),
  509. ('label', 'S%d' % chan_info['n_extra'])]
  510. elif chan_info['kind'] == 9: # Float signal
  511. dt = 'float32'
  512. dt = np.dtype(dt)
  513. return dt
  514. def get_sample_interval(info, chan_info):
  515. """
  516. Get sample interval for one channel
  517. """
  518. if info['system_id'] in [1, 2, 3, 4, 5]: # Before version 5
  519. sample_interval = (int(chan_info['divide']) * info['us_per_time'] *
  520. info['time_per_adc']) * 1e-6
  521. else:
  522. sample_interval = (int(chan_info['l_chan_dvd']) *
  523. info['us_per_time'] * info['dtime_base'])
  524. return sample_interval
  525. # headers structures :
  526. headerDescription = [
  527. ('system_id', 'i2'),
  528. ('copyright', 'S10'),
  529. ('creator', 'S8'),
  530. ('us_per_time', 'i2'),
  531. ('time_per_adc', 'i2'),
  532. ('filestate', 'i2'),
  533. ('first_data', 'i4'), # i8
  534. ('channels', 'i2'),
  535. ('chan_size', 'i2'),
  536. ('extra_data', 'i2'),
  537. ('buffersize', 'i2'),
  538. ('os_format', 'i2'),
  539. ('max_ftime', 'i4'), # i8
  540. ('dtime_base', 'f8'),
  541. ('datetime_detail', 'u1'),
  542. ('datetime_year', 'i2'),
  543. ('pad', 'S52'),
  544. ('comment1', 'S80'),
  545. ('comment2', 'S80'),
  546. ('comment3', 'S80'),
  547. ('comment4', 'S80'),
  548. ('comment5', 'S80'),
  549. ]
  550. channelHeaderDesciption1 = [
  551. ('del_size', 'i2'),
  552. ('next_del_block', 'i4'), # i8
  553. ('firstblock', 'i4'), # i8
  554. ('lastblock', 'i4'), # i8
  555. ('blocks', 'i2'),
  556. ('n_extra', 'i2'),
  557. ('pre_trig', 'i2'),
  558. ('free0', 'i2'),
  559. ('py_sz', 'i2'),
  560. ('max_data', 'i2'),
  561. ('comment', 'S72'),
  562. ('max_chan_time', 'i4'), # i8
  563. ('l_chan_dvd', 'i4'), # i8
  564. ('phy_chan', 'i2'),
  565. ('title', 'S10'),
  566. ('ideal_rate', 'f4'),
  567. ('kind', 'u1'),
  568. ('unused1', 'i1'),
  569. ]
  570. blockHeaderDesciption = [
  571. ('pred_block', 'i4'), # i8
  572. ('succ_block', 'i4'), # i8
  573. ('start_time', 'i4'), # i8
  574. ('end_time', 'i4'), # i8
  575. ('channel_num', 'i2'),
  576. ('items', 'i2'),
  577. ]
  578. dict_kind = {
  579. 0: 'empty',
  580. 1: 'Adc',
  581. 2: 'EventFall',
  582. 3: 'EventRise',
  583. 4: 'EventBoth',
  584. 5: 'Marker',
  585. 6: 'AdcMark',
  586. 7: 'RealMark',
  587. 8: 'TextMark',
  588. 9: 'RealWave',
  589. }
  590. DISKBLOCK = 512