fsl/data/melodicimage.py

@@ -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]
 

fsl/data/mesh.py

@@ -41,6 +41,13 @@ import fsl.transform.affine as affine
 log = logging.getLogger(__name__)
 
 
+class IncompatibleVerticesError(ValueError):
+    """``ValueError`` raised by the :meth:`Mesh.addVertices` method if
+    an attempt is made to add a vertex set with the wrong number of
+    vertices.
+    """
+
+
 class Mesh(notifier.Notifier, meta.Meta):
     """The ``Mesh`` class represents a 3D model. A mesh is defined by a
     collection of ``N`` vertices, and ``M`` triangles.  The triangles are
@@ -154,16 +161,8 @@ class Mesh(notifier.Notifier, meta.Meta):
     """
 
 
-    def __new__(cls, *args, **kwargs):
-        """Create a ``Mesh``. We must override ``__new__``, otherwise the
-        :class:`Meta` and :class:`Notifier` ``__new__`` methods will not be
-        called correctly.
-        """
-        return super(Mesh, cls).__new__(cls, *args, **kwargs)
-
-
     def __init__(self,
-                 indices,
+                 indices=None,
                  name='mesh',
                  dataSource=None,
                  vertices=None,
@@ -174,7 +173,8 @@ class Mesh(notifier.Notifier, meta.Meta):
         :meth:`addVertices` method.
 
         :arg indices:    A list of indices into the vertex data, defining the
-                         mesh triangles.
+                         mesh triangles. If not provided, must be provided
+                         after creation via the :meth:`indices` setter method.
 
         :arg name:       A name for this ``Mesh``.
 
@@ -187,22 +187,31 @@ class Mesh(notifier.Notifier, meta.Meta):
                          :meth:`addVertices` method along with ``vertices``.
         """
 
+        if indices is None and vertices is not None:
+            raise ValueError('Indices must be provided '
+                             'if vertices are provided')
+
         self.__name       = name
         self.__dataSource = dataSource
-        self.__nvertices  = indices.max() + 1
-        self.__selected   = None
+
+        # nvertices/indices are assigned in the
+        # indices setter method.
 
         # We potentially store two copies of
-        # the indices, with opposite unwinding
-        # orders. The vindices dict stores refs
-        # to one or the other for each vertex
-        # set.
-        self.__indices      = np.asarray(indices).reshape((-1, 3))
+        # the indices, - one set (__indices)
+        # as provided, and the other
+        # (__fixedIndices) with opposite
+        # unwinding orders. The vindices dict
+        # stores refs to one or the other for
+        # each vertex set.
+        self.__nvertices    = None
+        self.__indices      = None
         self.__fixedIndices = None
         self.__vindices     = collections.OrderedDict()
 
         # All of these are populated
         # in the addVertices method
+        self.__selected = None
         self.__vertices = collections.OrderedDict()
         self.__loBounds = collections.OrderedDict()
         self.__hiBounds = collections.OrderedDict()
@@ -220,8 +229,9 @@ class Mesh(notifier.Notifier, meta.Meta):
         # in the trimesh method
         self.__trimesh = collections.OrderedDict()
 
-        # Add initial vertex
-        # set if provided
+        # Add initial indices/vertices if provided
+        if indices is not None:
+            self.indices = indices
         if vertices is not None:
             self.addVertices(vertices, fixWinding=fixWinding)
 
@@ -290,10 +300,25 @@ class Mesh(notifier.Notifier, meta.Meta):
 
     @property
     def indices(self):
-        """The ``(M, 3)`` triangles of this mesh. """
+        """The ``(M, 3)`` triangles of this mesh. Returns ``None`` if
+        indices have not yet been assigned.
+        """
+        if self.__indices is None:
+            return None
         return self.__vindices[self.__selected]
 
 
+    @indices.setter
+    def indices(self, indices):
+        """Set the indices for this mesh. """
+        if self.__indices is not None:
+            raise ValueError('Indices are already set')
+
+        indices          = np.asarray(indices, dtype=np.int32)
+        self.__nvertices = int(indices.max()) + 1
+        self.__indices   = indices.reshape((-1, 3))
+
+
     @property
     def normals(self):
         """A ``(M, 3)`` array containing surface normals for every
@@ -337,7 +362,13 @@ class Mesh(notifier.Notifier, meta.Meta):
         ``Mesh`` instance. The bounding box is arranged like so:
 
             ``((xlow, ylow, zlow), (xhigh, yhigh, zhigh))``
+
+        Returns ``None`` if indices or vertices have not yet been assigned.
         """
+
+        if self.__indices is None or len(self.__vertices) == 0:
+            return None
+
         lo = self.__loBounds[self.__selected]
         hi = self.__hiBounds[self.__selected]
         return lo, hi
@@ -388,10 +419,13 @@ class Mesh(notifier.Notifier, meta.Meta):
 
         :returns:        The vertices, possibly reshaped
 
-        :raises:         ``ValueError`` if the provided ``vertices`` array
-                         has the wrong number of vertices.
+        :raises:         ``IncompatibleVerticesError`` if the provided
+                         ``vertices`` array has the wrong number of vertices.
         """
 
+        if self.__indices is None:
+            raise ValueError('Mesh indices have not yet been set')
+
         if key is None:
             key = 'default'
 
@@ -407,10 +441,9 @@ class Mesh(notifier.Notifier, meta.Meta):
         # reshape raised an error -
         # wrong number of vertices
         except ValueError:
-            raise ValueError('{}: invalid number of vertices: '
-                             '{} != ({}, 3)'.format(key,
-                                                    vertices.shape,
-                                                    self.nvertices))
+            raise IncompatibleVerticesError(
+                f'{key}: invalid number of vertices: '
+                f'{vertices.shape} != ({self.nvertices}, 3)')
 
         self.__vertices[key] = vertices
         self.__vindices[key] = self.__indices
@@ -581,7 +614,7 @@ class Mesh(notifier.Notifier, meta.Meta):
 
         # sort by ray. I'm Not sure if this is
         # needed - does trimesh do it for us?
-        rayIdxs = np.asarray(np.argsort(rays), np.int)
+        rayIdxs = np.asarray(np.argsort(rays))
         locs    = locs[rayIdxs]
         tris    = tris[rayIdxs]
 
@@ -695,7 +728,7 @@ def calcFaceNormals(vertices, indices):
     fnormals = np.cross((v1 - v0), (v2 - v0))
     fnormals = affine.normalise(fnormals)
 
-    return fnormals
+    return np.atleast_2d(fnormals)
 
 
 def calcVertexNormals(vertices, indices, fnormals):
@@ -709,7 +742,7 @@ def calcVertexNormals(vertices, indices, fnormals):
                    the mesh.
     """
 
-    vnormals = np.zeros((vertices.shape[0], 3), dtype=np.float)
+    vnormals = np.zeros((vertices.shape[0], 3), dtype=float)
 
     # TODO make fast. I can't figure
     # out how to use np.add.at to
@@ -758,15 +791,19 @@ def needsFixing(vertices, indices, fnormals, loBounds, hiBounds):
     ivert = np.argmin(dists)
     vert  = vertices[ivert]
 
-    # Pick a triangle that
-    # this vertex is in and
-    # ges its face normal
-    itri = np.where(indices == ivert)[0][0]
-    n    = fnormals[itri, :]
+    # Get all the triangles
+    # that this vertex is in
+    # and their face normals
+    itris = np.where(indices == ivert)[0]
+    norms = fnormals[itris, :]
 
-    # Make sure the angle between the
+    # Calculate the angle between each
     # normal, and a vector from the
-    # vertex to the camera is positive
-    # If it isn't, we need to flip the
-    # triangle winding order.
-    return np.dot(n, affine.normalise(camera - vert)) < 0
+    # vertex to the camera. If more than
+    # 50% of the angles are negative
+    # (== more than 90 degrees == the
+    # face is facing away from the
+    # camera), assume that we need to
+    # flip the triangle winding order.
+    angles = np.dot(norms, affine.normalise(camera - vert))
+    return ((angles >= 0).sum() / len(itris)) < 0.5

fsl/data/mghimage.py

@@ -10,8 +10,9 @@ Freesurfer ``mgh``/``mgz`` image files.
 
 
 import os.path as op
+import            pathlib
 
-import            six
+import numpy   as np
 import nibabel as nib
 
 import fsl.utils.path       as fslpath
@@ -37,7 +38,7 @@ class MGHImage(fslimage.Image):
      - http://nipy.org/nibabel/reference/nibabel.freesurfer.html
     """
 
-    def __init__(self, image, *args, **kwargs):
+    def __init__(self, image, **kwargs):
         """Create a ``MGHImage``.
 
         :arg image: Name of MGH file, or a
@@ -46,7 +47,7 @@ class MGHImage(fslimage.Image):
         All other arguments are passed through to :meth:`Image.__init__`
         """
 
-        if isinstance(image, six.string_types):
+        if isinstance(image, (str, pathlib.Path)):
             filename = op.abspath(image)
             name     = op.basename(filename)
             image    = nib.load(image)
@@ -54,15 +55,34 @@ class MGHImage(fslimage.Image):
             name     = 'MGH image'
             filename = None
 
-        data     = image.get_data()
+        data     = np.asanyarray(image.dataobj)
         xform    = image.affine
+        pixdim   = image.header.get_zooms()
         vox2surf = image.header.get_vox2ras_tkr()
 
+        # the image may have an affine which
+        # transforms the data into some space
+        # with a scaling that is different to
+        # the pixdims. So we create a header
+        # object with both the affine and the
+        # pixdims, so they are both preserved.
+        #
+        # Note that we have to set the zooms
+        # after the s/qform, otherwise nibabel
+        # will clobber them with zooms gleaned
+        # fron the affine.
+        header = nib.nifti1.Nifti1Header()
+        header.set_data_shape(data.shape)
+        header.set_sform(xform)
+        header.set_qform(xform)
+        header.set_zooms(pixdim)
+
         fslimage.Image.__init__(self,
                                 data,
-                                xform=xform,
+                                header=header,
                                 name=name,
-                                dataSource=filename)
+                                dataSource=filename,
+                                **kwargs)
 
         if filename is not None:
             self.setMeta('mghImageFile', filename)
