Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • paulmc/fslpy
  • ndcn0236/fslpy
  • seanf/fslpy
3 results
Show changes
Showing
with 2168 additions and 814 deletions
``fsl.wrappers.tbss``
=====================
.. automodule:: fsl.wrappers.tbss
:members:
:undoc-members:
:show-inheritance:
......@@ -8,24 +8,53 @@ within `FSL <https://fsl.fmrib.ox.ac.uk/fsl/fslwiki>`_ and by
|fsleyes_apidoc|_.
The top-level Python package for ``fslpy`` is called ``fsl``. It is broadly
split into the following sub-packages:
The top-level Python package for ``fslpy`` is called :mod:`fsl`. It is
broadly split into the following sub-packages:
+----------------------+-----------------------------------------------------+
| :mod:`fsl.data` | contains data abstractions and I/O routines for a |
| | range of FSL and neuroimaging file types. Most I/O |
| | routines use `nibabel <https://nipy.org/nibabel/>`_ |
| | extensively. |
+----------------------+-----------------------------------------------------+
| :mod:`fsl.utils` | contains a range of miscellaneous utilities, |
| | including :mod:`fsl.utils.path`, |
| | :mod:`fsl.utils.run`, and :mod:`fsl.utils.bids` |
+----------------------+-----------------------------------------------------+
| :mod:`fsl.scripts` | contains a range of scripts which are installed as |
| | FSL commands. |
+----------------------+-----------------------------------------------------+
| :mod:`fsl.transform` | contains functions and classes for working with |
| | FSL-style linear and non-linear transformations. |
+----------------------+-----------------------------------------------------+
| :mod:`fsl.version` | simply contains the ``fslpy`` version number. |
+----------------------+-----------------------------------------------------+
| :mod:`fsl.wrappers` | contains Python functions which can be used to |
| | invoke FSL commands. |
+----------------------+-----------------------------------------------------+
The :mod:`fsl` package provides the top-level Python package namespace for
``fslpy``, and for other FSL python libaries. It is a `native namespace
package <https://packaging.python.org/guides/packaging-namespace-packages/>`_,
which means that there is no ``fsl/__init__.py`` file.
Other libraries can use the ``fsl`` package namepace simply by also omitting a
``fsl/__init__.py`` file, and by ensuring that there are no naming conflicts
with any sub-packages of ``fslpy`` or any other projects which use the ``fsl``
package namespace.
.. autosummary::
fsl.data
fsl.utils
fsl.scripts
fsl.transform
fsl.version
fsl.wrappers
.. toctree::
:hidden:
self
fsl
fsl.data
fsl.scripts
fsl.transform
fsl.utils
fsl.wrappers
fsl.version
contributing
changelog
deprecation
dill
h5py
nibabel
nibabel.cifti2
nibabel.fileslice
nibabel.freesurfer
numpy
......@@ -8,3 +10,4 @@ numpy.linalg
scipy
scipy.ndimage
scipy.ndimage.interpolation
six
#!/usr/bin/env python
#
# __init__.py - The fslpy library.
#
# Author: Paul McCarthy <pauldmccarthy@gmail.com>
#
"""The :mod:`fsl` package is a library which contains convenience classes
and functions for use by FSL python tools. It is broadly split into the
following sub-packages:
.. autosummary::
fsl.data
fsl.utils
fsl.scripts
fsl.transform
fsl.version
fsl.wrappers
.. note:: The ``fsl`` namespace is a ``pkgutil``-style *namespace package* -
it can be used across different projects - see
https://packaging.python.org/guides/packaging-namespace-packages/
for details.
"""
__path__ = __import__('pkgutil').extend_path(__path__, __name__) # noqa
......@@ -58,7 +58,6 @@ import fsl.utils.image.resample as resample
import fsl.transform.affine as affine
import fsl.utils.notifier as notifier
import fsl.utils.settings as fslsettings
import fsl.utils.deprecated as deprecated
log = logging.getLogger(__name__)
......@@ -378,7 +377,7 @@ class AtlasLabel(object):
)
class AtlasDescription(object):
class AtlasDescription:
"""An ``AtlasDescription`` instance parses and stores the information
stored in the FSL XML file that describes a single FSL atlas. An XML
atlas specification file is assumed to have a structure that looks like
......@@ -561,7 +560,7 @@ class AtlasDescription(object):
imagefile = op.normpath(atlasDir + imagefile)
summaryimagefile = op.normpath(atlasDir + summaryimagefile)
i = fslimage.Image(imagefile, loadData=False, calcRange=False)
i = fslimage.Image(imagefile)
self.images .append(imagefile)
self.summaryImages.append(summaryimagefile)
......@@ -881,10 +880,17 @@ class LabelAtlas(Atlas):
of each present value. The proportions are returned as
values between 0 and 100.
.. note:: Calling this method will cause the atlas image data to be
loaded into memory.
.. note:: Use the :meth:`find` method to retrieve the ``AtlasLabel``
associated with each returned value.
"""
# Mask-based indexing requires the image
# data to be loaded into memory
self.data
# Extract the values that are in
# the mask, and their corresponding
# mask weights
......@@ -1085,24 +1091,6 @@ class StatisticAtlas(Atlas):
return avgvals
@deprecated.deprecated('2.6.0', '3.0.0', 'Use values instead')
def proportions(self, *args, **kwargs):
"""Deprecated - use :meth:`values` instead. """
return self.values(*args, **kwargs)
@deprecated.deprecated('2.6.0', '3.0.0', 'Use coordValues instead')
def coordProportions(self, *args, **kwargs):
"""Deprecated - use :meth:`coordValues` instead. """
return self.coordValues(*args, **kwargs)
@deprecated.deprecated('2.6.0', '3.0.0', 'Use maskValues instead')
def maskProportions(self, *args, **kwargs):
"""Deprecated - use :meth:`maskValues` instead. """
return self.maskValues(*args, **kwargs)
class ProbabilisticAtlas(StatisticAtlas):
"""A 4D atlas which contains one volume for each region. Each volume
contains probabiliy values for one region, between 0 and 100.
......
......@@ -9,20 +9,21 @@ files. Pillow is required to use the ``Bitmap`` class.
"""
import os.path as op
import logging
import six
import os.path as op
import pathlib
import logging
import numpy as np
import numpy as np
from . import image as fslimage
import fsl.data.image as fslimage
log = logging.getLogger(__name__)
BITMAP_EXTENSIONS = ['.bmp', '.png', '.jpg', '.jpeg',
'.tif', '.tiff', '.gif', '.rgba']
'.tif', '.tiff', '.gif', '.rgba',
'.jp2', '.jpg2', '.jp2k']
"""File extensions we understand. """
......@@ -34,7 +35,10 @@ BITMAP_DESCRIPTIONS = [
'TIFF',
'TIFF',
'Graphics Interchange Format',
'Raw RGBA']
'Raw RGBA',
'JPEG 2000',
'JPEG 2000',
'JPEG 2000']
"""A description for each :attr:`BITMAP_EXTENSION`. """
......@@ -51,17 +55,19 @@ class Bitmap(object):
data.
"""
if isinstance(bmp, six.string_types):
if isinstance(bmp, (pathlib.Path, str)):
try:
# Allow big images
import PIL.Image as Image
Image.MAX_IMAGE_PIXELS = 1e9
# Allow big/truncated images
import PIL.Image as Image
import PIL.ImageFile as ImageFile
Image .MAX_IMAGE_PIXELS = None
ImageFile.LOAD_TRUNCATED_IMAGES = True
except ImportError:
raise RuntimeError('Install Pillow to use the Bitmap class')
src = bmp
src = str(bmp)
img = Image.open(src)
# If this is a palette/LUT
......@@ -173,7 +179,7 @@ class Bitmap(object):
for ci, ch in enumerate(dtype.names):
data[ch] = self.data[..., ci]
data = np.array(data, order='F', copy=False)
data = np.asarray(data, order='F')
return fslimage.Image(data,
name=self.name,
......
This diff is collapsed.
......@@ -30,6 +30,7 @@ specification:
NIFTI_XFORM_ALIGNED_ANAT
NIFTI_XFORM_TALAIRACH
NIFTI_XFORM_MNI_152
NIFTI_XFORM_TEMPLATE_OTHER
"""
......@@ -81,7 +82,14 @@ NIFTI_XFORM_MNI_152 = 4
"""MNI 152 normalized coordinates."""
NIFTI_XFORM_ANALYZE = 5
NIFTI_XFORM_TEMPLATE_OTHER = 5
"""Coordinates aligned to some template that is not MNI152 or Talairach.
See https://www.nitrc.org/forum/message.php?msg_id=26394 for details.
"""
NIFTI_XFORM_ANALYZE = 6
"""Code which indicates that this is an ANALYZE image, not a NIFTI image. """
......@@ -98,6 +106,36 @@ NIFTI_UNITS_PPM = 40
NIFTI_UNITS_RADS = 48
# NIFTI datatype codes
NIFTI_DT_NONE = 0
NIFTI_DT_UNKNOWN = 0
NIFTI_DT_BINARY = 1
NIFTI_DT_UNSIGNED_CHAR = 2
NIFTI_DT_SIGNED_SHORT = 4
NIFTI_DT_SIGNED_INT = 8
NIFTI_DT_FLOAT = 16
NIFTI_DT_COMPLEX = 32
NIFTI_DT_DOUBLE = 64
NIFTI_DT_RGB = 128
NIFTI_DT_ALL = 255
NIFTI_DT_UINT8 = 2
NIFTI_DT_INT16 = 4
NIFTI_DT_INT32 = 8
NIFTI_DT_FLOAT32 = 16
NIFTI_DT_COMPLEX64 = 32
NIFTI_DT_FLOAT64 = 64
NIFTI_DT_RGB24 = 128
NIFTI_DT_INT8 = 256
NIFTI_DT_UINT16 = 512
NIFTI_DT_UINT32 = 768
NIFTI_DT_INT64 = 1024
NIFTI_DT_UINT64 = 1280
NIFTI_DT_FLOAT128 = 1536
NIFTI_DT_COMPLEX128 = 1792
NIFTI_DT_COMPLEX256 = 2048
NIFTI_DT_RGBA32 = 2304
# NIFTI file intent codes
NIFTI_INTENT_NONE = 0
NIFTI_INTENT_CORREL = 2
......
......@@ -33,15 +33,17 @@ import sys
import glob
import json
import shlex
import shutil
import logging
import binascii
import numpy as np
import nibabel as nib
import fsl.utils.tempdir as tempdir
import fsl.utils.memoize as memoize
import fsl.data.image as fslimage
import fsl.utils.tempdir as tempdir
import fsl.utils.memoize as memoize
import fsl.utils.platform as fslplatform
import fsl.data.image as fslimage
log = logging.getLogger(__name__)
......@@ -60,6 +62,25 @@ function). Versions prior to this require the series number to be passed.
"""
def dcm2niix() -> str:
"""Tries to find an absolute path to the ``dcm2niix`` command. Returns
``'dcm2niix'`` (unqualified) if a specific executable cannot be found.
"""
fsldir = fslplatform.platform.fsldir
candidates = [
shutil.which('dcm2niix')
]
if fsldir is not None:
candidates.insert(0, op.join(fsldir, 'bin', 'dcm2niix'))
for c in candidates:
if c is not None and op.exists(c):
return c
return 'dcm2niix'
class DicomImage(fslimage.Image):
"""The ``DicomImage`` is a volumetric :class:`.Image` with some associated
DICOM metadata.
......@@ -105,7 +126,7 @@ def installedVersion():
- Day
"""
cmd = 'dcm2niix -h'
cmd = f'{dcm2niix()} -h'
versionPattern = re.compile(r'v'
r'(?P<major>[0-9]+)\.'
r'(?P<minor>[0-9]+)\.'
......@@ -130,7 +151,7 @@ def installedVersion():
int(match.group('day')))
except Exception as e:
log.debug('Error parsing dcm2niix version string: {}'.format(e))
log.debug(f'Error parsing dcm2niix version string: {e}')
return None
......@@ -177,7 +198,7 @@ def scanDir(dcmdir):
raise RuntimeError('dcm2niix is not available or is too old')
dcmdir = op.abspath(dcmdir)
cmd = 'dcm2niix -b o -ba n -f %s -o . "{}"'.format(dcmdir)
cmd = f'{dcm2niix()} -b o -ba n -f %s -o . "{dcmdir}"'
series = []
with tempdir.tempdir() as td:
......@@ -194,6 +215,10 @@ def scanDir(dcmdir):
with open(fn, 'rt') as f:
meta = json.load(f)
meta['DicomDir'] = dcmdir
# SeriesDescription is not
# guaranteed to be present
if 'SeriesDescription' not in meta:
meta['SeriesDescription'] = meta['SeriesNumber']
series.append(meta)
# sort by series number
......@@ -233,7 +258,7 @@ def seriesCRC(series):
crc32 = str(binascii.crc32(uid.encode()))
if echo is not None and echo > 1:
crc32 = '{}.{}'.format(crc32, echo)
crc32 = f'{crc32}.{echo}'
return crc32
......@@ -268,14 +293,14 @@ def loadSeries(series):
else:
ident = snum
cmd = 'dcm2niix -b n -f %s -z n -o . -n "{}" "{}"'.format(ident, dcmdir)
cmd = f'{dcm2niix()} -b n -f %s -z n -o . -n "{ident}" "{dcmdir}"'
with tempdir.tempdir() as td:
with open(os.devnull, 'wb') as devnull:
sp.call(shlex.split(cmd), stdout=devnull, stderr=devnull)
files = glob.glob(op.join(td, '{}*.nii'.format(snum)))
files = glob.glob(op.join(td, f'{snum}*.nii'))
images = [nib.load(f, mmap=False) for f in files]
# copy images so nibabel no longer
......
......@@ -22,9 +22,12 @@ following functions are provided:
isFirstLevelAnalysis
loadDesign
loadContrasts
loadFTests
loadFsf
loadSettings
getThresholds
loadClusterResults
loadFEATDesignFile
The following functions return the names of various files of interest:
......@@ -38,11 +41,14 @@ The following functions return the names of various files of interest:
getPEFile
getCOPEFile
getZStatFile
getZFStatFile
getClusterMaskFile
getFClusterMaskFile
"""
import collections
import io
import logging
import os.path as op
import numpy as np
......@@ -165,70 +171,83 @@ def loadContrasts(featdir):
:arg featdir: A FEAT directory.
"""
matrix = None
numContrasts = 0
names = {}
designcon = op.join(featdir, 'design.con')
filename = op.join(featdir, 'design.con')
log.debug('Loading FEAT contrasts from {}'.format(designcon))
log.debug('Loading FEAT contrasts from %s', filename)
with open(designcon, 'rt') as f:
try:
designcon = loadFEATDesignFile(filename)
contrasts = np.genfromtxt(io.StringIO(designcon['Matrix']), ndmin=2)
numContrasts = int(designcon['NumContrasts'])
names = []
while True:
line = f.readline().strip()
if numContrasts != contrasts.shape[0]:
raise RuntimeError(f'Matrix shape {contrasts.shape} does not '
f'match number of contrasts {numContrasts}')
if line.startswith('/ContrastName'):
tkns = line.split(None, 1)
num = [c for c in tkns[0] if c.isdigit()]
num = int(''.join(num))
contrasts = [list(row) for row in contrasts]
# The /ContrastName field may not
# actually have a name specified
if len(tkns) > 1:
name = tkns[1].strip()
names[num] = name
for i in range(numContrasts):
cname = designcon.get(f'ContrastName{i + 1}', '')
if cname == '':
cname = f'{i + 1}'
names.append(cname)
elif line.startswith('/NumContrasts'):
numContrasts = int(line.split()[1])
except Exception as e:
log.debug('Error reading %s: %s', filename, e, exc_info=True)
raise RuntimeError(f'{filename} does not appear '
'to be a valid design.con file') from e
elif line == '/Matrix':
break
return names, contrasts
matrix = np.loadtxt(f, ndmin=2)
if matrix is None or \
numContrasts != matrix.shape[0]:
raise RuntimeError('{} does not appear to be a '
'valid design.con file'.format(designcon))
def loadFTests(featdir):
"""Loads F-tests from a FEAT directory. Returns a list of f-test vectors
(each of which is a list itself), where each vector contains a 1 or a 0
denoting the contrasts included in the F-test.
# Fill in any missing contrast names
if len(names) != numContrasts:
for i in range(numContrasts):
if i + 1 not in names:
names[i + 1] = str(i + 1)
:arg featdir: A FEAT directory.
"""
names = [names[c + 1] for c in range(numContrasts)]
contrasts = []
filename = op.join(featdir, 'design.fts')
for row in matrix:
contrasts.append(list(row))
if not op.exists(filename):
return []
return names, contrasts
log.debug('Loading FEAT F-tests from %s', filename)
try:
desfts = loadFEATDesignFile(filename)
ftests = np.genfromtxt(io.StringIO(desfts['Matrix']), ndmin=2)
ncols = int(desfts['NumWaves'])
nrows = int(desfts['NumContrasts'])
if ftests.shape != (nrows, ncols):
raise RuntimeError(f'Matrix shape {ftests.shape} does not match '
f'number of EVs/FTests ({ncols}, {nrows})')
def loadSettings(featdir):
"""Loads the analysis settings from a FEAT directory.
ftests = [list(row) for row in ftests]
Returns a dict containing the settings specified in the ``design.fsf``
file within the directory
except Exception as e:
log.debug('Error reading %s: %s', filename, e, exc_info=True)
raise RuntimeError(f'{filename} does not appear '
'to be a valid design.fts file') from e
:arg featdir: A FEAT directory.
return ftests
def loadFsf(designfsf):
"""Loads the analysis settings from a text file (.fsf) used to configure
FEAT.
Returns a dict containing the settings specified in the file
:arg designfsf: A .fsf file.
"""
settings = collections.OrderedDict()
designfsf = op.join(featdir, 'design.fsf')
log.debug('Loading FEAT settings from {}'.format(designfsf))
log.debug('Loading FEAT settings from %s', designfsf)
with open(designfsf, 'rt') as f:
......@@ -251,6 +270,20 @@ def loadSettings(featdir):
return settings
def loadSettings(featdir):
"""Loads the analysis settings from a FEAT directory.
Returns a dict containing the settings specified in the ``design.fsf``
file within the directory
:arg featdir: A FEAT directory.
"""
designfsf = op.join(featdir, 'design.fsf')
return loadFsf(designfsf)
def loadDesign(featdir, settings):
"""Loads the design matrix from a FEAT directory.
......@@ -296,19 +329,22 @@ def isFirstLevelAnalysis(settings):
return settings['level'] == '1'
def loadClusterResults(featdir, settings, contrast):
def loadClusterResults(featdir, settings, contrast, ftest=False):
"""If cluster thresholding was used in the FEAT analysis, this function
will load and return the cluster results for the specified (0-indexed)
contrast number.
contrast or f-test.
If there are no cluster results for the given contrast, ``None`` is
returned.
If there are no cluster results for the given contrast/f-test, ``None``
is returned.
An error will be raised if the cluster file cannot be parsed.
:arg featdir: A FEAT directory.
:arg settings: A FEAT settings dictionary.
:arg contrast: 0-indexed contrast number.
:arg contrast: 0-indexed contrast or f-test number.
:arg ftest: If ``False`` (default), return cluster results for
the contrast numbered ``contrast``. Otherwise, return
cluster results for the f-test numbered ``contrast``.
:returns: A list of ``Cluster`` instances, each of which contains
information about one cluster. A ``Cluster`` object has the
......@@ -329,11 +365,16 @@ def loadClusterResults(featdir, settings, contrast):
gravity.
``zcogz`` Z voxel coordinate of cluster centre of
gravity.
``copemax`` Maximum COPE value in cluster.
``copemaxx`` X voxel coordinate of maximum COPE value.
``copemax`` Maximum COPE value in cluster (not
present for f-tests).
``copemaxx`` X voxel coordinate of maximum COPE value
(not present for f-tests).
``copemaxy`` Y voxel coordinate of maximum COPE value.
(not present for f-tests).
``copemaxz`` Z voxel coordinate of maximum COPE value.
(not present for f-tests).
``copemean`` Mean COPE of all voxels in the cluster.
(not present for f-tests).
============ =========================================
"""
......@@ -343,8 +384,11 @@ def loadClusterResults(featdir, settings, contrast):
# the ZMax/COG etc coordinates
# are usually in voxel coordinates
coordXform = np.eye(4)
clusterFile = op.join(
featdir, 'cluster_zstat{}.txt'.format(contrast + 1))
if ftest: prefix = 'cluster_zfstat'
else: prefix = 'cluster_zstat'
clusterFile = op.join(featdir, f'{prefix}{contrast + 1}.txt')
if not op.exists(clusterFile):
......@@ -353,22 +397,16 @@ def loadClusterResults(featdir, settings, contrast):
# the cluster file will instead be called
# 'cluster_zstatX_std.txt', so we'd better
# check for that too.
clusterFile = op.join(
featdir, 'cluster_zstat{}_std.txt'.format(contrast + 1))
clusterFile = op.join(featdir, f'{prefix}{contrast + 1}_std.txt')
if not op.exists(clusterFile):
return None
# In higher levle analysis run in some standard
# In higher level analysis run in some standard
# space, the cluster coordinates are in standard
# space. We transform them to voxel coordinates.
# later on.
coordXform = fslimage.Image(
getDataFile(featdir),
loadData=False).worldToVoxMat
log.debug('Loading cluster results for contrast {} from {}'.format(
contrast, clusterFile))
coordXform = fslimage.Image(getDataFile(featdir)).worldToVoxMat
# The cluster.txt file is converted
# into a list of Cluster objects,
......@@ -386,10 +424,18 @@ def loadClusterResults(featdir, settings, contrast):
# if cluster thresholding was not used,
# the cluster table will not contain
# P valuse.
# P values.
if not hasattr(self, 'p'): self.p = 1.0
if not hasattr(self, 'logp'): self.logp = 0.0
# F-test cluster results will not have
# COPE-* results
if not hasattr(self, 'copemax'): self.copemax = np.nan
if not hasattr(self, 'copemaxx'): self.copemaxx = np.nan
if not hasattr(self, 'copemaxy'): self.copemaxy = np.nan
if not hasattr(self, 'copemaxz'): self.copemaxz = np.nan
if not hasattr(self, 'copemean'): self.copemean = np.nan
# This dict provides a mapping between
# Cluster object attribute names, and
# the corresponding column name in the
......@@ -421,10 +467,9 @@ def loadClusterResults(featdir, settings, contrast):
'COPE-MAX Z (mm)' : 'copemaxz',
'COPE-MEAN' : 'copemean'}
# An error will be raised if the
# cluster file does not exist (e.g.
# if the specified contrast index
# is invalid)
log.debug('Loading cluster results for contrast %s from %s',
contrast, clusterFile)
with open(clusterFile, 'rt') as f:
# Get every line in the file,
......@@ -446,12 +491,11 @@ def loadClusterResults(featdir, settings, contrast):
colNames = colNames.split('\t')
clusterLines = [cl .split('\t') for cl in clusterLines]
# Turn each cluster line into a
# Cluster instance. An error will
# be raised if the columm names
# are unrecognised (i.e. not in
# the colmap above), or if the
# file is poorly formed.
# Turn each cluster line into a Cluster
# instance. An error will be raised if the
# columm names are unrecognised (i.e. not
# in the colmap above), or if the file is
# poorly formed.
clusters = [Cluster(**dict(zip(colNames, cl))) for cl in clusterLines]
# Make sure all coordinates are in voxels -
......@@ -466,17 +510,51 @@ def loadClusterResults(featdir, settings, contrast):
zcog = [c.zcogx, c.zcogy, c.zcogz]
copemax = [c.copemaxx, c.copemaxy, c.copemaxz]
zmax = affine.transform([zmax], coordXform)[0].round()
zcog = affine.transform([zcog], coordXform)[0].round()
copemax = affine.transform([copemax], coordXform)[0].round()
zmax = affine.transform([zmax], coordXform)[0]
zcog = affine.transform([zcog], coordXform)[0]
copemax = affine.transform([copemax], coordXform)[0]
c.zmaxx, c.zmaxy, c.zmaxz = zmax
c.zcogx, c.zcogy, c.zcogz = zcog
c.copemax, c.copemaxy, c.copemaxz = copemax
c.zmaxx, c.zmaxy, c.zmaxz = zmax
c.zcogx, c.zcogy, c.zcogz = zcog
c.copemaxx, c.copemaxy, c.copemaxz = copemax
return clusters
def loadFEATDesignFile(filename):
"""Load a FEAT design file, e.g. ``design.mat``, ``design.con``, ``design.fts``.
These files contain key-value pairs, and are formatted according to an
undocumented structure where each key is of the form "/KeyName", and is
followed immediately by a whitespace character, and then the value.
:arg filename: File to load
:returns: A dictionary of key-value pairs. The values are all left
as strings.
"""
fields = {}
with open(filename, 'rt') as f:
content = f.read()
content = content.split('/')
for line in content:
line = line.strip()
if line == '':
continue
tokens = line.split(maxsplit=1)
if len(tokens) == 1:
name, value = tokens[0], ''
else:
name, value = tokens
fields[name] = value
return fields
def getDataFile(featdir):
"""Returns the name of the file in the FEAT directory which contains
the model input data (typically called ``filtered_func_data.nii.gz``).
......@@ -520,7 +598,7 @@ def getPEFile(featdir, ev):
:arg featdir: A FEAT directory.
:arg ev: The EV number (0-indexed).
"""
pefile = op.join(featdir, 'stats', 'pe{}'.format(ev + 1))
pefile = op.join(featdir, 'stats', f'pe{ev + 1}')
return fslimage.addExt(pefile, mustExist=True)
......@@ -532,7 +610,7 @@ def getCOPEFile(featdir, contrast):
:arg featdir: A FEAT directory.
:arg contrast: The contrast number (0-indexed).
"""
copefile = op.join(featdir, 'stats', 'cope{}'.format(contrast + 1))
copefile = op.join(featdir, 'stats', f'cope{contrast + 1}')
return fslimage.addExt(copefile, mustExist=True)
......@@ -544,10 +622,22 @@ def getZStatFile(featdir, contrast):
:arg featdir: A FEAT directory.
:arg contrast: The contrast number (0-indexed).
"""
zfile = op.join(featdir, 'stats', 'zstat{}'.format(contrast + 1))
zfile = op.join(featdir, 'stats', f'zstat{contrast + 1}')
return fslimage.addExt(zfile, mustExist=True)
def getZFStatFile(featdir, ftest):
"""Returns the path of the Z-statistic file for the specified F-test.
Raises a :exc:`~fsl.utils.path.PathError` if the file does not exist.
:arg featdir: A FEAT directory.
:arg ftest: The F-test number (0-indexed).
"""
zffile = op.join(featdir, 'stats', f'zfstat{ftest + 1}')
return fslimage.addExt(zffile, mustExist=True)
def getClusterMaskFile(featdir, contrast):
"""Returns the path of the cluster mask file for the specified contrast.
......@@ -556,5 +646,17 @@ def getClusterMaskFile(featdir, contrast):
:arg featdir: A FEAT directory.
:arg contrast: The contrast number (0-indexed).
"""
mfile = op.join(featdir, 'cluster_mask_zstat{}'.format(contrast + 1))
mfile = op.join(featdir, f'cluster_mask_zstat{contrast + 1}')
return fslimage.addExt(mfile, mustExist=True)
def getFClusterMaskFile(featdir, ftest):
"""Returns the path of the cluster mask file for the specified f-test.
Raises a :exc:`~fsl.utils.path.PathError` if the file does not exist.
:arg featdir: A FEAT directory.
:arg contrast: The f-test number (0-indexed).
"""
mfile = op.join(featdir, f'cluster_mask_zfstat{ftest + 1}')
return fslimage.addExt(mfile, mustExist=True)
......@@ -160,7 +160,7 @@ class FEATFSFDesign(object):
# Print a warning if we're
# using an old version of FEAT
if version < 6:
log.warning('Unsupported FEAT version: {}'.format(version))
log.warning('Unsupported FEAT version: %s', version)
# We need to parse the EVS a bit
# differently depending on whether
......@@ -210,8 +210,7 @@ class FEATFSFDesign(object):
continue
if (not self.__loadVoxEVs) or (ev.filename is None):
log.warning('Voxel EV image missing '
'for ev {}'.format(ev.index))
log.warning('Voxel EV image missing for ev %s', ev.index)
continue
design[:, ev.index] = ev.getData(x, y, z)
......@@ -250,8 +249,7 @@ class VoxelwiseEVMixin(object):
if op.exists(filename):
self.__filename = filename
else:
log.warning('Voxelwise EV file does not '
'exist: {}'.format(filename))
log.warning('Voxelwise EV file does not exist: %s', filename)
self.__filename = None
self.__image = None
......@@ -502,11 +500,11 @@ def getFirstLevelEVs(featDir, settings, designMat):
# - voxelwise EVs
for origIdx in range(origEVs):
title = settings[ 'evtitle{}' .format(origIdx + 1)]
shape = int(settings[ 'shape{}' .format(origIdx + 1)])
convolve = int(settings[ 'convolve{}' .format(origIdx + 1)])
deriv = int(settings[ 'deriv_yn{}' .format(origIdx + 1)])
basis = int(settings.get('basisfnum{}'.format(origIdx + 1), -1))
title = settings[ f'evtitle{origIdx + 1}']
shape = int(settings[ f'shape{origIdx + 1}'])
convolve = int(settings[ f'convolve{origIdx + 1}'])
deriv = int(settings[ f'deriv_yn{origIdx + 1}'])
basis = int(settings.get(f'basisfnum{origIdx + 1}', -1))
# Normal EV. This is just a column
# in the design matrix, defined by
......@@ -525,8 +523,7 @@ def getFirstLevelEVs(featDir, settings, designMat):
# The addExt function will raise an
# error if the file does not exist.
filename = op.join(
featDir, 'designVoxelwiseEV{}'.format(origIdx + 1))
filename = op.join(featDir, f'designVoxelwiseEV{origIdx + 1}')
filename = fslimage.addExt(filename, True)
evs.append(VoxelwiseEV(len(evs), origIdx, title, filename))
......@@ -607,7 +604,7 @@ def getFirstLevelEVs(featDir, settings, designMat):
startIdx = len(evs) + 1
if voxConfLocs != list(range(startIdx, startIdx + len(voxConfFiles))):
raise FSFError('Unsupported voxelwise confound ordering '
'({} -> {})'.format(startIdx, voxConfLocs))
f'({startIdx} -> {voxConfLocs})')
# Create the voxelwise confound EVs.
# We make a name for the EV from the
......@@ -680,7 +677,7 @@ def getHigherLevelEVs(featDir, settings, designMat):
for origIdx in range(origEVs):
# All we need is the title
title = settings['evtitle{}'.format(origIdx + 1)]
title = settings[f'evtitle{origIdx + 1}']
evs.append(NormalEV(len(evs), origIdx, title))
# Only the input file is specified for
......@@ -689,7 +686,7 @@ def getHigherLevelEVs(featDir, settings, designMat):
# name.
for origIdx in range(voxEVs):
filename = settings['evs_vox_{}'.format(origIdx + 1)]
filename = settings[f'evs_vox_{origIdx + 1}']
title = op.basename(fslimage.removeExt(filename))
evs.append(VoxelwiseEV(len(evs), origIdx, title, filename))
......@@ -705,12 +702,12 @@ def loadDesignMat(designmat):
:arg designmat: Path to the ``design.mat`` file.
"""
log.debug('Loading FEAT design matrix from {}'.format(designmat))
log.debug('Loading FEAT design matrix from %s', designmat)
matrix = np.loadtxt(designmat, comments='/', ndmin=2)
if matrix is None or matrix.size == 0 or len(matrix.shape) != 2:
raise FSFError('{} does not appear to be a '
'valid design.mat file'.format(designmat))
raise FSFError(f'{designmat} does not appear '
'to be a valid design.mat file')
return matrix
......@@ -63,8 +63,8 @@ class FEATImage(fslimage.Image):
path = op.join(path, 'filtered_func_data')
if not featanalysis.isFEATImage(path):
raise ValueError('{} does not appear to be data '
'from a FEAT analysis'.format(path))
raise ValueError(f'{path} does not appear to be '
'data from a FEAT analysis')
featDir = op.dirname(path)
settings = featanalysis.loadSettings( featDir)
......@@ -72,9 +72,11 @@ class FEATImage(fslimage.Image):
if featanalysis.hasStats(featDir):
design = featanalysis.loadDesign( featDir, settings)
names, cons = featanalysis.loadContrasts(featDir)
ftests = featanalysis.loadFTests( featDir)
else:
design = None
names, cons = [], []
ftests = []
fslimage.Image.__init__(self, path, **kwargs)
......@@ -83,26 +85,31 @@ class FEATImage(fslimage.Image):
self.__design = design
self.__contrastNames = names
self.__contrasts = cons
self.__ftests = ftests
self.__settings = settings
self.__residuals = None
self.__pes = [None] * self.numEVs()
self.__copes = [None] * self.numContrasts()
self.__zstats = [None] * self.numContrasts()
self.__zfstats = [None] * self.numFTests()
self.__clustMasks = [None] * self.numContrasts()
self.__fclustMasks = [None] * self.numFTests()
if 'name' not in kwargs:
self.name = '{}: {}'.format(self.__analysisName, self.name)
self.name = f'{self.__analysisName}: {self.name}'
def __del__(self):
"""Clears references to any loaded images."""
self.__design = None
self.__residuals = None
self.__pes = None
self.__copes = None
self.__zstats = None
self.__clustMasks = None
self.__design = None
self.__residuals = None
self.__pes = None
self.__copes = None
self.__zstats = None
self.__zfstats = None
self.__clustMasks = None
self.__fclustMasks = None
def getFEATDir(self):
......@@ -191,6 +198,11 @@ class FEATImage(fslimage.Image):
return len(self.__contrasts)
def numFTests(self):
"""Returns the number of f-tests in the analysis."""
return len(self.__ftests)
def contrastNames(self):
"""Returns a list containing the name of each contrast in the analysis.
"""
......@@ -206,6 +218,15 @@ class FEATImage(fslimage.Image):
return [list(c) for c in self.__contrasts]
def ftests(self):
"""Returns a list containing the analysis f-test vectors.
See :func:`.featanalysis.loadFTests`
"""
return [list(f) for f in self.__ftests]
def thresholds(self):
"""Returns the statistical thresholds used in the analysis.
......@@ -214,14 +235,16 @@ class FEATImage(fslimage.Image):
return featanalysis.getThresholds(self.__settings)
def clusterResults(self, contrast):
"""Returns the clusters found in the analysis.
def clusterResults(self, contrast, ftest=False):
"""Returns the clusters found in the analysis for the specified
contrast or f-test.
See :func:.featanalysis.loadClusterResults`
"""
return featanalysis.loadClusterResults(self.__featDir,
self.__settings,
contrast)
contrast,
ftest)
def getPE(self, ev):
......@@ -229,12 +252,10 @@ class FEATImage(fslimage.Image):
if self.__pes[ev] is None:
pefile = featanalysis.getPEFile(self.__featDir, ev)
evname = self.evNames()[ev]
self.__pes[ev] = fslimage.Image(
pefile,
name='{}: PE{} ({})'.format(
self.__analysisName,
ev + 1,
self.evNames()[ev]))
name=f'{self.__analysisName}: PE{ev + 1} ({evname})')
return self.__pes[ev]
......@@ -246,7 +267,7 @@ class FEATImage(fslimage.Image):
resfile = featanalysis.getResidualFile(self.__featDir)
self.__residuals = fslimage.Image(
resfile,
name='{}: residuals'.format(self.__analysisName))
name=f'{self.__analysisName}: residuals')
return self.__residuals
......@@ -256,12 +277,10 @@ class FEATImage(fslimage.Image):
if self.__copes[con] is None:
copefile = featanalysis.getCOPEFile(self.__featDir, con)
conname = self.contrastNames()[con]
self.__copes[con] = fslimage.Image(
copefile,
name='{}: COPE{} ({})'.format(
self.__analysisName,
con + 1,
self.contrastNames()[con]))
name=f'{self.__analysisName}: COPE{con + 1} ({conname})')
return self.__copes[con]
......@@ -270,35 +289,54 @@ class FEATImage(fslimage.Image):
"""
if self.__zstats[con] is None:
zfile = featanalysis.getZStatFile(self.__featDir, con)
zfile = featanalysis.getZStatFile(self.__featDir, con)
conname = self.contrastNames()[con]
self.__zstats[con] = fslimage.Image(
zfile,
name='{}: zstat{} ({})'.format(
self.__analysisName,
con + 1,
self.contrastNames()[con]))
name=f'{self.__analysisName}: zstat{con + 1} ({conname})')
return self.__zstats[con]
def getZFStats(self, ftest):
"""Returns the Z statistic image for the given f-test (0-indexed). """
if self.__zfstats[ftest] is None:
zfile = featanalysis.getZFStatFile(self.__featDir, ftest)
self.__zfstats[ftest] = fslimage.Image(
zfile,
name=f'{self.__analysisName}: zfstat{ftest + 1}')
return self.__zfstats[ftest]
def getClusterMask(self, con):
"""Returns the cluster mask image for the given contrast (0-indexed).
"""
if self.__clustMasks[con] is None:
mfile = featanalysis.getClusterMaskFile(self.__featDir, con)
mfile = featanalysis.getClusterMaskFile(self.__featDir, con)
conname = self.contrastNames()[con]
self.__clustMasks[con] = fslimage.Image(
mfile,
name='{}: cluster mask for zstat{} ({})'.format(
self.__analysisName,
con + 1,
self.contrastNames()[con]))
name=f'{self.__analysisName}: cluster mask '
f'for zstat{con + 1} ({conname})')
return self.__clustMasks[con]
def getFClusterMask(self, ftest):
"""Returns the cluster mask image for the given f-test (0-indexed).
"""
if self.__fclustMasks[ftest] is None:
mfile = featanalysis.getFClusterMaskFile(self.__featDir, ftest)
self.__fclustMasks[ftest] = fslimage.Image(
mfile,
name=f'{self.__analysisName}: cluster mask '
f'for zfstat{ftest + 1}')
return self.__fclustMasks[ftest]
def fit(self, contrast, xyz):
"""Calculates the model fit for the given contrast vector
at the given voxel. See the :func:`modelFit` function.
......
......@@ -16,18 +16,21 @@
"""
import os.path as op
import itertools as it
import math
import os.path as op
def loadLabelFile(filename,
includeLabel=None,
excludeLabel=None,
returnIndices=False):
"""Loads component labels from the specified file. The file is assuemd
returnIndices=False,
missingLabel='Unknown',
returnProbabilities=False):
"""Loads component labels from the specified file. The file is assumed
to be of the format generated by FIX, Melview or ICA-AROMA; such a file
should have a structure resembling the following::
filtered_func_data.ica
1, Signal, False
2, Unclassified Noise, True
......@@ -39,7 +42,6 @@ def loadLabelFile(filename,
8, Signal, False
[2, 5, 6, 7]
.. note:: This function will also parse files which only contain a
component list, e.g.::
......@@ -66,31 +68,46 @@ def loadLabelFile(filename,
- One or more labels for the component (multiple labels must be
comma-separated).
- ``'True'`` if the component has been classified as *bad*,
``'False'`` otherwise. This field is optional - if the last
comma-separated token on a line is not equal (case-insensitive)
to ``True`` or ``False``, it is interpreted as a component label.
- ``'True'`` if the component has been classified as *bad*, ``'False'``
otherwise. This field is optional - if the last non-numeric
comma-separated token on a line is not equal to ``True`` or ``False``
(case-insensitive) , it is interpreted as a component label.
- A value between 0 and 1, which gives the probability of the component
being signal, as generated by an automatic classifier (e.g. FIX). This
field is optional - it is output by some versions of FIX.
The last line of the file contains the index (starting from 1) of all
*bad* components, i.e. those components which are not classified as
signal or unknown.
:arg filename: Name of the label file to load.
:arg filename: Name of the label file to load.
:arg includeLabel: If the file contains a single line containing a
list component indices, this label will be used
for the components in the list. Defaults to
``'Unclassified noise'`` for FIX-like files, and
``'Movement'`` for ICA-AROMA-like files.
:arg includeLabel: If the file contains a single line containing a list
component indices, this label will be used for the
components in the list. Defaults to 'Unclassified
noise' for FIX-like files, and 'Movement' for
ICA-AROMA-like files.
:arg excludeLabel: If the file contains a single line containing
component indices, this label will be used for
the components that are not in the list.
Defaults to ``'Signal'`` for FIX-like files, and
``'Unknown'`` for ICA-AROMA-like files.
:arg excludeLabel: If the file contains a single line containing component
indices, this label will be used for the components
that are not in the list. Defaults to 'Signal' for
FIX-like files, and 'Unknown' for ICA-AROMA-like files.
:arg returnIndices: Defaults to ``False``. If ``True``, a list
containing the noisy component numbers that were
listed in the file is returned.
:arg returnIndices: Defaults to ``False``. If ``True``, a list containing
the noisy component numbers that were listed in the
file is returned.
:arg missingLabel: Label to use for any components which are not
present (only used for label files, not for noise
component files).
:arg returnProbabilities: Defaults to ``False``. If ``True``, a list
containing the component classification
probabilities is returned. If the file does not
contain probabilities, every value in this list
will be nan.
:returns: A tuple containing:
......@@ -102,72 +119,55 @@ def loadLabelFile(filename,
- If ``returnIndices is True``, a list of the noisy component
indices (starting from 1) that were specified in the file.
- If ``returnProbabilities is True``, a list of the component
classification probabilities that were specified in the
file (all nan if they are not in the file).
.. note:: Some label files generated by old versions of FIX/Melview do
not contain a line for every component (unknown/unlabelled
components may not be listed). For these files, and also for
files which only contain a component list, there is no way of
knowing how many components were in the data, so the returned
list may contain fewer entries than there are components.
"""
signalLabels = None
filename = op.abspath(filename)
filename = op.abspath(filename)
probabilities = None
signalLabels = None
with open(filename, 'rt') as f:
lines = f.readlines()
if len(lines) < 1:
raise InvalidLabelFileError('Invalid FIX classification '
raise InvalidLabelFileError(f'{filename}: Invalid FIX classification '
'file - not enough lines')
lines = [l.strip() for l in lines]
lines = [l for l in lines if l != '']
# If the file contains a single
# line, we assume that it is just
# a comma-separated list of noise
# components.
if len(lines) == 1:
line = lines[0]
# if the list is contained in
# square brackets, we assume
# that it is a FIX output file,
# where included components have
# been classified as noise, and
# excluded components as signal.
#
# Otherwise we assume that it
# is an AROMA file, where
# included components have
# been classified as being due
# to motion, and excluded
# components unclassified.
if includeLabel is None:
if line[0] == '[': includeLabel = 'Unclassified noise'
else: includeLabel = 'Movement'
if excludeLabel is None:
if line[0] == '[': excludeLabel = 'Signal'
else: excludeLabel = 'Unknown'
else:
signalLabels = [excludeLabel]
# Remove any leading/trailing
# whitespace or brackets.
line = lines[0].strip(' []')
melDir = None
noisyComps = [int(i) for i in line.split(',')]
allLabels = []
for i in range(max(noisyComps)):
if (i + 1) in noisyComps: allLabels.append([includeLabel])
else: allLabels.append([excludeLabel])
# Otherwise, we assume that
# it is a full label file.
else:
# If the file contains one or two lines, we
# assume that it is just a comma-separated list
# of noise components (possibly preceeded by
# the MELODIC directory path)
if len(lines) <= 2:
melDir, noisyComps, allLabels, signalLabels = \
_parseSingleLineLabelFile(lines, includeLabel, excludeLabel)
probabilities = [math.nan] * len(allLabels)
melDir = lines[0]
noisyComps = lines[-1].strip(' []').split(',')
noisyComps = [c for c in noisyComps if c != '']
noisyComps = [int(c) for c in noisyComps]
# Otherwise, we assume that it is a full label file.
else:
melDir, noisyComps, allLabels, probabilities = \
_parseFullLabelFile(filename, lines, missingLabel)
# There's no way to validate
# the melodic directory path,
# but let's try anyway.
if melDir is not None:
if len(melDir.split(',')) >= 3:
raise InvalidLabelFileError(
f'{filename}: First line does not look like '
f'a MELODIC directory path: {melDir}')
# The melodic directory path should
# either be an absolute path, or
......@@ -176,38 +176,6 @@ def loadLabelFile(filename,
if not op.isabs(melDir):
melDir = op.join(op.dirname(filename), melDir)
# Parse the labels for every component
allLabels = []
for i, compLine in enumerate(lines[1:-1]):
tokens = compLine.split(',')
tokens = [t.strip() for t in tokens]
if len(tokens) < 3:
raise InvalidLabelFileError(
'Invalid FIX classification file - '
'line {}: {}'.format(i + 1, compLine))
try:
compIdx = int(tokens[0])
except ValueError:
raise InvalidLabelFileError(
'Invalid FIX classification file - '
'line {}: {}'.format(i + 1, compLine))
if tokens[-1].lower() in ('true', 'false'):
compLabels = tokens[1:-1]
else:
compLabels = tokens[1:]
if compIdx != i + 1:
raise InvalidLabelFileError(
'Invalid FIX classification file - wrong component '
'number at line {}: {}'.format(i + 1, compLine))
allLabels.append(compLabels)
# Validate the labels against
# the noisy list - all components
# in the noisy list should not
......@@ -218,8 +186,8 @@ def loadLabelFile(filename,
noise = isNoisyComponent(labels, signalLabels)
if noise and (comp not in noisyComps):
raise InvalidLabelFileError('Noisy component {} has invalid '
'labels: {}'.format(comp, labels))
raise InvalidLabelFileError(f'{filename}: Noisy component {comp} '
f'has invalid labels: {labels}')
for comp in noisyComps:
......@@ -228,44 +196,187 @@ def loadLabelFile(filename,
noise = isNoisyComponent(labels, signalLabels)
if not noise:
raise InvalidLabelFileError('Noisy component {} is missing '
'a noise label'.format(comp))
raise InvalidLabelFileError(f'{filename}: Noisy component {comp} '
'is missing a noise label')
retval = [melDir, allLabels]
if returnIndices: return melDir, allLabels, noisyComps
else: return melDir, allLabels
if returnIndices: retval.append(noisyComps)
if returnProbabilities: retval.append(probabilities)
return tuple(retval)
def _parseSingleLineLabelFile(lines, includeLabel, excludeLabel):
"""Called by :func:`loadLabelFile`. Parses the contents of an
ICA-AROMA-style label file which just contains a list of noise
components (and possibly the MELODIC directory path), e.g.::
filtered_func_data.ica
[2, 5, 6, 7]
"""
signalLabels = None
noisyComps = lines[-1]
if len(lines) == 2: melDir = lines[0]
else: melDir = None
# if the list is contained in
# square brackets, we assume
# that it is a FIX output file,
# where included components have
# been classified as noise, and
# excluded components as signal.
#
# Otherwise we assume that it
# is an AROMA file, where
# included components have
# been classified as being due
# to motion, and excluded
# components unclassified.
if includeLabel is None:
if noisyComps[0] == '[': includeLabel = 'Unclassified noise'
else: includeLabel = 'Movement'
if excludeLabel is None:
if noisyComps[0] == '[': excludeLabel = 'Signal'
else: excludeLabel = 'Unknown'
else:
signalLabels = [excludeLabel]
# Remove any leading/trailing
# whitespace or brackets.
noisyComps = noisyComps.strip(' []')
noisyComps = [int(i) for i in noisyComps.split(',')]
allLabels = []
for i in range(max(noisyComps)):
if (i + 1) in noisyComps: allLabels.append([includeLabel])
else: allLabels.append([excludeLabel])
return melDir, noisyComps, allLabels, signalLabels
def _parseFullLabelFile(filename, lines, missingLabel):
"""Called by :func:`loadLabelFile`. Parses the contents of a
FIX/Melview-style label file which contains labels for each component,
e.g.:
filtered_func_data.ica
1, Signal, False
2, Unclassified Noise, True
3, Unknown, False
4, Signal, False
5, Unclassified Noise, True
6, Unclassified Noise, True
7, Unclassified Noise, True
8, Signal, False
[2, 5, 6, 7]
"""
melDir = lines[0]
noisyComps = lines[-1].strip(' []').split(',')
noisyComps = [c for c in noisyComps if c != '']
noisyComps = [int(c) for c in noisyComps]
# Parse the labels for every component.
# Initially store as a {comp : ([labels], probability)} dict.
allLabels = {}
for i, compLine in enumerate(lines[1:-1]):
tokens = compLine.split(',')
tokens = [t.strip() for t in tokens]
if len(tokens) < 3:
raise InvalidLabelFileError(
f'{filename}: Invalid FIX classification '
f'file - line: {i + 1}: {compLine}')
try:
compIdx = int(tokens[0])
if compIdx in allLabels:
raise ValueError()
except ValueError:
raise InvalidLabelFileError(
f'{filename}: Invalid FIX classification '
f'file - line {i + 1}: {compLine}')
tokens = tokens[1:]
probability = math.nan
# last token could be classification probability
if _isfloat(tokens[-1]):
probability = float(tokens[-1])
tokens = tokens[:-1]
# true/false is ignored as it is superfluous
if tokens[-1].lower() in ('true', 'false'):
tokens = tokens[:-1]
allLabels[compIdx] = tokens, probability
# Convert {comp : [labels]} into a list
# of lists, filling in missing components
allLabelsList = []
probabilities = []
for i in range(max(it.chain(allLabels.keys(), noisyComps))):
labels, prob = allLabels.get(i + 1, ([missingLabel], math.nan))
allLabelsList.append(labels)
probabilities.append(prob)
allLabels = allLabelsList
return melDir, noisyComps, allLabels, probabilities
def _isfloat(s):
"""Returns True if the given string appears to contain a floating
point number, False otherwise.
"""
try:
float(s)
return True
except Exception:
return False
def saveLabelFile(allLabels,
filename,
dirname=None,
listBad=True,
signalLabels=None):
signalLabels=None,
probabilities=None):
"""Saves the given classification labels to the specified file. The
classifications are saved in the format described in the
:func:`loadLabelFile` method.
:arg allLabels: A list of lists, one list for each component, where
each list contains the labels for the corresponding
component.
:arg allLabels: A list of lists, one list for each component, where
each list contains the labels for the corresponding
component.
:arg filename: Name of the file to which the labels should be saved.
:arg filename: Name of the file to which the labels should be saved.
:arg dirname: If provided, is output as the first line of the file.
Intended to be a relative path to the MELODIC analysis
directory with which this label file is associated. If
not provided, a ``'.'`` is output as the first line.
:arg dirname: If provided, is output as the first line of the file.
Intended to be a relative path to the MELODIC analysis
directory with which this label file is associated. If
not provided, a ``'.'`` is output as the first line.
:arg listBad: If ``True`` (the default), the last line of the file
will contain a comma separated list of components which
are deemed 'noisy' (see :func:`isNoisyComponent`).
:arg listBad: If ``True`` (the default), the last line of the file
will contain a comma separated list of components which
are deemed 'noisy' (see :func:`isNoisyComponent`).
:arg signalLabels: Labels which should be deemed 'signal' - see the
:func:`isNoisyComponent` function.
:arg signalLabels: Labels which should be deemed 'signal' - see the
:func:`isNoisyComponent` function.
:arg probabilities: Classification probabilities. If provided, the
probability for each component is saved to the file.
"""
lines = []
noisyComps = []
if probabilities is not None and len(probabilities) != len(allLabels):
raise ValueError('len(probabilities) != len(allLabels)')
# The first line - the melodic directory name
if dirname is None:
dirname = '.'
......@@ -283,6 +394,9 @@ def saveLabelFile(allLabels,
labels = [l.replace(',', '_') for l in labels]
tokens = [str(comp)] + labels + [str(noise)]
if probabilities is not None:
tokens.append(f'{probabilities[i]:0.6f}')
lines.append(', '.join(tokens))
if noise:
......@@ -318,4 +432,3 @@ class InvalidLabelFileError(Exception):
"""Exception raised by the :func:`loadLabelFile` function when an attempt
is made to load an invalid label file.
"""
pass
......@@ -67,7 +67,8 @@ CORE_GEOMETRY_FILES = ['?h.orig',
'?h.pial',
'?h.white',
'?h.inflated',
'?h.sphere']
'?h.sphere',
'?h.pial_semi_inflated']
"""File patterns for identifying the core Freesurfer geometry files. """
......@@ -76,7 +77,8 @@ CORE_GEOMETRY_DESCRIPTIONS = [
"Freesurfer surface (pial)",
"Freesurfer surface (white matter)",
"Freesurfer surface (inflated)",
"Freesurfer surface (sphere)"]
"Freesurfer surface (sphere)",
"Freesurfer surface (pial semi-inflated)"]
"""A description for each extension in :attr:`GEOMETRY_EXTENSIONS`. """
......
......@@ -25,12 +25,14 @@ are available:
import glob
import re
import os.path as op
import numpy as np
import nibabel as nib
import fsl.utils.path as fslpath
import fsl.utils.bids as bids
import fsl.data.constants as constants
import fsl.data.mesh as fslmesh
......@@ -45,6 +47,15 @@ EXTENSION_DESCRIPTIONS = ['GIFTI surface file', 'GIFTI file']
"""A description for each of the :data:`ALLOWED_EXTENSIONS`. """
VERTEX_DATA_EXTENSIONS = ['.func.gii',
'.shape.gii',
'.label.gii',
'.time.gii']
"""File suffixes which are interpreted as GIFTI vertex data files,
containing data values for every vertex in the mesh.
"""
class GiftiMesh(fslmesh.Mesh):
"""Class which represents a GIFTI surface image. This is essentially
just a 3D model made of triangles.
......@@ -90,9 +101,24 @@ class GiftiMesh(fslmesh.Mesh):
for i, v in enumerate(vertices):
if i == 0: key = infile
else: key = '{}_{}'.format(infile, i)
else: key = f'{infile}_{i}'
self.addVertices(v, key, select=(i == 0), fixWinding=fixWinding)
self.setMeta(infile, surfimg)
self.meta[infile] = surfimg
# Copy all metadata entries for the GIFTI image
for k, v in surfimg.meta.items():
self.meta[k] = v
# and also for each GIFTI data array - triangles
# are stored under "faces", and pointsets are
# stored under "vertices"/[0,1,2...] (as there may
# be multiple pointsets in a file)
self.meta['vertices'] = {}
for i, arr in enumerate(surfimg.darrays):
if arr.intent == constants.NIFTI_INTENT_POINTSET:
self.meta['vertices'][i] = dict(arr.meta)
elif arr.intent == constants.NIFTI_INTENT_TRIANGLE:
self.meta['faces'] = dict(arr.meta)
if vdata is not None:
self.addVertexData(infile, vdata)
......@@ -102,8 +128,11 @@ class GiftiMesh(fslmesh.Mesh):
# as the specfiied one.
if loadAll:
# Only attempt to auto-load sensibly
# named gifti files (i.e. *.surf.gii,
# rather than *.gii).
surfFiles = relatedFiles(infile, [ALLOWED_EXTENSIONS[0]])
nvertices = vertices[0].shape[0]
surfFiles = relatedFiles(infile, ALLOWED_EXTENSIONS)
for sfile in surfFiles:
......@@ -116,7 +145,7 @@ class GiftiMesh(fslmesh.Mesh):
continue
self.addVertices(vertices[0], sfile, select=False)
self.setMeta(sfile, surfimg)
self.meta[sfile] = surfimg
def loadVertices(self, infile, key=None, *args, **kwargs):
......@@ -140,10 +169,10 @@ class GiftiMesh(fslmesh.Mesh):
for i, v in enumerate(vertices):
if i == 0: key = infile
else: key = '{}_{}'.format(infile, i)
else: key = f'{infile}_{i}'
vertices[i] = self.addVertices(v, key, *args, **kwargs)
self.setMeta(infile, surfimg)
self.meta[infile] = surfimg
return vertices
......@@ -207,15 +236,15 @@ def loadGiftiMesh(filename):
vdata = [d for d in gimg.darrays if d.intent not in (pscode, tricode)]
if len(triangles) != 1:
raise ValueError('{}: GIFTI surface files must contain '
'exactly one triangle array'.format(filename))
raise ValueError(f'{filename}: GIFTI surface files must '
'contain exactly one triangle array')
if len(pointsets) == 0:
raise ValueError('{}: GIFTI surface files must contain '
'at least one pointset array'.format(filename))
raise ValueError(f'{filename}: GIFTI surface files must '
'contain at least one pointset array')
vertices = [ps.data for ps in pointsets]
indices = triangles[0].data
indices = np.atleast_2d(triangles[0].data)
if len(vdata) == 0: vdata = None
else: vdata = prepareGiftiVertexData(vdata, filename)
......@@ -259,17 +288,17 @@ def prepareGiftiVertexData(darrays, filename=None):
vertices.
"""
intents = set([d.intent for d in darrays])
intents = {d.intent for d in darrays}
if len(intents) != 1:
raise ValueError('{} contains multiple (or no) intents'
': {}'.format(filename, intents))
raise ValueError(f'{filename} contains multiple '
f'(or no) intents: {intents}')
intent = intents.pop()
if intent in (constants.NIFTI_INTENT_POINTSET,
constants.NIFTI_INTENT_TRIANGLE):
raise ValueError('{} contains surface data'.format(filename))
raise ValueError(f'{filename} contains surface data')
# Just a single array - return it as-is.
# n.b. Storing (M, N) data in a single
......@@ -284,8 +313,8 @@ def prepareGiftiVertexData(darrays, filename=None):
vdata = [d.data for d in darrays]
if any([len(d.shape) != 1 for d in vdata]):
raise ValueError('{} contains one or more non-vector '
'darrays'.format(filename))
raise ValueError(f'{filename} contains one or '
'more non-vector darrays')
vdata = np.vstack(vdata).T
vdata = vdata.reshape(vdata.shape[0], -1)
......@@ -298,45 +327,110 @@ def relatedFiles(fname, ftypes=None):
directory which appear to be related with the given one. Files which
share the same prefix are assumed to be related to the given file.
This function assumes that the GIFTI files are named according to a
standard convention - the following conventions are supported:
- HCP-style, i.e.: ``<subject>.<hemi>.<type>.<space>.<ftype>.gii``
- BIDS-style, i.e.:
``<source_prefix>_hemi-<hemi>[_space-<space>]*_<suffix>.<ftype>.gii``
If the files are not named according to one of these conventions, this
function will return an empty list.
:arg fname: Name of the file to search for related files for
:arg ftype: If provided, only files with suffixes in this list are
searched for. Defaults to files which contain vertex data.
searched for. Defaults to :attr:`VERTEX_DATA_EXTENSIONS`.
"""
if ftypes is None:
ftypes = ['.func.gii', '.shape.gii', '.label.gii', '.time.gii']
ftypes = VERTEX_DATA_EXTENSIONS
# We want to return all files in the same
# directory which have the following name:
path = op.abspath(fname)
dirname, fname = op.split(path)
# We want to identify all files in the same
# directory which are associated with the
# given file. We assume that the files are
# named according to one of the following
# conventions:
#
# [subj].[hemi].[type].*.[ftype]
# - HCP style:
# <subject>.<hemi>.<type>.<space>.<ftype>.gii
#
# where
# - [subj] is the subject ID, and matches fname
# - BIDS style:
# <source_prefix>_hemi-<hemi>[_space-<space>]*.<ftype>.gii
#
# - [hemi] is the hemisphere, and matches fname
# We cannot assume consistent ordering of
# the entities (key-value pairs) within a
# BIDS style filename, so we cannot simply
# use a regular expression or glob pattern.
# Instead, for each style we define:
#
# - [type] defines the file contents
# - a "matcher" function, which tests
# whether the file matches the style,
# and returns the important elements
# from the file name.
#
# - suffix is func, shape, label, time, or `ftype`
path = op.abspath(fname)
dirname, fname = op.split(path)
# get the [subj].[hemi] prefix
try:
subj, hemi, _ = fname.split('.', 2)
prefix = '.'.join((subj, hemi))
except Exception:
# - a "searcher" function, which takes
# the elements of the input file
# that were extracted by the matcher,
# and searches for other related files
# HCP style - extract "<subject>.<hemi>"
# and "<space>".
def matchhcp(f):
pat = r'^(.*\.[LR])\..*\.(.*)\..*\.gii$'
match = re.match(pat, f)
if match:
return match.groups()
else:
return None
def searchhcp(match, ftype):
prefix, space = match
template = f'{prefix}.*.{space}{ftype}'
return glob.glob(op.join(dirname, template))
# BIDS style - extract all entities (kv
# pairs), ignoring specific irrelevant
# ones.
def matchbids(f):
try: match = bids.BIDSFile(f)
except ValueError: return None
match.entities.pop('desc', None)
return match
def searchbids(match, ftype):
allfiles = glob.glob(op.join(dirname, '*{}'.format(ftype)))
for f in allfiles:
try: bf = bids.BIDSFile(f)
except ValueError: continue
if bf.match(match, False):
yield f
# find the first style that matches
matchers = [matchhcp, matchbids]
searchers = [searchhcp, searchbids]
for matcher, searcher in zip(matchers, searchers):
match = matcher(fname)
if match:
break
# Give up if the file does
# not match any known style.
else:
return []
# Build a list of files in the same
# directory and matching the template
related = []
for ftype in ftypes:
hits = glob.glob(op.join(dirname, '{}*{}'.format(prefix, ftype)))
hits = searcher(match, ftype)
# eliminate dupes
related.extend([h for h in hits if h not in related])
# exclude the file itself
return [r for r in related if r != path]
This diff is collapsed.
......@@ -7,6 +7,9 @@
"""This module provides the :class:`ImageWrapper` class, which can be used
to manage data access to ``nibabel`` NIFTI images.
.. note:: This module is deprecated - it is being moved to the FSLeyes project,
and will be removed in a future version of fslpy.
Terminology
-----------
......@@ -42,9 +45,10 @@ import collections
import collections.abc as abc
import itertools as it
import numpy as np
import nibabel as nib
import numpy as np
import fsl.data.image as fslimage
import fsl.utils.deprecated as deprecated
import fsl.utils.notifier as notifier
import fsl.utils.naninfrange as nir
import fsl.utils.idle as idle
......@@ -148,6 +152,8 @@ class ImageWrapper(notifier.Notifier):
"""
@deprecated.deprecated('3.9.0', '4.0.0',
'The ImageWrapper has been migrated to FSLeyes')
def __init__(self,
image,
name=None,
......@@ -175,8 +181,6 @@ class ImageWrapper(notifier.Notifier):
data range is updated directly on reads/writes.
"""
import fsl.data.image as fslimage
self.__image = image
self.__name = name
self.__taskThread = None
......@@ -388,6 +392,14 @@ class ImageWrapper(notifier.Notifier):
self.__data = np.asanyarray(self.__image.dataobj)
@property
def dataIsLoaded(self):
"""Return true if the image data has been loaded into memory, ``False``
otherwise.
"""
return self.__data is not None
def __getData(self, sliceobj, isTuple=False):
"""Retrieves the image data at the location specified by ``sliceobj``.
......@@ -717,104 +729,29 @@ class ImageWrapper(notifier.Notifier):
self.__updateDataRangeOnWrite(slices, values)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to fsl.data.image')
def isValidFancySliceObj(sliceobj, shape):
"""Returns ``True`` if the given ``sliceobj`` is a valid and fancy slice
object.
``nibabel`` refers to slice objects as "fancy" if they comprise anything
but tuples of integers and simple ``slice`` objects. The ``ImageWrapper``
class supports one type of "fancy" slicing, where the ``sliceobj`` is a
boolean ``numpy`` array of the same shape as the image.
This function returns ``True`` if the given ``sliceobj`` adheres to these
requirements, ``False`` otherwise.
"""
# We only support boolean numpy arrays
# which have the same shape as the image
return (isinstance(sliceobj, np.ndarray) and
sliceobj.dtype == np.bool and
np.prod(sliceobj.shape) == np.prod(shape))
"""Deprecated - moved to :mod:`fsl.data.image`."""
return fslimage.isValidFancySliceObj(sliceobj, shape)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to fsl.data.image')
def canonicalSliceObj(sliceobj, shape):
"""Returns a canonical version of the given ``sliceobj``. See the
``nibabel.fileslice.canonical_slicers`` function.
"""
# Fancy slice objects must have
# the same shape as the data
if isValidFancySliceObj(sliceobj, shape):
return sliceobj.reshape(shape)
else:
if not isinstance(sliceobj, tuple):
sliceobj = (sliceobj,)
if len(sliceobj) > len(shape):
sliceobj = sliceobj[:len(shape)]
return nib.fileslice.canonical_slicers(sliceobj, shape)
"""Deprecated - moved to :mod:`fsl.data.image`."""
return fslimage.canonicalSliceObj(sliceobj, shape)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to fsl.data.image')
def expectedShape(sliceobj, shape):
"""Given a slice object, and the shape of an array to which
that slice object is going to be applied, returns the expected
shape of the result.
.. note:: It is assumed that the ``sliceobj`` has been passed through
the :func:`canonicalSliceObj` function.
:arg sliceobj: Something which can be used to slice an array
of shape ``shape``.
:arg shape: Shape of the array being sliced.
:returns: A tuple containing:
- Expected number of dimensions of the result
- Expected shape of the result (or ``None`` if
``sliceobj`` is fancy).
"""
if isValidFancySliceObj(sliceobj, shape):
return 1, None
# Truncate some dimensions from the
# slice object if it has too many
# (e.g. trailing dims of length 1).
elif len(sliceobj) > len(shape):
sliceobj = sliceobj[:len(shape)]
# Figure out the number of dimensions
# that the result should have, given
# this slice object.
expShape = []
for i in range(len(sliceobj)):
# Each dimension which has an
# int slice will be collapsed
if isinstance(sliceobj[i], int):
continue
start = sliceobj[i].start
stop = sliceobj[i].stop
if start is None: start = 0
if stop is None: stop = shape[i]
stop = min(stop, shape[i])
expShape.append(stop - start)
return len(expShape), expShape
"""Deprecated - moved to :mod:`fsl.data.image`."""
return fslimage.expectedShape(sliceobj, shape)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def sliceObjToSliceTuple(sliceobj, shape):
"""Turns an array slice object into a tuple of (low, high) index
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Turns an array slice object into a tuple of (low, high) index
pairs, one pair for each dimension in the given shape
:arg sliceobj: Something which can be used to slice an array of shape
......@@ -853,8 +790,11 @@ def sliceObjToSliceTuple(sliceobj, shape):
return tuple(indices)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def sliceTupleToSliceObj(slices):
"""Turns a sequence of (low, high) index pairs into a tuple of array
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Turns a sequence of (low, high) index pairs into a tuple of array
``slice`` objects.
:arg slices: A sequence of (low, high) index pairs.
......@@ -868,8 +808,11 @@ def sliceTupleToSliceObj(slices):
return tuple(sliceobj)
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def adjustCoverage(oldCoverage, slices):
"""Adjusts/expands the given ``oldCoverage`` so that it covers the
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Adjusts/expands the given ``oldCoverage`` so that it covers the
given set of ``slices``.
:arg oldCoverage: A ``numpy`` array of shape ``(2, n)`` containing
......@@ -916,8 +859,11 @@ return code for the :func:`sliceOverlap` function.
"""
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def sliceOverlap(slices, coverage):
"""Determines whether the given ``slices`` overlap with the given
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Determines whether the given ``slices`` overlap with the given
``coverage``.
:arg slices: A sequence of (low, high) index pairs, assumed to cover
......@@ -983,8 +929,11 @@ def sliceOverlap(slices, coverage):
elif np.all(overlapStates == OVERLAP_ALL): return OVERLAP_ALL
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def sliceCovered(slices, coverage):
"""Returns ``True`` if the portion of the image data calculated by
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Returns ``True`` if the portion of the image data calculated by
the given ``slices` has already been calculated, ``False`` otherwise.
:arg slices: A sequence of (low, high) index pairs, assumed to cover
......@@ -1015,8 +964,11 @@ def sliceCovered(slices, coverage):
return True
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def calcExpansion(slices, coverage):
"""Calculates a series of *expansion* slices, which can be used to expand
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Calculates a series of *expansion* slices, which can be used to expand
the given ``coverage`` so that it includes the given ``slices``.
:arg slices: Slices that the coverage needs to be expanded to cover.
......@@ -1185,8 +1137,11 @@ def calcExpansion(slices, coverage):
return volumes, expansions
@deprecated.deprecated('3.9.0', '4.0.0', 'Moved to FSLeyes')
def collapseExpansions(expansions, numDims):
"""Scans through the given list of expansions (each assumed to pertain
"""Deprecated - the imagewrapper has been moved to FSLeyes.
Scans through the given list of expansions (each assumed to pertain
to a single 3D image), and combines any which cover the same
image area, and cover adjacent volumes.
......
......@@ -33,7 +33,6 @@ import logging
import os.path as op
import numpy as np
import fsl.utils.path as fslpath
import fsl.data.image as fslimage
import fsl.data.featanalysis as featanalysis
......@@ -63,10 +62,9 @@ def isMelodicImage(path):
def isMelodicDir(path):
"""Returns ``True`` if the given path looks like it is contained within
a MELODIC directory, ``False`` otherwise. A melodic directory:
"""Returns ``True`` if the given path looks like it is a MELODIC directory,
``False`` otherwise. A MELODIC directory:
- Must be named ``*.ica``.
- Must contain a file called ``melodic_IC.nii.gz`` or
``melodic_oIC.nii.gz``.
- Must contain a file called ``melodic_mix``.
......@@ -75,12 +73,7 @@ def isMelodicDir(path):
path = op.abspath(path)
if op.isdir(path): dirname = path
else: dirname = op.dirname(path)
sufs = ['.ica']
if not any([dirname.endswith(suf) for suf in sufs]):
if not op.isdir(path):
return False
# Must contain an image file called
......@@ -88,7 +81,7 @@ def isMelodicDir(path):
prefixes = ['melodic_IC', 'melodic_oIC']
for p in prefixes:
try:
fslimage.addExt(op.join(dirname, p))
fslimage.addExt(op.join(path, p))
break
except fslimage.PathError:
pass
......@@ -97,8 +90,8 @@ def isMelodicDir(path):
# Must contain files called
# melodic_mix and melodic_FTmix
if not op.exists(op.join(dirname, 'melodic_mix')): return False
if not op.exists(op.join(dirname, 'melodic_FTmix')): return False
if not op.exists(op.join(path, 'melodic_mix')): return False
if not op.exists(op.join(path, 'melodic_FTmix')): return False
return True
......@@ -108,10 +101,13 @@ def getAnalysisDir(path):
to that MELODIC directory is returned. Otherwise, ``None`` is returned.
"""
meldir = fslpath.deepest(path, ['.ica'])
if not op.isdir(path):
path = op.dirname(path)
if meldir is not None and isMelodicDir(meldir):
return meldir
while path not in (op.sep, ''):
if isMelodicDir(path):
return path
path = op.dirname(path)
return None
......@@ -137,10 +133,18 @@ def getDataFile(meldir):
if topDir is None:
return None
dataFile = op.join(topDir, 'filtered_func_data')
# People often rename filtered_func_data.nii.gz
# to something like filtered_func_data_clean.nii.gz,
# because that is the recommended approach when
# performing ICA-based denoising). So we try both.
candidates = ['filtered_func_data', 'filtered_func_data_clean']
try: return fslimage.addExt(dataFile)
except fslimage.PathError: return None
for candidate in candidates:
dataFile = op.join(topDir, candidate)
try: return fslimage.addExt(dataFile)
except fslimage.PathError: continue
return None
def getMeanFile(meldir):
......@@ -187,7 +191,7 @@ def getNumComponents(meldir):
contained in the given directrory.
"""
icImg = fslimage.Image(getICFile(meldir), loadData=False, calcRange=False)
icImg = fslimage.Image(getICFile(meldir))
return icImg.shape[3]
......
......@@ -74,9 +74,7 @@ class MelodicImage(fslimage.Image):
dataFile = self.getDataFile()
if dataFile is not None:
dataImage = fslimage.Image(dataFile,
loadData=False,
calcRange=False)
dataImage = fslimage.Image(dataFile)
if dataImage.ndim >= 4:
self.__tr = dataImage.pixdim[3]
......
This diff is collapsed.