[med-svn] [python-mne] 05/376: first working version of setup_read_raw
Yaroslav Halchenko
debian at onerussian.com
Fri Nov 27 17:21:55 UTC 2015
This is an automated email from the git hooks/post-receive script.
yoh pushed a commit to annotated tag v0.1
in repository python-mne.
commit 621ae3168b64f22220ea0996b398df91f834fe9b
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date: Tue Dec 28 12:45:35 2010 -0500
first working version of setup_read_raw
---
examples/read_evoked.py | 8 +-
examples/read_raw.py | 14 ++
fiff/__init__.py | 4 +-
fiff/evoked.py | 244 +---------------------------------
fiff/{evoked.py => meas_info.py} | 211 +----------------------------
fiff/open.py | 6 +-
fiff/proj.py | 17 +--
fiff/raw.py | 280 +++++++++++++++++++++++++++++++++++++++
fiff/tag.py | 21 ++-
fiff/tree.py | 4 +-
10 files changed, 335 insertions(+), 474 deletions(-)
diff --git a/examples/read_evoked.py b/examples/read_evoked.py
index 6c2c402..95c48c4 100644
--- a/examples/read_evoked.py
+++ b/examples/read_evoked.py
@@ -20,6 +20,8 @@ fname = 'sm02a1-ave.fif'
data = fiff.read_evoked(fname)
-# import pylab as pl
-# pl.plot(data['evoked']['times'], data['evoked']['epochs'][:306,:].T)
-# pl.show()
+import pylab as pl
+pl.plot(data['evoked']['times'], data['evoked']['epochs'][:306,:].T)
+pl.xlabel('time (ms)')
+pl.ylabel('MEG data (T)')
+pl.show()
diff --git a/examples/read_raw.py b/examples/read_raw.py
new file mode 100644
index 0000000..c868296
--- /dev/null
+++ b/examples/read_raw.py
@@ -0,0 +1,14 @@
+import fiff
+
+fname = 'sm02a5_raw.fif'
+
+# fid, tree, directory = fiff.fiff_open(fname, verbose=True)
+
+raw = fiff.setup_read_raw(fname)
+# data, times = fiff.read_raw_segment(raw, from_=None, to=None, sel=None)
+
+# import pylab as pl
+# pl.plot(data['evoked']['times'], data['evoked']['epochs'][:306,:].T)
+# pl.show()
+
+
diff --git a/fiff/__init__.py b/fiff/__init__.py
index b32fefb..1984433 100644
--- a/fiff/__init__.py
+++ b/fiff/__init__.py
@@ -1,4 +1,6 @@
+from .constants import FIFF
from .open import fiff_open
from .evoked import read_evoked
from .cov import read_cov
-from .constants import FIFF
+from .raw import setup_read_raw, read_raw_segment
+
diff --git a/fiff/evoked.py b/fiff/evoked.py
index d9048dd..09cd9f3 100644
--- a/fiff/evoked.py
+++ b/fiff/evoked.py
@@ -1,250 +1,10 @@
import numpy as np
-from .bunch import Bunch
from .constants import FIFF
from .open import fiff_open
-from .ctf import read_ctf_comp
-from .tag import read_tag, find_tag
+from .tag import read_tag
from .tree import dir_tree_find
-from .proj import read_proj
-from .channels import read_bad_channels
-
-
-def read_meas_info(source, tree=None):
- """[info,meas] = fiff_read_meas_info(source,tree)
-
- Read the measurement info
-
- If tree is specified, source is assumed to be an open file id,
- otherwise a the name of the file to read. If tree is missing, the
- meas output argument should not be specified.
- """
- if tree is None:
- fid, tree, _ = fiff_open(source)
- open_here = True
- else:
- fid = source
- open_here = False
-
- # Find the desired blocks
- meas = dir_tree_find(tree, FIFF.FIFFB_MEAS)
- if len(meas) == 0:
- if open_here:
- fid.close()
- raise ValueError, 'Could not find measurement data'
- if len(meas) > 1:
- if open_here:
- fid.close()
- raise ValueError, 'Cannot read more that 1 measurement data'
- meas = meas[0]
-
- meas_info = dir_tree_find(meas, FIFF.FIFFB_MEAS_INFO)
- if len(meas_info) == 0:
- if open_here:
- fid.close()
- raise ValueError, 'Could not find measurement info'
- if len(meas_info) > 1:
- if open_here:
- fid.close()
- raise ValueError, 'Cannot read more that 1 measurement info'
- meas_info = meas_info[0]
-
- # Read measurement info
- dev_head_t = None
- ctf_head_t = None
- meas_date = None
- highpass = None
- lowpass = None
- nchan = None
- sfreq = None
- chs = []
- p = 0
- for k in range(meas_info.nent):
- kind = meas_info.directory[k].kind
- pos = meas_info.directory[k].pos
- if kind == FIFF.FIFF_NCHAN:
- tag = read_tag(fid, pos)
- nchan = tag.data
- elif kind == FIFF.FIFF_SFREQ:
- tag = read_tag(fid, pos)
- sfreq = tag.data
- elif kind == FIFF.FIFF_CH_INFO:
- tag = read_tag(fid, pos)
- chs.append(tag.data)
- p += 1
- elif kind == FIFF.FIFF_LOWPASS:
- tag = read_tag(fid, pos)
- lowpass = tag.data
- elif kind == FIFF.FIFF_HIGHPASS:
- tag = read_tag(fid, pos)
- highpass = tag.data
- elif kind == FIFF.FIFF_MEAS_DATE:
- tag = read_tag(fid, pos)
- meas_date = tag.data
- elif kind == FIFF.FIFF_COORD_TRANS:
- tag = read_tag(fid, pos)
- cand = tag.data
- if cand.from_ == FIFF.FIFFV_COORD_DEVICE and \
- cand.to == FIFF.FIFFV_COORD_HEAD: # XXX : from
- dev_head_t = cand
- elif cand.from_ == FIFF.FIFFV_MNE_COORD_CTF_HEAD and \
- cand.to == FIFF.FIFFV_COORD_HEAD:
- ctf_head_t = cand
-
- # XXX : fix
- # Check that we have everything we need
- # if ~exist('nchan','var')
- # if open_here
- # fclose(fid);
- # end
- # error(me,'Number of channels in not defined');
- # end
- # if ~exist('sfreq','var')
- # if open_here
- # fclose(fid);
- # end
- # error(me,'Sampling frequency is not defined');
- # end
- # if ~exist('chs','var')
- # if open_here
- # fclose(fid);
- # end
- # error(me,'Channel information not defined');
- # end
- # if length(chs) ~= nchan
- # if open_here
- # fclose(fid);
- # end
- # error(me,'Incorrect number of channel definitions found');
- # end
-
- if dev_head_t is None or ctf_head_t is None:
- hpi_result = dir_tree_find(meas_info, FIFF.FIFFB_HPI_RESULT)
- if len(hpi_result) == 1:
- hpi_result = hpi_result[0]
- for k in range(hpi_result.nent):
- kind = hpi_result.directory[k].kind
- pos = hpi_result.directory[k].pos
- if kind == FIFF.FIFF_COORD_TRANS:
- tag = read_tag(fid, pos)
- cand = tag.data;
- if cand.from_ == FIFF.FIFFV_COORD_DEVICE and \
- cand.to == FIFF.FIFFV_COORD_HEAD: # XXX: from
- dev_head_t = cand;
- elif cand.from_ == FIFF.FIFFV_MNE_COORD_CTF_HEAD and \
- cand.to == FIFF.FIFFV_COORD_HEAD:
- ctf_head_t = cand;
-
- # Locate the Polhemus data
- isotrak = dir_tree_find(meas_info,FIFF.FIFFB_ISOTRAK)
- if len(isotrak):
- isotrak = isotrak[0]
- else:
- if len(isotrak) == 0:
- if open_here:
- fid.close()
- raise ValueError, 'Isotrak not found'
- if len(isotrak) > 1:
- if open_here:
- fid.close()
- raise ValueError, 'Multiple Isotrak found'
-
- dig = []
- if len(isotrak) == 1:
- for k in range(isotrak.nent):
- kind = isotrak.directory[k].kind;
- pos = isotrak.directory[k].pos;
- if kind == FIFF.FIFF_DIG_POINT:
- tag = read_tag(fid,pos);
- dig.append(tag.data)
- dig[-1]['coord_frame'] = FIFF.FIFFV_COORD_HEAD
-
- # Locate the acquisition information
- acqpars = dir_tree_find(meas_info, FIFF.FIFFB_DACQ_PARS);
- acq_pars = []
- acq_stim = []
- if len(acqpars) == 1:
- for k in range(acqpars.nent):
- kind = acqpars.directory[k].kind
- pos = acqpars.directory[k].pos
- if kind == FIFF.FIFF_DACQ_PARS:
- tag = read_tag(fid, pos)
- acq_pars = tag.data
- elif kind == FIFF.FIFF_DACQ_STIM:
- tag = read_tag(fid, pos)
- acq_stim = tag.data
-
- # Load the SSP data
- projs = read_proj(fid, meas_info)
-
- # Load the CTF compensation data
- comps = read_ctf_comp(fid, meas_info, chs)
-
- # Load the bad channel list
- bads = read_bad_channels(fid, meas_info)
-
- #
- # Put the data together
- #
- if tree.id is not None:
- info = dict(file_id=tree.id)
- else:
- info = dict(file_id=None)
-
- # Make the most appropriate selection for the measurement id
- if meas_info.parent_id is None:
- if meas_info.id is None:
- if meas.id is None:
- if meas.parent_id is None:
- info['meas_id'] = info.file_id
- else:
- info['meas_id'] = meas.parent_id
- else:
- info['meas_id'] = meas.id
- else:
- info['meas_id'] = meas_info.id
- else:
- info['meas_id'] = meas_info.parent_id;
-
- if meas_date is None:
- info['meas_date'] = [info['meas_id']['secs'], info['meas_id']['usecs']]
- else:
- info['meas_date'] = meas_date
-
- info['nchan'] = nchan
- info['sfreq'] = sfreq
- info['highpass'] = highpass if highpass is not None else 0
- info['lowpass'] = lowpass if lowpass is not None else info['sfreq']/2.0
-
- # Add the channel information and make a list of channel names
- # for convenience
- info['chs'] = chs;
- info['ch_names'] = [ch.ch_name for ch in chs]
-
- #
- # Add the coordinate transformations
- #
- info['dev_head_t'] = dev_head_t
- info['ctf_head_t'] = ctf_head_t
- if dev_head_t is not None and ctf_head_t is not None:
- info['dev_ctf_t'] = info['dev_head_t']
- info['dev_ctf_t'].to = info['ctf_head_t'].from_ # XXX : see if better name
- info['dev_ctf_t'].trans = np.dot(np.inv(ctf_head_t.trans), info.dev_ctf_t.trans)
- else:
- info['dev_ctf_t'] = []
-
- # All kinds of auxliary stuff
- info['dig'] = dig
- info['bads'] = bads
- info['projs'] = projs
- info['comps'] = comps
- info['acq_pars'] = acq_pars
- info['acq_stim'] = acq_stim
-
- if open_here:
- fid.close()
-
- return info, meas
+from .meas_info import read_meas_info
def read_evoked(fname, setno=0):
diff --git a/fiff/evoked.py b/fiff/meas_info.py
similarity index 54%
copy from fiff/evoked.py
copy to fiff/meas_info.py
index d9048dd..2f6aa9d 100644
--- a/fiff/evoked.py
+++ b/fiff/meas_info.py
@@ -1,12 +1,11 @@
import numpy as np
-from .bunch import Bunch
-from .constants import FIFF
-from .open import fiff_open
-from .ctf import read_ctf_comp
-from .tag import read_tag, find_tag
from .tree import dir_tree_find
+from .open import fiff_open
+from .constants import FIFF
+from .tag import read_tag
from .proj import read_proj
+from .ctf import read_ctf_comp
from .channels import read_bad_channels
@@ -64,7 +63,7 @@ def read_meas_info(source, tree=None):
pos = meas_info.directory[k].pos
if kind == FIFF.FIFF_NCHAN:
tag = read_tag(fid, pos)
- nchan = tag.data
+ nchan = int(tag.data)
elif kind == FIFF.FIFF_SFREQ:
tag = read_tag(fid, pos)
sfreq = tag.data
@@ -164,6 +163,7 @@ def read_meas_info(source, tree=None):
acq_pars = []
acq_stim = []
if len(acqpars) == 1:
+ acqpars = acqpars[0]
for k in range(acqpars.nent):
kind = acqpars.directory[k].kind
pos = acqpars.directory[k].pos
@@ -245,202 +245,3 @@ def read_meas_info(source, tree=None):
fid.close()
return info, meas
-
-
-def read_evoked(fname, setno=0):
- """
- [data] = fiff_read_evoked(fname,setno)
-
- Read one evoked data set
-
- """
-
- if setno < 0:
- raise ValueError, 'Data set selector must be positive'
-
- print 'Reading %s ...\n' % fname
- fid, tree, _ = fiff_open(fname);
-
- # Read the measurement info
- info, meas = read_meas_info(fid, tree)
- info['filename'] = fname
-
- # Locate the data of interest
- processed = dir_tree_find(meas,FIFF.FIFFB_PROCESSED_DATA);
- if len(processed) == 0:
- fid.close()
- raise ValueError, 'Could not find processed data'
-
- evoked = dir_tree_find(meas,FIFF.FIFFB_EVOKED)
- if len(evoked) == 0:
- fid.close(fid)
- raise ValueError, 'Could not find evoked data'
-
- # Identify the aspects
- #
- naspect = 0
- is_smsh = None
- sets = []
- for k in range(len(evoked)):
- aspects_k = dir_tree_find(evoked[k], FIFF.FIFFB_ASPECT)
- set_k = dict(aspects=aspects_k,
- naspect=len(aspects_k))
- sets.append(set_k)
-
- if sets[k]['naspect'] > 0:
- if is_smsh is None:
- is_smsh = np.zeros((1, sets[k]['naspect']))
- else:
- is_smsh = np.c_[is_smsh, np.zeros((1, sets[k]['naspect']))]
- naspect += sets[k]['naspect']
-
- saspects = dir_tree_find(evoked[k], FIFF.FIFFB_SMSH_ASPECT)
- nsaspects = len(saspects)
- if nsaspects > 0:
- sets[k]['naspect'] += nsaspects
- sets[k]['naspect'] = [sets[k]['naspect'], saspects] # XXX : potential bug
- is_smsh = np.c_[is_smsh, np.ones(1, sets[k]['naspect'])]
- naspect += nsaspects;
-
- print '\t%d evoked data sets containing a total of %d data aspects' \
- ' in %s\n' % (len(evoked), naspect, fname)
-
- if setno >= naspect or setno < 0:
- fid.close()
- raise ValueError, 'Data set selector out of range'
-
- # Next locate the evoked data set
- #
- p = 0
- goon = True
- for k in range(len(evoked)):
- for a in range(sets[k]['naspect']):
- if p == setno:
- my_evoked = evoked[k]
- my_aspect = sets[k]['aspects'][a]
- goon = False
- break
- p += 1
- if not goon:
- break
- else:
- # The desired data should have been found but better to check
- fid.close(fid)
- raise ValueError, 'Desired data set not found'
-
- # Now find the data in the evoked block
- nchan = 0
- sfreq = -1
- q = 0
- chs = []
- comment = None
- for k in range(my_evoked.nent):
- kind = my_evoked.directory[k].kind
- pos = my_evoked.directory[k].pos
- if kind == FIFF.FIFF_COMMENT:
- tag = read_tag(fid, pos)
- comment = tag.data
- elif kind == FIFF.FIFF_FIRST_SAMPLE:
- tag = read_tag(fid, pos)
- first = tag.data
- elif kind == FIFF.FIFF_LAST_SAMPLE:
- tag = read_tag(fid, pos)
- last = tag.data
- elif kind == FIFF.FIFF_NCHAN:
- tag = read_tag(fid, pos)
- nchan = tag.data
- elif kind == FIFF.FIFF_SFREQ:
- tag = read_tag(fid, pos)
- sfreq = tag.data
- elif kind ==FIFF.FIFF_CH_INFO:
- tag = read_tag(fid, pos)
- chs.append(tag.data)
- q += 1
-
- if comment is None:
- comment = 'No comment'
-
- # Local channel information?
- if nchan > 0:
- if chs is None:
- fid.close()
- raise ValueError, 'Local channel information was not found when it was expected.'
-
- if len(chs) != nchan:
- fid.close()
- raise ValueError, 'Number of channels and number of channel definitions are different'
-
- info.chs = chs
- info.nchan = nchan
- print '\tFound channel information in evoked data. nchan = %d\n' % nchan
- if sfreq > 0:
- info['sfreq'] = sfreq
-
- nsamp = last - first + 1
- print '\tFound the data of interest:\n'
- print '\t\tt = %10.2f ... %10.2f ms (%s)\n' % (
- 1000*first/info['sfreq'], 1000*last/info['sfreq'], comment)
- if info['comps'] is not None:
- print '\t\t%d CTF compensation matrices available\n' % len(info['comps'])
-
- # Read the data in the aspect block
- nave = 1
- epoch = []
- for k in range(my_aspect.nent):
- kind = my_aspect.directory[k].kind
- pos = my_aspect.directory[k].pos
- if kind == FIFF.FIFF_COMMENT:
- tag = read_tag(fid, pos)
- comment = tag.data
- elif kind == FIFF.FIFF_ASPECT_KIND:
- tag = read_tag(fid, pos)
- aspect_kind = tag.data
- elif kind == FIFF.FIFF_NAVE:
- tag = read_tag(fid, pos)
- nave = tag.data
- elif kind == FIFF.FIFF_EPOCH:
- tag = read_tag(fid, pos)
- epoch.append(tag)
-
- print '\t\tnave = %d aspect type = %d\n' % (nave, aspect_kind)
-
- nepoch = len(epoch)
- if nepoch != 1 and nepoch != info.nchan:
- fid.close()
- raise ValueError, 'Number of epoch tags is unreasonable '\
- '(nepoch = %d nchan = %d)' % (nepoch, info.nchan)
-
- if nepoch == 1:
- # Only one epoch
- all_data = epoch[0].data
- # May need a transpose if the number of channels is one
- if all_data.shape[1] == 1 and info.nchan == 1:
- all_data = all_data.T
-
- else:
- # Put the old style epochs together
- all_data = epoch[0].data.T
- for k in range(1, nepoch):
- all_data = np.r_[all_data, epoch[k].data.T]
-
- if all_data.shape[1] != nsamp:
- fid.close()
- raise ValueError, 'Incorrect number of samples (%d instead of %d)' % (
- all_data.shape[1], nsamp)
-
- # Calibrate
- cals = np.array([info['chs'][k].cal for k in range(info['nchan'])])
- all_data = np.dot(np.diag(cals.ravel()), all_data)
-
- # Put it all together
- data = dict(info=info, evoked=dict(aspect_kind=aspect_kind,
- is_smsh=is_smsh[setno],
- nave=nave, first=first,
- last=last, comment=comment,
- times=np.arange(first, last+1,
- dtype=np.float) / info['sfreq'],
- epochs=all_data))
-
- fid.close()
-
- return data
diff --git a/fiff/open.py b/fiff/open.py
index 6da7123..a704419 100644
--- a/fiff/open.py
+++ b/fiff/open.py
@@ -1,4 +1,4 @@
-from .read_tag import read_tag_info, read_tag
+from .tag import read_tag_info, read_tag
from .tree import make_dir_tree
from .constants import FIFF
@@ -40,7 +40,9 @@ def fiff_open(fname, verbose=False):
directory = list()
while tag.next >= 0:
pos = fid.tell()
- directory.append(read_tag_info(fid))
+ tag = read_tag_info(fid)
+ tag.pos = pos
+ directory.append(tag)
tree, _ = make_dir_tree(fid, directory)
diff --git a/fiff/proj.py b/fiff/proj.py
index 089d07b..1cdfcb2 100644
--- a/fiff/proj.py
+++ b/fiff/proj.py
@@ -21,7 +21,7 @@ def read_proj(fid, node):
tag = find_tag(fid, nodes[0], FIFF.FIFF_NCHAN)
if tag is not None:
- global_nchan = tag.data
+ global_nchan = int(tag.data)
items = dir_tree_find(nodes[0], FIFF.FIFFB_PROJ_ITEM)
for i in range(len(items)):
@@ -30,7 +30,7 @@ def read_proj(fid, node):
item = items[i]
tag = find_tag(fid, item, FIFF.FIFF_NCHAN)
if tag is not None:
- nchan = tag.data
+ nchan = int(tag.data)
else:
nchan = global_nchan
@@ -91,14 +91,15 @@ def read_proj(fid, node):
projdata.append(one)
if len(projdata) > 0:
- print '\tRead a total of %d projection items:\n', len(projdata)
+ print '\tRead a total of %d projection items:' % len(projdata)
for k in range(len(projdata)):
- print '\t\t%s (%d x %d)' % (projdata[k].desc,
- projdata[k].data.nrow,
- projdata[k].data.ncol)
if projdata[k].active:
- print ' active\n'
+ misc = 'active'
else:
- print ' idle\n'
+ misc = ' idle'
+ print '\t\t%s (%d x %d) %s' % (projdata[k].desc,
+ projdata[k].data.nrow,
+ projdata[k].data.ncol,
+ misc)
return projdata
diff --git a/fiff/raw.py b/fiff/raw.py
new file mode 100644
index 0000000..22f4ccb
--- /dev/null
+++ b/fiff/raw.py
@@ -0,0 +1,280 @@
+import numpy as np
+
+from .constants import FIFF
+from .open import fiff_open
+from .evoked import read_meas_info
+from .tree import dir_tree_find
+from .tag import read_tag
+
+
+def setup_read_raw(fname, allow_maxshield=False):
+ """
+ %
+ % [data] = setup_read_raw(fname,allow_maxshield)
+ %
+ % Read information about raw data file
+ %
+ % fname Name of the file to read
+ % allow_maxshield Accept unprocessed MaxShield data
+ %
+ """
+
+ # Open the file
+ print 'Opening raw data file %s...' % fname
+ fid, tree, _ = fiff_open(fname)
+
+ # Read the measurement info
+ info, meas = read_meas_info(fid, tree)
+
+ # Locate the data of interest
+ raw = dir_tree_find(meas, FIFF.FIFFB_RAW_DATA)
+ if raw is None:
+ raw = dir_tree_find(meas, FIFF.FIFFB_CONTINUOUS_DATA)
+ if allow_maxshield:
+ raw = dir_tree_find(meas, FIFF.FIFFB_SMSH_RAW_DATA)
+ if raw is None:
+ raise ValueError, 'No raw data in %s' % fname
+ else:
+ if raw is None:
+ raise ValueError, 'No raw data in %s' % fname
+
+ if len(raw) == 1:
+ raw = raw[0]
+
+ # Set up the output structure
+ info['filename'] = fname
+
+ data = dict(fid=fid, info=info, first_samp=0, last_samp=0)
+
+ # Process the directory
+ directory = raw['directory']
+ nent = raw['nent']
+ nchan = int(info['nchan'])
+ first = 0
+ first_samp = 0
+ first_skip = 0
+
+ # Get first sample tag if it is there
+ if directory[first].kind == FIFF.FIFF_FIRST_SAMPLE:
+ tag = read_tag(fid, directory[first].pos)
+ first_samp = int(tag.data)
+ first += 1
+
+ # Omit initial skip
+ if directory[first].kind == FIFF.FIFF_DATA_SKIP:
+ # This first skip can be applied only after we know the buffer size
+ tag = read_tag(fid, directory[first].pos)
+ first_skip = int(tag.data)
+ first += first
+
+ data['first_samp'] = first_samp
+
+ # Go through the remaining tags in the directory
+ rawdir = list()
+ nskip = 0
+ for k in range(first, nent):
+ ent = directory[k]
+ if ent.kind == FIFF.FIFF_DATA_SKIP:
+ tag = read_tag(fid, ent.pos)
+ nskip = int(tag.data)
+ elif ent.kind == FIFF.FIFF_DATA_BUFFER:
+ # Figure out the number of samples in this buffer
+ if ent.type == FIFF.FIFFT_DAU_PACK16:
+ nsamp = ent.size / (2.0*nchan)
+ elif ent.type == FIFF.FIFFT_SHORT:
+ nsamp = ent.size / (2.0*nchan)
+ elif ent.type == FIFF.FIFFT_FLOAT:
+ nsamp = ent.size / (4.0*nchan)
+ elif ent.type == FIFF.FIFFT_INT:
+ nsamp = ent.size / (4.0*nchan)
+ else:
+ fid.close()
+ raise ValueError, 'Cannot handle data buffers of type %d' % ent.type
+
+ # Do we have an initial skip pending?
+ if first_skip > 0:
+ first_samp += nsamp*first_skip
+ data['first_samp'] = first_samp
+ first_skip = 0
+
+ # Do we have a skip pending?
+ if nskip > 0:
+ rawdir.append(dict(ent=None, first=first_samp,
+ last=first_samp + nskip*nsamp - 1,
+ nsamp=nskip*nsamp))
+ first_samp += nskip*nsamp
+ nskip = 0
+
+ # Add a data buffer
+ rawdir.append(dict(ent=ent, first=first_samp,
+ last=first_samp + nsamp -1,
+ nsamp=nsamp))
+ first_samp += nsamp
+
+ data['last_samp'] = first_samp - 1
+
+ # Add the calibration factors
+ cals = np.zeros(data['info']['nchan'])
+ for k in range(data['info']['nchan']):
+ cals[k] = data['info']['chs'][k]['range']*data['info']['chs'][k]['cal']
+
+ data['cals'] = cals
+ data['rawdir'] = rawdir
+ data['proj'] = None
+ data['comp'] = None
+ print '\tRange : %d ... %d = %9.3f ... %9.3f secs' % (
+ data['first_samp'], data['last_samp'],
+ float(data['first_samp']) / data['info']['sfreq'],
+ float(data['last_samp']) / data['info']['sfreq'])
+ print 'Ready.\n'
+
+ return data
+
+
+def read_raw_segment(raw, from_=None, to=None, sel=None):
+ """
+ %
+ % [data,times] = fiff_read_raw_segment(raw,from_,to,sel)
+ %
+ % Read a specific raw data segment
+ %
+ % raw - structure returned by fiff_setup_read_raw
+ % from_ - first sample to include. If omitted, defaults to the
+ % first sample in data
+ % to - last sample to include. If omitted, defaults to the last
+ % sample in data
+ % sel - optional channel selection vector
+ %
+ % data - returns the data matrix (channels x samples)
+ % times - returns the time values corresponding to the samples (optional)
+ %
+ """
+
+ if to is None:
+ to = raw['last_samp']
+ if from_ is None:
+ from_ = raw['first_samp']
+
+ # Initial checks
+ from_ = float(from_)
+ to = float(to)
+ if from_ < raw['first_samp']:
+ from_ = raw['first_samp']
+
+ if to > raw['last_samp']:
+ to = raw['last_samp']
+
+ if from_ > to:
+ raise ValueError, 'No data in this range'
+
+ print 'Reading %d ... %d = %9.3f ... %9.3f secs...' % (
+ from_, to, from_/raw['info']['sfreq'], to/raw['info']['sfreq'])
+
+ # Initialize the data and calibration vector
+ nchan = raw['info']['nchan']
+ dest = 1
+ cal = np.diag(raw['cals'].ravel())
+
+ if sel is None:
+ data = np.zeros((nchan, to - from_ + 1))
+ if raw['proj'] is None and raw['comp'] is None:
+ mult = None
+ else:
+ if raw['proj'] is None:
+ mult = raw['comp'] * cal
+ elif raw['comp'] is None:
+ mult = raw['proj'] * cal
+ else:
+ mult = raw['proj'] * raw['comp'] * cal
+
+ else:
+ data = np.zeros((len(sel), to - from_ + 1))
+ if raw['proj'] is None and raw['comp'] is None:
+ mult = None
+ cal = np.diag(raw['cals'][sel].ravel())
+ else:
+ if raw['proj'] is None:
+ mult = raw['comp'][sel,:] * cal
+ elif raw['comp'] is None:
+ mult = raw['proj'][sel,:]*cal
+ else:
+ mult = raw['proj'][sel,:] * raw['comp'] * cal
+
+ do_debug = False
+ if cal is not None:
+ from scipy import sparse
+ cal = sparse.csr_matrix(cal)
+
+ if mult is not None:
+ from scipy import sparse
+ mult = sparse.csr_matrix(sparse(mult))
+
+ for k in range(len(raw['rawdir'])):
+ this = raw['rawdir'][k]
+
+ # Do we need this buffer
+ if this['last'] > from_:
+ if this['ent'] is None:
+ # Take the easy route: skip is translated to zeros
+ if do_debug:
+ print 'S'
+ if sel is None:
+ one = np.zeros((nchan, this['nsamp']))
+ else:
+ one = np.zeros((len(sel), this['nsamp']))
+ else:
+ tag = read_tag(raw['fid'], this['ent'].pos)
+
+ # Depending on the state of the projection and selection
+ # we proceed a little bit differently
+ if mult is None:
+ if sel is None:
+ one = cal * tag.data.reshape(nchan, this['nsamp']).astype(np.float)
+ else:
+ one = tag.data.reshape(nchan, this['nsamp']).astype(np.float)
+ one = cal * one[sel,:]
+ else:
+ one = mult * tag.data.reshape(tag.data,nchan,this['nsamp']).astype(np.float)
+
+ # The picking logic is a bit complicated
+ if to >= this['last'] and from_ <= this['first']:
+ # We need the whole buffer
+ first_pick = 0
+ last_pick = this['nsamp']
+ if do_debug:
+ print 'W'
+
+ elif from_ > this['first']:
+ first_pick = from_ - this['first'] + 1
+ if to < this['last']:
+ # Something from the middle
+ last_pick = this['nsamp'] + to - this['last']
+ if do_debug:
+ print 'M'
+
+ else:
+ # From the middle to the end
+ last_pick = this['nsamp']
+ if do_debug:
+ print 'E'
+ else:
+ # From the beginning to the middle
+ first_pick = 1
+ last_pick = to - this['first'] + 1
+ if do_debug:
+ print 'B'
+
+ # Now we are ready to pick
+ picksamp = last_pick - first_pick + 1
+ if picksamp > 0:
+ data[:, dest:dest+picksamp-1] = one[:, first_pick:last_pick]
+ dest += picksamp
+
+ # Done?
+ if this['last'] >= to:
+ print ' [done]\n'
+ break
+
+ times = np.range(from_, to) / raw['info']['sfreq']
+
+ return data, times
diff --git a/fiff/tag.py b/fiff/tag.py
index e0144bd..991b41a 100644
--- a/fiff/tag.py
+++ b/fiff/tag.py
@@ -6,31 +6,30 @@ from .constants import FIFF
class Tag(object):
"""docstring for Tag"""
- def __init__(self, kind, type, size, next):
- self.kind = kind
- self.type = type
- self.size = size
- self.next = next
+ def __init__(self, kind, type_, size, next, pos=None):
+ self.kind = int(kind)
+ self.type = int(type_)
+ self.size = int(size)
+ self.next = int(next)
+ self.pos = pos if pos is not None else next
+ self.pos = int(self.pos)
def __repr__(self):
- out = "kind: %s - type: %s - size: %s - next: %s" % (
- self.kind, self.type, self.size, self.next)
+ out = "kind: %s - type: %s - size: %s - next: %s - pos: %s" % (
+ self.kind, self.type, self.size, self.next, self.pos)
if hasattr(self, 'data'):
out += " - data: %s\n" % self.data
else:
out += "\n"
return out
- @property
- def pos(self):
- return self.next
def read_tag_info(fid):
s = fid.read(4*4)
tag = Tag(*struct.unpack(">iiii", s))
if tag.next == 0:
fid.seek(tag.size, 1)
- else:
+ elif tag.next > 0:
fid.seek(tag.next, 0)
return tag
diff --git a/fiff/tree.py b/fiff/tree.py
index c211df2..e3ef2a0 100644
--- a/fiff/tree.py
+++ b/fiff/tree.py
@@ -1,5 +1,5 @@
from .bunch import Bunch
-from .read_tag import read_tag
+from .tag import read_tag
def dir_tree_find(tree, kind):
@@ -25,7 +25,7 @@ def dir_tree_find(tree, kind):
return nodes
-def make_dir_tree(fid, directory, start=0, indent=0, verbose=True):
+def make_dir_tree(fid, directory, start=0, indent=0, verbose=False):
"""Create the directory tree structure
"""
FIFF_BLOCK_START = 104
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-mne.git
More information about the debian-med-commit
mailing list