@@ -127,3 +147,30 @@ class MGHImage(fslimage.Image):
         coordinates into the surface coordinate system for this image.
         """
         return self.__worldToSurfMat
+
+
+def voxToSurfMat(img):
+    """Generate an affine which can transform the voxel coordinates of
+    the given image into a corresponding Freesurfer surface coordinate
+    system (known as "Torig", or "vox2ras-tkr").
+
+    See https://surfer.nmr.mgh.harvard.edu/fswiki/CoordinateSystems
+
+    :arg img: An :class:`.Image` object.
+
+    :return:  A ``(4, 4)`` matrix encoding an affine transformation from the
+              image voxel coordinate system to the corresponding Freesurfer
+              surface coordinate system.
+    """
+
+    zooms = np.array(img.pixdim[:3])
+    dims  = img.shape[ :3] * zooms / 2
+
+    xform        = np.zeros((4, 4), dtype=np.float32)
+    xform[ 0, 0] = -zooms[0]
+    xform[ 1, 2] =  zooms[2]
+    xform[ 2, 1] = -zooms[1]
+    xform[ 3, 3] = 1
+    xform[:3, 3] = [dims[0], -dims[2], dims[1]]
+
+    return xform

fsl/data/vest.py

@@ -11,9 +11,14 @@
 
    looksLikeVestLutFile
    loadVestLutFile
+   loadVestFile
+   generateVest
 """
 
 
+import textwrap as tw
+import             io
+
 import numpy as np
 
 
@@ -76,3 +81,70 @@ def loadVestLutFile(path, normalise=True):
 
     else:
         return colours
+
+
+def loadVestFile(path, ignoreHeader=True):
+    """Loads numeric data from a VEST file, returning it as a ``numpy`` array.
+
+    :arg ignoreHeader: if ``True`` (the default), the matrix shape specified
+                       in the VEST header information is ignored, and the shape
+                       inferred from the data. Otherwise, if the number of
+                       rows/columns specified in the VEST header information
+                       does not match the matrix shape, a ``ValueError`` is
+                       raised.
+
+    :returns:          a ``numpy`` array containing the matrix data in the
+                       VEST file.
+    """
+
+    data = np.loadtxt(path, comments=['#', '/'])
+
+    if not ignoreHeader:
+        nrows, ncols = None, None
+        with open(path, 'rt') as f:
+            for line in f:
+                if   'NumWaves'  in line: ncols = int(line.split()[1])
+                elif 'NumPoints' in line: nrows = int(line.split()[1])
+                else: continue
+
+                if (ncols is not None) and (nrows is not None):
+                    break
+
+        if tuple(data.shape) != (nrows, ncols):
+            raise ValueError(f'Invalid VEST file ({path}) - data shape '
+                             f'({data.shape}) does not match header '
+                             f'({nrows}, {ncols})')
+
+    return data
+
+
+def generateVest(data):
+    """Generates VEST-formatted text for the given ``numpy`` array.
+
+    :arg data: A 1D or 2D numpy array.
+    :returns:  A string containing a VEST header, and the ``data``.
+    """
+
+    data = np.asanyarray(data)
+
+    if len(data.shape) not in (1, 2):
+        raise ValueError(f'unsupported number of dimensions: {data.shape}')
+
+    data = np.atleast_2d(data)
+
+    if np.issubdtype(data.dtype, np.integer): fmt = '%d'
+    else:                                     fmt = '%0.12f'
+
+    sdata = io.StringIO()
+    np.savetxt(sdata, data, fmt=fmt)
+    sdata = sdata.getvalue()
+
+    nrows, ncols = data.shape
+
+    vest = tw.dedent(f"""
+    /NumWaves {ncols}
+    /NumPoints {nrows}
+    /Matrix
+    """).strip() + '\n' + sdata
+
+    return vest.strip()

fsl/scripts/Text2Vest.py

+#!/usr/bin/env python
+#
+# Text2Vest.py - Convert an ASCII text matrix file into a VEST file.
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""``Text2Vest`` simply takes a plain text ASCII text matrix file, and
+adds a VEST header.
+"""
+
+
+import sys
+
+import numpy as np
+
+import fsl.data.vest as fslvest
+
+
+usage = "Usage: Text2Vest <text_file> <vest_file>"
+
+
+def main(argv=None):
+    """Convert a plain text file to a VEST file. """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) != 2:
+        print(usage)
+        return 0
+
+    infile, outfile = argv
+    data            = np.loadtxt(infile, ndmin=2)
+    vest            = fslvest.generateVest(data)
+
+    with open(outfile, 'wt') as f:
+        f.write(vest)
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/Vest2Text.py

+#!/usr/bin/env python
+#
+# Vest2Text.py - Convert a VEST matrix file into a plain text ASCII file.
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""``Vest2Text`` takes a VEST file containing a 2D matrix, and converts it
+into a plain-text ASCII file.
+"""
+
+
+import sys
+
+import numpy as np
+
+import fsl.data.vest as fslvest
+
+
+usage = "Usage: Vest2Text <vest_file> <text_file>"
+
+
+def main(argv=None):
+    """Convert a VEST file to a plain text file. """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) != 2:
+        print(usage)
+        return 0
+
+    infile, outfile = argv
+    data            = fslvest.loadVestFile(infile)
+
+    if np.issubdtype(data.dtype, np.integer): fmt = '%d'
+    else:                                     fmt = '%0.12f'
+
+    np.savetxt(outfile, data, fmt=fmt)
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/__init__.py

-#!/usr/bin/env python
-#
-# __init__.py - The fsl.scripts package.
-#
-# Author: Paul McCarthy <pauldmccarthy@gmail.com>
-#
-"""The ``fsl.scripts`` package contains all of the executable scripts provided
-by ``fslpy``.
-"""

fsl/scripts/atlasq.py

@@ -19,20 +19,9 @@ import              warnings
 import              logging
 import numpy     as np
 
-# if h5py <= 2.7.1 is installed,
-# it will be imported via nibabel,
-# and will cause a numpy warning
-# to be emitted.
-with warnings.catch_warnings():
-    warnings.filterwarnings("ignore", category=FutureWarning)
-    import fsl.data.image as fslimage
-
-# If wx is not present, then fsl.utils.platform
-# will complain that it is not present.
-logging.getLogger('fsl.utils.platform').setLevel(logging.ERROR)
-
-import fsl.data.atlases as fslatlases  # noqa
-import fsl.version      as fslversion  # noqa
+import fsl.data.image   as fslimage
+import fsl.data.atlases as fslatlases
+import fsl.version      as fslversion
 
 
 log = logging.getLogger(__name__)
@@ -381,7 +370,7 @@ def maskQuery(atlas, masks, *args, **kwargs):
 
             labels = []
             props  = []
-            zprops = atlas.maskProportions(mask)
+            zprops = atlas.maskValues(mask)
 
             for i in range(len(zprops)):
                 if zprops[i] > 0:
@@ -405,7 +394,7 @@ def coordQuery(atlas, coords, voxel, *args, **kwargs):
 
         if isinstance(atlas, fslatlases.ProbabilisticAtlas):
 
-            props   = atlas.proportions(coord, voxel=voxel)
+            props   = atlas.values(coord, voxel=voxel)
             labels  = []
             nzprops = []
 
@@ -569,10 +558,10 @@ def parseArgs(args):
     usages = {
         'main'    : 'usage: atlasq [-h] command [options]',
         'ohi'     : textwrap.dedent("""
-                     usage: atlasq ohi -h
-                            atlasq ohi --dumpatlases
-                            atlasq ohi -a atlas -c X,Y,Z
-                            atlasq ohi -a atlas -m mask
+                     usage: atlasquery -h
+                            atlasquery --dumpatlases
+                            atlasquery -a atlas -c X,Y,Z
+                            atlasquery -a atlas -m mask
         """).strip(),
         'list'    : 'usage: atlasq list [-e]',
         'summary' : 'usage: atlasq summary atlas',

fsl/scripts/fsl_abspath.py

+#!/usr/bin/env python
+#
+# fsl_abspath.py - Make a path absolute
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""The fsl_abspath command - makes relative paths absolute.
+"""
+
+
+import os.path as op
+import            sys
+
+
+usage = """
+usage: fsl_abspath path
+""".strip()
+
+
+def main(argv=None):
+    """fsl_abspath - make a relative path absolute. """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) != 1:
+        print(usage)
+        return 1
+
+    print(op.realpath(argv[0]))
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/fsl_apply_x5.py

@@ -55,7 +55,7 @@ def parseArgs(args):
                         choices=('nearest', 'linear', 'cubic'),
                         default='linear'),
         'ref'    : dict(help=helps['ref'],
-                        type=ft.partial(parse_data.Image, loadData=False)),
+                        type=parse_data.Image),
     }
 
     parser.add_argument(*flags['input'],  **opts['input'])

fsl/scripts/fsl_convert_x5.py

@@ -56,7 +56,7 @@ def parseArgs(args):
     subparsers = parser.add_subparsers(dest='ctype')
     flirt      = subparsers.add_parser('flirt', epilog=epilog)
     fnirt      = subparsers.add_parser('fnirt', epilog=epilog)
-    imgtype    = ft.partial(parse_data.Image, loadData=False)
+    imgtype    = parse_data.Image
 
     flirt.add_argument('input',                  help=helps['input'])
     flirt.add_argument('output',                 help=helps['output'])

fsl/scripts/fsl_ents.py

@@ -9,8 +9,6 @@ time series from a MELODIC ``.ica`` directory.
 """
 
 
-from __future__ import print_function
-
 import os.path as op
 import            sys
 import            argparse
@@ -18,12 +16,8 @@ import            warnings
 
 import numpy   as np
 
-# See atlasq.py for explanation
-with warnings.catch_warnings():
-    warnings.filterwarnings("ignore", category=FutureWarning)
-
-    import fsl.data.fixlabels       as fixlabels
-    import fsl.data.melodicanalysis as melanalysis
+import fsl.data.fixlabels       as fixlabels
+import fsl.data.melodicanalysis as melanalysis
 
 
 DTYPE = np.float64

fsl/scripts/imcp.py

@@ -12,19 +12,13 @@ The :func:`main` function is essentially a wrapper around the
 """
 
 
-from __future__ import print_function
-
 import os.path        as op
 import                   sys
-import                   warnings
+import                   logging
 
 import fsl.utils.path as fslpath
-
-# See atlasq.py for explanation
-with warnings.catch_warnings():
-    warnings.filterwarnings("ignore", category=FutureWarning)
-    import fsl.utils.imcp as imcp
-    import fsl.data.image as fslimage
+import fsl.utils.imcp as imcp
+import fsl.data.image as fslimage
 
 
 usage = """Usage:
@@ -59,6 +53,11 @@ def main(argv=None):
         print(usage)
         return 1
 
+    # When converting to NIFTI2, nibabel
+    # emits an annoying message via log.warning:
+    #   sizeof_hdr should be 540; set sizeof_hdr to 540
+    logging.getLogger('nibabel').setLevel(logging.ERROR)
+
     try:
         srcs = [fslimage.fixExt(s) for s in srcs]
         srcs = fslpath.removeDuplicates(

fsl/scripts/imglob.py

@@ -9,17 +9,11 @@ NIFTI/ANALYZE image files.
 """
 
 
-from __future__ import print_function
-
+import itertools      as it
+import                   glob
 import                   sys
-import                   warnings
 import fsl.utils.path as fslpath
 
-# See atlasq.py for explanation
-with warnings.catch_warnings():
-    warnings.filterwarnings("ignore", category=FutureWarning)
-    import fsl.data.image as fslimage
-
 
 usage = """
 Usage: imglob [-extension/extensions] <list of names>
@@ -27,8 +21,17 @@ Usage: imglob [-extension/extensions] <list of names>
        -extensions for image list with full extensions
 """.strip()
 
-exts   = fslimage.ALLOWED_EXTENSIONS
-groups = fslimage.FILE_GROUPS
+
+# The lists below are defined in the
+# fsl.data.image class, but are duplicated
+# here for performance (to avoid import of
+# nibabel/numpy/etc).
+exts   = ['.nii.gz', '.nii', '.img', '.hdr', '.img.gz', '.hdr.gz']
+"""List of supported image file extensions. """
+
+
+groups = [('.hdr', '.img'), ('.hdr.gz', '.img.gz')]
+"""List of known image file groups (image/header file pairs). """
 
 
 def imglob(paths, output=None):
@@ -61,12 +64,38 @@ def imglob(paths, output=None):
 
     imgfiles = []
 
+    # Expand any wildcard paths if provided.
+    # Depending on the way that imglob is
+    # invoked, this may not get done by the
+    # calling shell.
+    #
+    # We also have to handle incomplete
+    # wildcards, e.g. if the user provides
+    # "img_??", we need to add possible
+    # file suffixes before it can be
+    # expanded.
+    expanded = []
+    for path in paths:
+        if any(c in path for c in '*?[]'):
+            if fslpath.hasExt(path, exts):
+                globs = [path]
+            else:
+                globs = [f'{path}{ext}' for ext in exts]
+            globs = [glob.glob(g) for g in globs]
+            expanded.extend(it.chain(*globs))
+        else:
+            expanded.append(path)
+
+    paths = expanded
+
     # Build a list of all image files (both
     # hdr and img and otherwise) that match
     for path in paths:
         try:
-            path = fslimage.removeExt(path)
-            imgfiles.extend(fslimage.addExt(path, unambiguous=False))
+            path = fslpath.removeExt(path, allowedExts=exts)
+            imgfiles.extend(fslpath.addExt(path,
+                                           allowedExts=exts,
+                                           unambiguous=False))
         except fslpath.PathError:
             continue
 

fsl/scripts/imln.py

+#!/usr/bin/env python
+#
+# imln.py - Create symbolic links to image files.
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""This module defines the ``imln`` application, for creating sym-links
+to NIFTI image files.
+
+.. note:: When creating links to relative paths, ln requires that the path is
+          relative to the link location, rather than the invocation
+          location. This is *not* currently supported by imln, and possibly
+          never will be.
+"""
+
+
+import os.path        as op
+import                   os
+import                   sys
+import fsl.utils.path as fslpath
+
+
+# The lists below are defined in the
+# fsl.data.image class, but are duplicated
+# here for performance (to avoid import of
+# nibabel/numpy/etc).
+exts = ['.nii.gz', '.nii',
+        '.img',    '.hdr',
+        '.img.gz', '.hdr.gz',
+        '.mnc',    '.mnc.gz']
+"""List of file extensions that are supported by ``imtest``.
+"""
+
+groups = [('.hdr', '.img'), ('.hdr.gz', '.img.gz')]
+"""List of known image file groups (image/header file pairs). """
+
+
+usage = """
+Usage: imln <file1> <file2>
+  Makes a link (called file2) to file1
+  NB: filenames can be basenames or include an extension
+""".strip()
+
+
+def main(argv=None):
+    """``imln`` - create sym-links to images. """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) != 2:
+        print(usage)
+        return 1
+
+    target, linkbase = argv
+    target           = fslpath.removeExt(target,   exts)
+    linkbase         = fslpath.removeExt(linkbase, exts)
+
+    # Target must exist, so we can
+    # infer the correct extension(s).
+    # Error on incomplete file groups
+    # (e.g. a.img without a.hdr).
+    try:
+        targets = fslpath.getFileGroup(target,
+                                       allowedExts=exts,
+                                       fileGroups=groups,
+                                       unambiguous=True)
+    except Exception as e:
+        print(f'Error: {e}')
+        return 1
+
+    for target in targets:
+        if not op.exists(target):
+            continue
+
+        ext  = fslpath.getExt(target, exts)
+        link = f'{linkbase}{ext}'
+
+        try:
+
+            # emulate old imln behaviour - if
+            # link already exists, it is removed
+            if op.exists(link):
+                os.remove(link)
+
+            os.symlink(target, link)
+
+        except Exception as e:
+            print(f'Error: {e}')
+            return 1
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/immv.py

@@ -13,19 +13,13 @@ The :func:`main` function is essentially a wrapper around the
 """
 
 
-from __future__ import print_function
-
 import os.path        as op
 import                   sys
-import                   warnings
+import                   logging
 
 import fsl.utils.path as fslpath
-
-# See atlasq.py for explanation
-with warnings.catch_warnings():
-    warnings.filterwarnings("ignore", category=FutureWarning)
-    import fsl.utils.imcp as imcp
-    import fsl.data.image as fslimage
+import fsl.utils.imcp as imcp
+import fsl.data.image as fslimage
 
 
 usage = """Usage:
@@ -60,6 +54,11 @@ def main(argv=None):
         print(usage)
         return 1
 
+    # When converting to NIFTI2, nibabel
+    # emits an annoying message via log.warning:
+    #   sizeof_hdr should be 540; set sizeof_hdr to 540
+    logging.getLogger('nibabel').setLevel(logging.ERROR)
+
     try:
         srcs = [fslimage.fixExt(s) for s in srcs]
         srcs = fslpath.removeDuplicates(

fsl/scripts/imrm.py

+#!/usr/bin/env python
+#
+# imrm.py - Remove image files.
+#
+# Author: Paul McCarthy <paulmc@fmrib.ox.ac.uk>
+#
+"""This module defines the ``imrm`` application, for removing NIFTI image
+files.
+"""
+
+
+import os.path            as op
+import                       os
+import                       sys
+import fsl.scripts.imglob as imglob
+
+
+usage = """Usage: imrm <list of image names to remove>
+NB: filenames can be basenames or not
+""".strip()
+
+
+def main(argv=None):
+    """Removes all images which are specified on the command line. """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) < 1:
+        print(usage)
+        return 1
+
+    paths = imglob.imglob(argv, 'all')
+
+    for path in paths:
+        if op.exists(path):
+            os.remove(path)
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/imtest.py

+#!/usr/bin/env python
+#
+# imtest.py - Test whether an image file exists or not.
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+"""The ``imtest`` script can be used to test whether an image file exists or
+not, without having to know the file suffix (.nii, .nii.gz, etc).
+"""
+
+
+import os.path        as op
+import                   sys
+import fsl.utils.path as fslpath
+
+
+# The lists below are defined in the
+# fsl.data.image class, but are duplicated
+# here for performance (to avoid import of
+# nibabel/numpy/etc).
+exts = ['.nii.gz', '.nii',
+        '.img',    '.hdr',
+        '.img.gz', '.hdr.gz',
+        '.mnc',    '.mnc.gz']
+"""List of file extensions that are supported by ``imtest``.
+"""
+
+groups = [('.hdr', '.img'), ('.hdr.gz', '.img.gz')]
+"""List of known image file groups (image/header file pairs). """
+
+
+def imtest(path):
+    """Returns ``True`` if the given image path exists, False otherwise. """
+
+    path = fslpath.removeExt(path, exts)
+    path = op.realpath(path)
+
+    # getFileGroup will raise an error
+    # if the image (including all
+    # components - i.e. header and
+    # image) does not exist
+    try:
+        fslpath.getFileGroup(path,
+                             allowedExts=exts,
+                             fileGroups=groups,
+                             unambiguous=True)
+        return True
+
+    except fslpath.PathError:
+        return False
+
+
+def main(argv=None):
+    """Test if an image path exists, and prints ``'1'`` if it does or ``'0'``
+    if it doesn't.
+    """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    # emulate old fslio/imtest - always return 0
+    if len(argv) != 1:
+        print('0')
+        return 0
+
+    if imtest(argv[0]):
+        print('1')
+    else:
+        print('0')
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

fsl/scripts/remove_ext.py

+#!/usr/bin/env python
+#
+# remove_ext.py - Remove file extensions from NIFTI image paths
+#
+# Author: Paul McCarthy <pauldmccarthy@gmail.com>
+#
+
+
+import                   sys
+import fsl.utils.path as fslpath
+
+
+usage = """Usage: remove_ext <list of image paths to remove extension from>
+""".strip()
+
+
+# This list is defined in the
+# fsl.data.image class, but are duplicated
+# here for performance (to avoid import of
+# nibabel/numpy/etc).
+exts = ['.nii.gz', '.nii',
+        '.img',    '.hdr',
+        '.img.gz', '.hdr.gz',
+        '.mnc',    '.mnc.gz']
+"""List of file extensions that are removed by ``remove_ext``. """
+
+
+def main(argv=None):
+    """Removes file extensions from all paths which are specified on the
+    command line.
+    """
+
+    if argv is None:
+        argv = sys.argv[1:]
+
+    if len(argv) < 1:
+        print(usage)
+        return 1
+
+    removed = []
+
+    for path in argv:
+        removed.append(fslpath.removeExt(path, exts))
+
+    print(' '.join(removed))
+
+    return 0
+
+
+if __name__ == '__main__':
+    sys.exit(main())

tests/__init__.py → fsl/tests/__init__.py

@@ -10,6 +10,7 @@
 import              os
 import              sys
 import              glob
+import              hashlib
 import              shutil
 import              fnmatch
 import              logging
@@ -20,14 +21,12 @@ import os.path   as op
 import numpy     as np
 import nibabel   as nib
 
-from six import StringIO
+from io import StringIO
 
-
-try: from unittest import mock
-except ImportError: import mock
+from unittest import mock
 
 import fsl.data.image                     as fslimage
-from   fsl.utils.tempdir import              tempdir
+from   fsl.utils.tempdir  import             tempdir
 from   fsl.utils.platform import platform as fslplatform
 
 
@@ -50,7 +49,10 @@ def mockFSLDIR(**kwargs):
                 if not op.isdir(subdir):
                     os.makedirs(subdir)
                 for fname in files:
-                    touch(op.join(subdir, fname))
+                    fname = op.join(subdir, fname)
+                    touch(fname)
+                    if subdir == bindir:
+                        os.chmod(fname, 0o755)
             fslplatform.fsldir = fsldir
             fslplatform.fsldevdir = None
 
@@ -68,7 +70,7 @@ def touch(fname):
         pass
 
 
-class CaptureStdout(object):
+class CaptureStdout:
     """Context manager which captures stdout and stderr. """
 
     def __init__(self):
@@ -87,6 +89,7 @@ class CaptureStdout(object):
 
         sys.stdout = self.__mock_stdout
         sys.stderr = self.__mock_stderr
+        return self
 
 
     def __exit__(self, *args, **kwargs):
@@ -145,6 +148,8 @@ def testdir(contents=None, suffix=""):
             shutil.rmtree(self.testdir)
 
     return ctx(contents)
+testdir.__test__ = False
+
 
 def make_dummy_files(paths):
     """Creates dummy files for all of the given paths. """
@@ -284,7 +289,8 @@ def make_mock_feat_analysis(featdir,
                             copes=True,
                             zstats=True,
                             residuals=True,
-                            clustMasks=True):
+                            clustMasks=True,
+                            zfstats=True):
 
     if xform is None:
         xform = np.eye(4)
@@ -317,6 +323,7 @@ def make_mock_feat_analysis(featdir,
             data = np.ravel_multi_index(data, shape)
             data = data.reshape(list(shape) + [1]).repeat(timepoints, axis=3)
             data[..., :] += range(i, i + timepoints)
+            data = data.astype(np.int32)
 
             img = nib.nifti1.Nifti1Image(data, xform)
 
@@ -341,6 +348,11 @@ def make_mock_feat_analysis(featdir,
         otherFiles .extend(files)
         otherShapes.extend([shape] * len(files))
 
+    if zfstats:
+        files = glob.glob(op.join(featdir, 'stats', 'zfstat*nii.gz'))
+        otherFiles .extend(files)
+        otherShapes.extend([shape] * len(files))
+
     if residuals:
         files = glob.glob(op.join(featdir, 'stats', 'res4d.nii.gz'))
         otherFiles .extend(files)
@@ -428,3 +440,10 @@ def make_random_mask(filename, shape, xform, premask=None, minones=1):
     img.save(filename)
 
     return img
+
+
+def sha256(filename):
+    hashobj = hashlib.sha256()
+    with open(filename, 'rb') as f:
+        hashobj.update(f.read())
+    return hashobj.hexdigest()