From c6c6a1de03d8838868af5c0540af9ef37c842284 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Tue, 24 Feb 2015 01:34:18 -0800
Subject: [PATCH 01/16] Add an 'axis' parameter to concat_images, plus two
 tests.

---
 nibabel/funcs.py            | 17 ++++++++++++-----
 nibabel/tests/test_funcs.py | 17 +++++++++++++++++
 2 files changed, 29 insertions(+), 5 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 645fe09b2b..5a927ffc1b 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -88,8 +88,8 @@ def squeeze_image(img):
                  img.extra)
 
 
-def concat_images(images, check_affines=True):
-    ''' Concatenate images in list to single image, along last dimension
+def concat_images(images, check_affines=True, axis=None):
+    ''' Concatenate images in list to single image, along specified dimension
 
     Parameters
     ----------
@@ -98,7 +98,9 @@ def concat_images(images, check_affines=True):
     check_affines : {True, False}, optional
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
-
+    axis : int, optional
+        If None, concatenates on the last dimension.
+        If not None, concatenates on the specified dimension.
     Returns
     -------
     concat_img : ``SpatialImage``
@@ -122,8 +124,13 @@ def concat_images(images, check_affines=True):
         if check_affines:
             if not np.all(img.affine == affine):
                 raise ValueError('Affines do not match')
-        out_data[i] = img.get_data()
-    out_data = np.rollaxis(out_data, 0, len(i0shape)+1)
+        out_data[i] = img.get_data().copy()
+    if axis is not None:
+        out_data = np.concatenate(out_data, axis=axis)
+    elif np.all([d.shape[-1] == 1 for d in out_data]):
+        out_data = np.concatenate(out_data, axis=d.ndim-1)
+    else:
+        out_data = np.rollaxis(out_data, 0, len(i0shape)+1)
     klass = img0.__class__
     return klass(out_data, affine, header)
 
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 8a18afd739..d80b77b1dc 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -59,6 +59,23 @@ def test_concat():
         for img in imgs:
             del(img)
 
+    # Test axis parameter and trailing unary dimension
+    shape_4D = np.asarray(shape + (1,))
+    data0 = np.arange(10).reshape(shape_4D)
+    affine = np.eye(4)
+    img0_mem = Nifti1Image(data0, affine)
+    img1_mem = Nifti1Image(data0 - 10, affine)
+
+    concat_img1 = concat_images([img0_mem, img1_mem])
+    expected_shape1 = shape_4D.copy()
+    expected_shape1[-1] *= 2
+    assert_array_equal(concat_img1.shape, expected_shape1)
+
+    concat_img2 = concat_images([img0_mem, img1_mem], axis=0)
+    expected_shape2 = shape_4D.copy()
+    expected_shape2[0] *= 2
+    assert_array_equal(concat_img2.shape, expected_shape2)
+
 
 def test_closest_canonical():
     arr = np.arange(24).reshape((2,3,4,1))

From 2626884a21f6071a693886abdec0c7f1393cde3a Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Tue, 24 Feb 2015 15:42:05 -0800
Subject: [PATCH 02/16] Try again, this time with lists and more tests...

---
 nibabel/funcs.py            |  7 +++----
 nibabel/tests/test_funcs.py | 21 ++++++++++++++++++++-
 2 files changed, 23 insertions(+), 5 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 5a927ffc1b..da99543e4f 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -117,19 +117,18 @@ def concat_images(images, check_affines=True, axis=None):
     affine = img0.affine
     header = img0.header
     out_shape = (n_imgs, ) + i0shape
-    out_data = np.empty(out_shape)
+    out_data = []
     for i, img in enumerate(images):
         if is_filename:
             img = load(img)
         if check_affines:
             if not np.all(img.affine == affine):
                 raise ValueError('Affines do not match')
-        out_data[i] = img.get_data().copy()
+        out_data.append(img.get_data())
     if axis is not None:
         out_data = np.concatenate(out_data, axis=axis)
-    elif np.all([d.shape[-1] == 1 for d in out_data]):
-        out_data = np.concatenate(out_data, axis=d.ndim-1)
     else:
+        out_data = np.asarray(out_data)
         out_data = np.rollaxis(out_data, 0, len(i0shape)+1)
     klass = img0.__class__
     return klass(out_data, affine, header)
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index d80b77b1dc..3949b64e52 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -66,16 +66,35 @@ def test_concat():
     img0_mem = Nifti1Image(data0, affine)
     img1_mem = Nifti1Image(data0 - 10, affine)
 
-    concat_img1 = concat_images([img0_mem, img1_mem])
+    # 4d, same shape, append on axis 3
+    concat_img1 = concat_images([img0_mem, img1_mem], axis=3)
     expected_shape1 = shape_4D.copy()
     expected_shape1[-1] *= 2
     assert_array_equal(concat_img1.shape, expected_shape1)
 
+    # 4d, same shape, append on axis 0
     concat_img2 = concat_images([img0_mem, img1_mem], axis=0)
     expected_shape2 = shape_4D.copy()
     expected_shape2[0] *= 2
     assert_array_equal(concat_img2.shape, expected_shape2)
 
+    # 4d, same shape, append on axis -1
+    concat_img3 = concat_images([img0_mem, img1_mem], axis=-1)
+    expected_shape3 = shape_4D.copy()
+    expected_shape3[-1] *= 2
+    assert_array_equal(concat_img3.shape, expected_shape3)
+
+    # 4d, different shape, append on axis that's different
+    print('%s %s' % (str(concat_img3.shape), str(img1_mem.shape)))
+    concat_img4 = concat_images([concat_img3, img1_mem], axis=-1)
+    expected_shape4 = shape_4D.copy()
+    expected_shape4[-1] *= 3
+    assert_array_equal(concat_img4.shape, expected_shape4)
+
+    # 4d, different shape, append on axis that's not different...
+    # Doesn't work!
+    assert_raises(ValueError, concat_images, [concat_img3, img1_mem], axis=1)
+
 
 def test_closest_canonical():
     arr = np.arange(24).reshape((2,3,4,1))

From c39caac72baa017bfe8b47f0226d931a32b52b03 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Tue, 24 Feb 2015 15:43:43 -0800
Subject: [PATCH 03/16] Add greater coverage of different shapes.

---
 nibabel/tests/test_funcs.py | 117 ++++++++++++++++++------------------
 1 file changed, 59 insertions(+), 58 deletions(-)

diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 3949b64e52..22608240d2 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -30,70 +30,71 @@ def _as_fname(img):
 
 
 def test_concat():
-    shape = (1,2,5)
-    data0 = np.arange(10).reshape(shape)
-    affine = np.eye(4)
-    img0_mem = Nifti1Image(data0, affine)
-    data1 = data0 - 10
-    img1_mem = Nifti1Image(data1, affine)
-    img2_mem = Nifti1Image(data1, affine+1)
-    img3_mem = Nifti1Image(data1.T, affine)
-    all_data = np.concatenate(
-        [data0[:,:,:,np.newaxis],data1[:,:,:,np.newaxis]],3)
-    # Check filenames and in-memory images work
-    with InTemporaryDirectory():
-        imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
-        img_files = [_as_fname(img) for img in imgs]
-        for img0, img1, img2, img3 in (imgs, img_files):
-            all_imgs = concat_images([img0, img1])
-            assert_array_equal(all_imgs.get_data(), all_data)
-            assert_array_equal(all_imgs.affine, affine)
-            # check that not-matching affines raise error
-            assert_raises(ValueError, concat_images, [img0, img2])
-            assert_raises(ValueError, concat_images, [img0, img3])
-            # except if check_affines is False
-            all_imgs = concat_images([img0, img1])
-            assert_array_equal(all_imgs.get_data(), all_data)
-            assert_array_equal(all_imgs.affine, affine)
-        # Delete images as prophylaxis for windows access errors
-        for img in imgs:
-            del(img)
+    for shape in ((1,2,5), (7,3,1), (13,11,11), (0,1,1)):
+        numel = np.asarray(shape).prod()
+        data0 = np.arange(numel).reshape(shape)
+        affine = np.eye(4)
+        img0_mem = Nifti1Image(data0, affine)
+        data1 = data0 - 10
+        img1_mem = Nifti1Image(data1, affine)
+        img2_mem = Nifti1Image(data1, affine+1)
+        img3_mem = Nifti1Image(data1.T, affine)
+        all_data = np.concatenate(
+            [data0[:,:,:,np.newaxis],data1[:,:,:,np.newaxis]],3)
+        # Check filenames and in-memory images work
+        with InTemporaryDirectory():
+            imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
+            img_files = [_as_fname(img) for img in imgs]
+            for img0, img1, img2, img3 in (imgs, img_files):
+                all_imgs = concat_images([img0, img1])
+                assert_array_equal(all_imgs.get_data(), all_data)
+                assert_array_equal(all_imgs.affine, affine)
+                # check that not-matching affines raise error
+                assert_raises(ValueError, concat_images, [img0, img2])
+                assert_raises(ValueError, concat_images, [img0, img3])
+                # except if check_affines is False
+                all_imgs = concat_images([img0, img1])
+                assert_array_equal(all_imgs.get_data(), all_data)
+                assert_array_equal(all_imgs.affine, affine)
+            # Delete images as prophylaxis for windows access errors
+            for img in imgs:
+                del(img)
 
-    # Test axis parameter and trailing unary dimension
-    shape_4D = np.asarray(shape + (1,))
-    data0 = np.arange(10).reshape(shape_4D)
-    affine = np.eye(4)
-    img0_mem = Nifti1Image(data0, affine)
-    img1_mem = Nifti1Image(data0 - 10, affine)
+        # Test axis parameter and trailing unary dimension
+        shape_4D = np.asarray(shape + (1,))
+        data0 = np.arange(numel).reshape(shape_4D)
+        affine = np.eye(4)
+        img0_mem = Nifti1Image(data0, affine)
+        img1_mem = Nifti1Image(data0 - 10, affine)
 
-    # 4d, same shape, append on axis 3
-    concat_img1 = concat_images([img0_mem, img1_mem], axis=3)
-    expected_shape1 = shape_4D.copy()
-    expected_shape1[-1] *= 2
-    assert_array_equal(concat_img1.shape, expected_shape1)
+        # 4d, same shape, append on axis 3
+        concat_img1 = concat_images([img0_mem, img1_mem], axis=3)
+        expected_shape1 = shape_4D.copy()
+        expected_shape1[-1] *= 2
+        assert_array_equal(concat_img1.shape, expected_shape1)
 
-    # 4d, same shape, append on axis 0
-    concat_img2 = concat_images([img0_mem, img1_mem], axis=0)
-    expected_shape2 = shape_4D.copy()
-    expected_shape2[0] *= 2
-    assert_array_equal(concat_img2.shape, expected_shape2)
+        # 4d, same shape, append on axis 0
+        concat_img2 = concat_images([img0_mem, img1_mem], axis=0)
+        expected_shape2 = shape_4D.copy()
+        expected_shape2[0] *= 2
+        assert_array_equal(concat_img2.shape, expected_shape2)
 
-    # 4d, same shape, append on axis -1
-    concat_img3 = concat_images([img0_mem, img1_mem], axis=-1)
-    expected_shape3 = shape_4D.copy()
-    expected_shape3[-1] *= 2
-    assert_array_equal(concat_img3.shape, expected_shape3)
+        # 4d, same shape, append on axis -1
+        concat_img3 = concat_images([img0_mem, img1_mem], axis=-1)
+        expected_shape3 = shape_4D.copy()
+        expected_shape3[-1] *= 2
+        assert_array_equal(concat_img3.shape, expected_shape3)
 
-    # 4d, different shape, append on axis that's different
-    print('%s %s' % (str(concat_img3.shape), str(img1_mem.shape)))
-    concat_img4 = concat_images([concat_img3, img1_mem], axis=-1)
-    expected_shape4 = shape_4D.copy()
-    expected_shape4[-1] *= 3
-    assert_array_equal(concat_img4.shape, expected_shape4)
+        # 4d, different shape, append on axis that's different
+        print('%s %s' % (str(concat_img3.shape), str(img1_mem.shape)))
+        concat_img4 = concat_images([concat_img3, img1_mem], axis=-1)
+        expected_shape4 = shape_4D.copy()
+        expected_shape4[-1] *= 3
+        assert_array_equal(concat_img4.shape, expected_shape4)
 
-    # 4d, different shape, append on axis that's not different...
-    # Doesn't work!
-    assert_raises(ValueError, concat_images, [concat_img3, img1_mem], axis=1)
+        # 4d, different shape, append on axis that's not different...
+        # Doesn't work!
+        assert_raises(ValueError, concat_images, [concat_img3, img1_mem], axis=1)
 
 
 def test_closest_canonical():

From afaa5fe9a789e212dad7519f0a43f417848bda48 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Wed, 11 Mar 2015 14:06:45 -0700
Subject: [PATCH 04/16] Make this work for all 3D and 4D combinations possible,
 across all axes possible.  Test extensively.

---
 nibabel/funcs.py            |  45 ++++++----
 nibabel/tests/test_funcs.py | 158 +++++++++++++++++++++---------------
 2 files changed, 125 insertions(+), 78 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index da99543e4f..f5be5856c4 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -98,9 +98,9 @@ def concat_images(images, check_affines=True, axis=None):
     check_affines : {True, False}, optional
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
-    axis : int, optional
-        If None, concatenates on the last dimension.
-        If not None, concatenates on the specified dimension.
+    axis : None or int, optional
+        If None, concatenates on the 4th dimension.
+        If not None, concatenates on the specified dimension [-2 to 3).
     Returns
     -------
     concat_img : ``SpatialImage``
@@ -113,23 +113,40 @@ def concat_images(images, check_affines=True, axis=None):
     if not hasattr(img0, 'get_data'):
         img0 = load(img0)
         is_filename = True
-    i0shape = img0.shape
     affine = img0.affine
     header = img0.header
-    out_shape = (n_imgs, ) + i0shape
-    out_data = []
+
+    if axis is None:  # collect images in output array for efficiency
+        out_shape = (n_imgs, ) + img0.shape[:3]
+        out_data = np.empty(out_shape)
+    else:  # collect images in list for use with np.concatenate
+        out_data = [None] * n_imgs
+
     for i, img in enumerate(images):
         if is_filename:
             img = load(img)
-        if check_affines:
-            if not np.all(img.affine == affine):
-                raise ValueError('Affines do not match')
-        out_data.append(img.get_data())
-    if axis is not None:
-        out_data = np.concatenate(out_data, axis=axis)
+        if check_affines and not np.all(img.affine == affine):
+            raise ValueError('Affines do not match')
+
+        if axis is None and img.get_data().ndim == 4 and img.get_data().shape[3] == 1:
+            out_data[i] = np.reshape(img.get_data(), img.get_data().shape[:-1])
+        else:
+            out_data[i] = img.get_data()
+
+        if is_filename:
+            del img
+
+    if axis is None:
+        out_data = np.rollaxis(out_data, 0, out_data.ndim)
     else:
-        out_data = np.asarray(out_data)
-        out_data = np.rollaxis(out_data, 0, len(i0shape)+1)
+        # Massage the output, to allow combining 3D and 4D images.
+        is_3D = [len(d.shape) == 3 for d in out_data]
+        is_4D = [len(d.shape) == 4 for d in out_data]
+        if np.any(is_3D) and np.any(is_4D):
+            out_data = [data if is_4D[di] else np.reshape(data, data.shape + (1,))
+                        for di, data in enumerate(out_data)]
+        out_data = np.concatenate(out_data, axis=axis)
+
     klass = img0.__class__
     return klass(out_data, affine, header)
 
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 22608240d2..c6780ef4ce 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -30,71 +30,101 @@ def _as_fname(img):
 
 
 def test_concat():
-    for shape in ((1,2,5), (7,3,1), (13,11,11), (0,1,1)):
-        numel = np.asarray(shape).prod()
-        data0 = np.arange(numel).reshape(shape)
-        affine = np.eye(4)
-        img0_mem = Nifti1Image(data0, affine)
-        data1 = data0 - 10
-        img1_mem = Nifti1Image(data1, affine)
-        img2_mem = Nifti1Image(data1, affine+1)
-        img3_mem = Nifti1Image(data1.T, affine)
-        all_data = np.concatenate(
-            [data0[:,:,:,np.newaxis],data1[:,:,:,np.newaxis]],3)
-        # Check filenames and in-memory images work
-        with InTemporaryDirectory():
-            imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
-            img_files = [_as_fname(img) for img in imgs]
-            for img0, img1, img2, img3 in (imgs, img_files):
-                all_imgs = concat_images([img0, img1])
-                assert_array_equal(all_imgs.get_data(), all_data)
-                assert_array_equal(all_imgs.affine, affine)
-                # check that not-matching affines raise error
-                assert_raises(ValueError, concat_images, [img0, img2])
-                assert_raises(ValueError, concat_images, [img0, img3])
-                # except if check_affines is False
-                all_imgs = concat_images([img0, img1])
-                assert_array_equal(all_imgs.get_data(), all_data)
-                assert_array_equal(all_imgs.affine, affine)
-            # Delete images as prophylaxis for windows access errors
-            for img in imgs:
-                del(img)
-
-        # Test axis parameter and trailing unary dimension
-        shape_4D = np.asarray(shape + (1,))
-        data0 = np.arange(numel).reshape(shape_4D)
-        affine = np.eye(4)
+
+    # Build combinations of 3D, 4D w/size[3] == 1, and 4D w/size[3] == 3
+    all_shapes_3D = ((1, 2, 5), (7, 3, 1), (13, 11, 11), (0, 1, 1))
+    all_shapes_4D_unary = tuple((shape + (1,) for shape in all_shapes_3D))
+    all_shapes_4D_multi = tuple((shape + (3,) for shape in all_shapes_3D))
+    all_shapes = all_shapes_3D + all_shapes_4D_unary + all_shapes_4D_multi
+
+    affine = np.eye(4)
+    # Loop over all possible combinations of images, in first and
+    #   second position.
+    for data0_shape in all_shapes:
+        data0_numel = np.asarray(data0_shape).prod()
+        data0 = np.arange(data0_numel).reshape(data0_shape)
         img0_mem = Nifti1Image(data0, affine)
-        img1_mem = Nifti1Image(data0 - 10, affine)
-
-        # 4d, same shape, append on axis 3
-        concat_img1 = concat_images([img0_mem, img1_mem], axis=3)
-        expected_shape1 = shape_4D.copy()
-        expected_shape1[-1] *= 2
-        assert_array_equal(concat_img1.shape, expected_shape1)
-
-        # 4d, same shape, append on axis 0
-        concat_img2 = concat_images([img0_mem, img1_mem], axis=0)
-        expected_shape2 = shape_4D.copy()
-        expected_shape2[0] *= 2
-        assert_array_equal(concat_img2.shape, expected_shape2)
-
-        # 4d, same shape, append on axis -1
-        concat_img3 = concat_images([img0_mem, img1_mem], axis=-1)
-        expected_shape3 = shape_4D.copy()
-        expected_shape3[-1] *= 2
-        assert_array_equal(concat_img3.shape, expected_shape3)
-
-        # 4d, different shape, append on axis that's different
-        print('%s %s' % (str(concat_img3.shape), str(img1_mem.shape)))
-        concat_img4 = concat_images([concat_img3, img1_mem], axis=-1)
-        expected_shape4 = shape_4D.copy()
-        expected_shape4[-1] *= 3
-        assert_array_equal(concat_img4.shape, expected_shape4)
-
-        # 4d, different shape, append on axis that's not different...
-        # Doesn't work!
-        assert_raises(ValueError, concat_images, [concat_img3, img1_mem], axis=1)
+
+        for data1_shape in all_shapes:
+            data1_numel = np.asarray(data1_shape).prod()
+            data1 = np.arange(data1_numel).reshape(data1_shape)
+            img1_mem = Nifti1Image(data1, affine)
+            img2_mem = Nifti1Image(data1, affine+1)  # bad affine
+            img3_mem = Nifti1Image(data1.T, affine)  # bad data shape
+
+            # Loop over every possible axis, including None (explicit and implied)
+            for axis in (list(range(-2, 3)) + [None, '__default__']):
+
+                # Allow testing default vs. passing explicit param
+                if axis == '__default__':
+                    np_concat_kwargs = dict(axis=-1)
+                    concat_imgs_kwargs = dict()
+                    axis = None  # Convert downstream
+                elif axis is None:
+                    np_concat_kwargs = dict(axis=-1)
+                    concat_imgs_kwargs = dict(axis=axis)
+                else:
+                    np_concat_kwargs = dict(axis=axis)
+                    concat_imgs_kwargs = dict(axis=axis)
+
+                # Create expected output
+                try:
+                    # Error will be thrown if the np.concatenate fails.
+                    #   However, when axis=None, the concatenate is possible
+                    #   but our efficient logic (where all images are
+                    #   3D and the same size) fails, so we also
+                    #   have to expect errors for those.
+                    expect_error = False
+                    if data0.ndim == 3 and data1.ndim == 4:
+                        expect_error = axis is None and data1.shape[3] != 1
+                        all_data = np.concatenate([data0[..., np.newaxis], data1],
+                                                  **np_concat_kwargs)
+                    elif data0.ndim == 4 and data1.ndim == 3:
+                        expect_error = axis is None and data0.shape[3] != 1
+                        all_data = np.concatenate([data0, data1[..., np.newaxis]],
+                                                  **np_concat_kwargs)
+                    elif data0.ndim == 4 and data1.ndim == 4:
+                        expect_error = axis is None and (data0.shape[3] != 1 or
+                                                         data1.shape[3] != 1)
+                        all_data = np.concatenate([data0, data1],
+                                                  **np_concat_kwargs)
+                    elif axis is None:  # 3D from here and below
+                        all_data = np.concatenate(
+                            [data0[..., np.newaxis], data1[..., np.newaxis]], 3)
+                    else:  # both 3D, appending on final axis
+                        all_data = np.concatenate([data0, data1],
+                                                  **np_concat_kwargs)
+                except ValueError:
+                    # Shapes are not combinable
+                    expect_error = True
+
+                # Check filenames and in-memory images work
+                with InTemporaryDirectory():
+                    # Try mem-based, file-based, and mixed
+                    imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
+                    img_files = [_as_fname(img) for img in imgs]
+                    for img0, img1, img2, img3 in (imgs, img_files):
+                        try:
+                            all_imgs = concat_images([img0, img1],
+                                                     **concat_imgs_kwargs)
+                            assert_array_equal(all_imgs.get_data(), all_data)
+                            assert_array_equal(all_imgs.affine, affine)
+                            assert_false(expect_error, "Expected a concatenation error, but got none.")
+                        except ValueError as ve:
+                            assert_true(expect_error, ve.message)
+
+                        # check that not-matching affines raise error
+                        assert_raises(ValueError, concat_images, [img0, img2], **concat_imgs_kwargs)
+                        assert_raises(ValueError, concat_images, [img0, img3], **concat_imgs_kwargs)
+
+                        # except if check_affines is False
+                        try:
+                            all_imgs = concat_images([img0, img1], **concat_imgs_kwargs)
+                            assert_array_equal(all_imgs.get_data(), all_data)
+                            assert_array_equal(all_imgs.affine, affine)
+                            assert_false(expect_error, "Expected a concatenation error, but got none.")
+                        except ValueError as ve:
+                            assert_true(expect_error, ve.message)
 
 
 def test_closest_canonical():

From 3331a51d928298e53eb799e37c66527e1f07f0f0 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Wed, 11 Mar 2015 14:07:28 -0700
Subject: [PATCH 05/16] Allow mixed files and objects.

---
 nibabel/funcs.py            | 15 ++++++++-------
 nibabel/tests/test_funcs.py |  3 ++-
 2 files changed, 10 insertions(+), 8 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index f5be5856c4..91d47db8ec 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -94,7 +94,7 @@ def concat_images(images, check_affines=True, axis=None):
     Parameters
     ----------
     images : sequence
-       sequence of ``SpatialImage`` or of filenames\s
+       sequence of ``SpatialImage`` or filenames\s
     check_affines : {True, False}, optional
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
@@ -109,32 +109,33 @@ def concat_images(images, check_affines=True, axis=None):
     '''
     n_imgs = len(images)
     img0 = images[0]
-    is_filename = False
     if not hasattr(img0, 'get_data'):
         img0 = load(img0)
-        is_filename = True
     affine = img0.affine
     header = img0.header
+    i0shape = img0.shape
+    del img0
 
     if axis is None:  # collect images in output array for efficiency
-        out_shape = (n_imgs, ) + img0.shape[:3]
+        out_shape = (n_imgs, ) + i0shape[:3]
         out_data = np.empty(out_shape)
     else:  # collect images in list for use with np.concatenate
         out_data = [None] * n_imgs
 
     for i, img in enumerate(images):
-        if is_filename:
+        if not hasattr(img, 'get_data'):
             img = load(img)
+
         if check_affines and not np.all(img.affine == affine):
             raise ValueError('Affines do not match')
 
+        # Special case for 4D image with size[3] == 1; reshape to work!
         if axis is None and img.get_data().ndim == 4 and img.get_data().shape[3] == 1:
             out_data[i] = np.reshape(img.get_data(), img.get_data().shape[:-1])
         else:
             out_data[i] = img.get_data()
 
-        if is_filename:
-            del img
+        del img
 
     if axis is None:
         out_data = np.rollaxis(out_data, 0, out_data.ndim)
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index c6780ef4ce..e4b344904e 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -103,7 +103,8 @@ def test_concat():
                     # Try mem-based, file-based, and mixed
                     imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
                     img_files = [_as_fname(img) for img in imgs]
-                    for img0, img1, img2, img3 in (imgs, img_files):
+                    imgs_mixed = [imgs[0], img_files[1], imgs[2], img_files[3]]
+                    for img0, img1, img2, img3 in (imgs, img_files, imgs_mixed):
                         try:
                             all_imgs = concat_images([img0, img1],
                                                      **concat_imgs_kwargs)

From 41732f447f9e9bdf68aa40a7f07ebcaa0e790d0a Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Wed, 11 Mar 2015 14:15:43 -0700
Subject: [PATCH 06/16] Improve efficiency: load img0 once, del reference, and
 don't check affine on first image.

---
 nibabel/funcs.py            | 30 +++++++++++++++---------------
 nibabel/tests/test_funcs.py |  3 +++
 2 files changed, 18 insertions(+), 15 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 91d47db8ec..b240af6728 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -107,26 +107,27 @@ def concat_images(images, check_affines=True, axis=None):
        New image resulting from concatenating `images` across last
        dimension
     '''
+
     n_imgs = len(images)
-    img0 = images[0]
-    if not hasattr(img0, 'get_data'):
-        img0 = load(img0)
-    affine = img0.affine
-    header = img0.header
-    i0shape = img0.shape
-    del img0
-
-    if axis is None:  # collect images in output array for efficiency
-        out_shape = (n_imgs, ) + i0shape[:3]
-        out_data = np.empty(out_shape)
-    else:  # collect images in list for use with np.concatenate
-        out_data = [None] * n_imgs
+    if n_imgs == 0:
+        raise ValueError('Cannot concatenate an empty list of images.')
 
     for i, img in enumerate(images):
         if not hasattr(img, 'get_data'):
             img = load(img)
 
-        if check_affines and not np.all(img.affine == affine):
+        if i == 0:  # first image, initialize data from loaded image
+            affine = img.affine
+            header = img.header
+            klass = img.__class__
+
+            if axis is None:  # collect images in output array for efficiency
+                out_shape = (n_imgs, ) + img.shape[:3]
+                out_data = np.empty(out_shape)
+            else:  # collect images in list for use with np.concatenate
+                out_data = [None] * n_imgs
+
+        elif check_affines and not np.all(img.affine == affine):
             raise ValueError('Affines do not match')
 
         # Special case for 4D image with size[3] == 1; reshape to work!
@@ -148,7 +149,6 @@ def concat_images(images, check_affines=True, axis=None):
                         for di, data in enumerate(out_data)]
         out_data = np.concatenate(out_data, axis=axis)
 
-    klass = img0.__class__
     return klass(out_data, affine, header)
 
 
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index e4b344904e..ba6b3eb5dd 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -30,6 +30,8 @@ def _as_fname(img):
 
 
 def test_concat():
+    # Smoke test: concat empty list.
+    assert_raises(ValueError, concat_images, [])
 
     # Build combinations of 3D, 4D w/size[3] == 1, and 4D w/size[3] == 3
     all_shapes_3D = ((1, 2, 5), (7, 3, 1), (13, 11, 11), (0, 1, 1))
@@ -128,6 +130,7 @@ def test_concat():
                             assert_true(expect_error, ve.message)
 
 
+
 def test_closest_canonical():
     arr = np.arange(24).reshape((2,3,4,1))
     # no funky stuff, returns same thing

From 47fc8f0dcdfa98467d0c172f0c0213b2812e76a1 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Wed, 11 Mar 2015 14:16:45 -0700
Subject: [PATCH 07/16] Add a final comment.

---
 nibabel/funcs.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index b240af6728..540b6d6355 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -144,7 +144,7 @@ def concat_images(images, check_affines=True, axis=None):
         # Massage the output, to allow combining 3D and 4D images.
         is_3D = [len(d.shape) == 3 for d in out_data]
         is_4D = [len(d.shape) == 4 for d in out_data]
-        if np.any(is_3D) and np.any(is_4D):
+        if np.any(is_3D) and np.any(is_4D):  # Convert all to 4D
             out_data = [data if is_4D[di] else np.reshape(data, data.shape + (1,))
                         for di, data in enumerate(out_data)]
         out_data = np.concatenate(out_data, axis=axis)

From f8ad0786b9d1c69291046562875e46ed97e7ac56 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 14:45:09 -0700
Subject: [PATCH 08/16] Code reorganization, comment editing.

---
 nibabel/funcs.py            | 11 +++++++----
 nibabel/tests/test_funcs.py | 15 ++++++++-------
 2 files changed, 15 insertions(+), 11 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 540b6d6355..44d4bbbdcc 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -94,13 +94,16 @@ def concat_images(images, check_affines=True, axis=None):
     Parameters
     ----------
     images : sequence
-       sequence of ``SpatialImage`` or filenames\s
+       sequence of ``SpatialImage`` or filenames of the same dimensionality\s
     check_affines : {True, False}, optional
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
     axis : None or int, optional
-        If None, concatenates on the 4th dimension.
-        If not None, concatenates on the specified dimension [-2 to 3).
+        If None, concatenates on a new dimension.  This rrequires all images
+           to be the same shape).
+        If not None, concatenates on the specified dimension.  This requires
+           all images to be the same shape, except on the specified dimension.
+           For 4D images, axis must be between -2 and 3.
     Returns
     -------
     concat_img : ``SpatialImage``
@@ -110,7 +113,7 @@ def concat_images(images, check_affines=True, axis=None):
 
     n_imgs = len(images)
     if n_imgs == 0:
-        raise ValueError('Cannot concatenate an empty list of images.')
+        raise ValueError("Cannot concatenate an empty list of images.")
 
     for i, img in enumerate(images):
         if not hasattr(img, 'get_data'):
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index ba6b3eb5dd..749c52a390 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -110,11 +110,12 @@ def test_concat():
                         try:
                             all_imgs = concat_images([img0, img1],
                                                      **concat_imgs_kwargs)
-                            assert_array_equal(all_imgs.get_data(), all_data)
-                            assert_array_equal(all_imgs.affine, affine)
-                            assert_false(expect_error, "Expected a concatenation error, but got none.")
                         except ValueError as ve:
                             assert_true(expect_error, ve.message)
+                        else:
+                            assert_false(expect_error, "Expected a concatenation error, but got none.")
+                            assert_array_equal(all_imgs.get_data(), all_data)
+                            assert_array_equal(all_imgs.affine, affine)
 
                         # check that not-matching affines raise error
                         assert_raises(ValueError, concat_images, [img0, img2], **concat_imgs_kwargs)
@@ -123,12 +124,12 @@ def test_concat():
                         # except if check_affines is False
                         try:
                             all_imgs = concat_images([img0, img1], **concat_imgs_kwargs)
-                            assert_array_equal(all_imgs.get_data(), all_data)
-                            assert_array_equal(all_imgs.affine, affine)
-                            assert_false(expect_error, "Expected a concatenation error, but got none.")
                         except ValueError as ve:
                             assert_true(expect_error, ve.message)
-
+                        else:
+                            assert_false(expect_error, "Expected a concatenation error, but got none.")
+                            assert_array_equal(all_imgs.get_data(), all_data)
+                            assert_array_equal(all_imgs.affine, affine)
 
 
 def test_closest_canonical():

From 84640b39b8f8a8804f8e86336978ad4fd4c44eb2 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 14:49:01 -0700
Subject: [PATCH 09/16] Remove 3D/4D special-case code.

---
 nibabel/funcs.py            | 14 ++------------
 nibabel/tests/test_funcs.py | 20 +++-----------------
 2 files changed, 5 insertions(+), 29 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 44d4bbbdcc..cd1d71f772 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -125,7 +125,7 @@ def concat_images(images, check_affines=True, axis=None):
             klass = img.__class__
 
             if axis is None:  # collect images in output array for efficiency
-                out_shape = (n_imgs, ) + img.shape[:3]
+                out_shape = (n_imgs, ) + img.shape
                 out_data = np.empty(out_shape)
             else:  # collect images in list for use with np.concatenate
                 out_data = [None] * n_imgs
@@ -133,23 +133,13 @@ def concat_images(images, check_affines=True, axis=None):
         elif check_affines and not np.all(img.affine == affine):
             raise ValueError('Affines do not match')
 
-        # Special case for 4D image with size[3] == 1; reshape to work!
-        if axis is None and img.get_data().ndim == 4 and img.get_data().shape[3] == 1:
-            out_data[i] = np.reshape(img.get_data(), img.get_data().shape[:-1])
-        else:
-            out_data[i] = img.get_data()
+        out_data[i] = img.get_data()
 
         del img
 
     if axis is None:
         out_data = np.rollaxis(out_data, 0, out_data.ndim)
     else:
-        # Massage the output, to allow combining 3D and 4D images.
-        is_3D = [len(d.shape) == 3 for d in out_data]
-        is_4D = [len(d.shape) == 4 for d in out_data]
-        if np.any(is_3D) and np.any(is_4D):  # Convert all to 4D
-            out_data = [data if is_4D[di] else np.reshape(data, data.shape + (1,))
-                        for di, data in enumerate(out_data)]
         out_data = np.concatenate(out_data, axis=axis)
 
     return klass(out_data, affine, header)
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 749c52a390..80beec0f84 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -76,23 +76,9 @@ def test_concat():
                     #   but our efficient logic (where all images are
                     #   3D and the same size) fails, so we also
                     #   have to expect errors for those.
-                    expect_error = False
-                    if data0.ndim == 3 and data1.ndim == 4:
-                        expect_error = axis is None and data1.shape[3] != 1
-                        all_data = np.concatenate([data0[..., np.newaxis], data1],
-                                                  **np_concat_kwargs)
-                    elif data0.ndim == 4 and data1.ndim == 3:
-                        expect_error = axis is None and data0.shape[3] != 1
-                        all_data = np.concatenate([data0, data1[..., np.newaxis]],
-                                                  **np_concat_kwargs)
-                    elif data0.ndim == 4 and data1.ndim == 4:
-                        expect_error = axis is None and (data0.shape[3] != 1 or
-                                                         data1.shape[3] != 1)
-                        all_data = np.concatenate([data0, data1],
-                                                  **np_concat_kwargs)
-                    elif axis is None:  # 3D from here and below
-                        all_data = np.concatenate(
-                            [data0[..., np.newaxis], data1[..., np.newaxis]], 3)
+                    expect_error = data0.ndim != data1.ndim
+                    if axis is None:  # 3D from here and below
+                        all_data = np.concatenate([data0[..., np.newaxis],  data1[..., np.newaxis]],**np_concat_kwargs)
                     else:  # both 3D, appending on final axis
                         all_data = np.concatenate([data0, data1],
                                                   **np_concat_kwargs)

From c836dcefbcab49b76a8caacf13f588af08407b98 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 14:50:47 -0700
Subject: [PATCH 10/16] Fix broadcasting bug (img1.shape == (1,2,3),
 img2.shape=(1,2,1), axis=None would broadcast rather than error)

---
 nibabel/funcs.py | 8 +++++++-
 1 file changed, 7 insertions(+), 1 deletion(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index cd1d71f772..a2705069d2 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -122,10 +122,11 @@ def concat_images(images, check_affines=True, axis=None):
         if i == 0:  # first image, initialize data from loaded image
             affine = img.affine
             header = img.header
+            shape = img.shape
             klass = img.__class__
 
             if axis is None:  # collect images in output array for efficiency
-                out_shape = (n_imgs, ) + img.shape
+                out_shape = (n_imgs, ) + shape
                 out_data = np.empty(out_shape)
             else:  # collect images in list for use with np.concatenate
                 out_data = [None] * n_imgs
@@ -133,6 +134,11 @@ def concat_images(images, check_affines=True, axis=None):
         elif check_affines and not np.all(img.affine == affine):
             raise ValueError('Affines do not match')
 
+        elif axis is None and not np.array_equal(shape, img.shape):
+            # shape mismatch; numpy broadcasting can hide these.
+            raise ValueError("Image %d (shape=%s) does not match first image "
+                             " shape (%s)." % (i, shape, img.shape))
+
         out_data[i] = img.get_data()
 
         del img

From 49b353a4c5964774667ea6852b399474cb17b2fc Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 15:59:35 -0700
Subject: [PATCH 11/16] Similar bug in axis=int pathway, due to np.concatenate
 "smartness".

---
 nibabel/funcs.py | 23 +++++++++++++++++++----
 1 file changed, 19 insertions(+), 4 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index a2705069d2..e6f6d4e0b9 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -88,6 +88,19 @@ def squeeze_image(img):
                  img.extra)
 
 
+def _shape_equal_excluding(shape1, shape2, exclude_axes=None):
+    """ Helper function to compare two array shapes, excluding any
+    axis specified."""
+
+    if len(shape1) != len(shape2):
+        return False
+
+    idx_mask = np.ones((len(shape1),), dtype=bool)
+    idx_mask[exclude_axes] = False
+    return np.array_equal(np.asarray(shape1)[idx_mask],
+                          np.asarray(shape2)[idx_mask])
+
+
 def concat_images(images, check_affines=True, axis=None):
     ''' Concatenate images in list to single image, along specified dimension
 
@@ -134,10 +147,12 @@ def concat_images(images, check_affines=True, axis=None):
         elif check_affines and not np.all(img.affine == affine):
             raise ValueError('Affines do not match')
 
-        elif axis is None and not np.array_equal(shape, img.shape):
-            # shape mismatch; numpy broadcasting can hide these.
-            raise ValueError("Image %d (shape=%s) does not match first image "
-                             " shape (%s)." % (i, shape, img.shape))
+        elif ((axis is None and not np.array_equal(shape, img.shape)) or
+              (axis is not None and not _shape_equal_excluding(shape, img.shape,
+                                                               exclude_axes=[axis]))):
+            # shape mismatch; numpy broadcast / concatenate can hide these.
+            raise ValueError("Image #%d (shape=%s) does not match the first "
+                             "image shape (%s)." % (i, shape, img.shape))
 
         out_data[i] = img.get_data()
 

From 99f116828e933260282ce10662ea4324c55f7f77 Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 15:59:53 -0700
Subject: [PATCH 12/16] Test 2D - 5D; remove some tests to increase speed.

---
 nibabel/tests/test_funcs.py | 164 ++++++++++++++++++------------------
 1 file changed, 84 insertions(+), 80 deletions(-)

diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 80beec0f84..2f486f4036 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -34,88 +34,92 @@ def test_concat():
     assert_raises(ValueError, concat_images, [])
 
     # Build combinations of 3D, 4D w/size[3] == 1, and 4D w/size[3] == 3
-    all_shapes_3D = ((1, 2, 5), (7, 3, 1), (13, 11, 11), (0, 1, 1))
-    all_shapes_4D_unary = tuple((shape + (1,) for shape in all_shapes_3D))
-    all_shapes_4D_multi = tuple((shape + (3,) for shape in all_shapes_3D))
-    all_shapes = all_shapes_3D + all_shapes_4D_unary + all_shapes_4D_multi
+    all_shapes_5D = ((1, 4, 5, 3, 3),
+                     (7, 3, 1, 4, 5),
+                     (0, 2, 1, 4, 5))
 
     affine = np.eye(4)
-    # Loop over all possible combinations of images, in first and
-    #   second position.
-    for data0_shape in all_shapes:
-        data0_numel = np.asarray(data0_shape).prod()
-        data0 = np.arange(data0_numel).reshape(data0_shape)
-        img0_mem = Nifti1Image(data0, affine)
-
-        for data1_shape in all_shapes:
-            data1_numel = np.asarray(data1_shape).prod()
-            data1 = np.arange(data1_numel).reshape(data1_shape)
-            img1_mem = Nifti1Image(data1, affine)
-            img2_mem = Nifti1Image(data1, affine+1)  # bad affine
-            img3_mem = Nifti1Image(data1.T, affine)  # bad data shape
-
-            # Loop over every possible axis, including None (explicit and implied)
-            for axis in (list(range(-2, 3)) + [None, '__default__']):
-
-                # Allow testing default vs. passing explicit param
-                if axis == '__default__':
-                    np_concat_kwargs = dict(axis=-1)
-                    concat_imgs_kwargs = dict()
-                    axis = None  # Convert downstream
-                elif axis is None:
-                    np_concat_kwargs = dict(axis=-1)
-                    concat_imgs_kwargs = dict(axis=axis)
-                else:
-                    np_concat_kwargs = dict(axis=axis)
-                    concat_imgs_kwargs = dict(axis=axis)
-
-                # Create expected output
-                try:
-                    # Error will be thrown if the np.concatenate fails.
-                    #   However, when axis=None, the concatenate is possible
-                    #   but our efficient logic (where all images are
-                    #   3D and the same size) fails, so we also
-                    #   have to expect errors for those.
-                    expect_error = data0.ndim != data1.ndim
-                    if axis is None:  # 3D from here and below
-                        all_data = np.concatenate([data0[..., np.newaxis],  data1[..., np.newaxis]],**np_concat_kwargs)
-                    else:  # both 3D, appending on final axis
-                        all_data = np.concatenate([data0, data1],
-                                                  **np_concat_kwargs)
-                except ValueError:
-                    # Shapes are not combinable
-                    expect_error = True
-
-                # Check filenames and in-memory images work
-                with InTemporaryDirectory():
-                    # Try mem-based, file-based, and mixed
-                    imgs = [img0_mem, img1_mem, img2_mem, img3_mem]
-                    img_files = [_as_fname(img) for img in imgs]
-                    imgs_mixed = [imgs[0], img_files[1], imgs[2], img_files[3]]
-                    for img0, img1, img2, img3 in (imgs, img_files, imgs_mixed):
-                        try:
-                            all_imgs = concat_images([img0, img1],
-                                                     **concat_imgs_kwargs)
-                        except ValueError as ve:
-                            assert_true(expect_error, ve.message)
-                        else:
-                            assert_false(expect_error, "Expected a concatenation error, but got none.")
-                            assert_array_equal(all_imgs.get_data(), all_data)
-                            assert_array_equal(all_imgs.affine, affine)
-
-                        # check that not-matching affines raise error
-                        assert_raises(ValueError, concat_images, [img0, img2], **concat_imgs_kwargs)
-                        assert_raises(ValueError, concat_images, [img0, img3], **concat_imgs_kwargs)
-
-                        # except if check_affines is False
-                        try:
-                            all_imgs = concat_images([img0, img1], **concat_imgs_kwargs)
-                        except ValueError as ve:
-                            assert_true(expect_error, ve.message)
-                        else:
-                            assert_false(expect_error, "Expected a concatenation error, but got none.")
-                            assert_array_equal(all_imgs.get_data(), all_data)
-                            assert_array_equal(all_imgs.affine, affine)
+    for dim in range(2, 6):
+        all_shapes_ND = tuple((shape[:dim] for shape in all_shapes_5D))
+        all_shapes_N1D_unary = tuple((shape + (1,) for shape in all_shapes_ND))
+        all_shapes = all_shapes_ND + all_shapes_N1D_unary
+
+        # Loop over all possible combinations of images, in first and
+        #   second position.
+        for data0_shape in all_shapes:
+            data0_numel = np.asarray(data0_shape).prod()
+            data0 = np.arange(data0_numel).reshape(data0_shape)
+            img0_mem = Nifti1Image(data0, affine)
+
+            for data1_shape in all_shapes:
+                data1_numel = np.asarray(data1_shape).prod()
+                data1 = np.arange(data1_numel).reshape(data1_shape)
+                img1_mem = Nifti1Image(data1, affine)
+                img2_mem = Nifti1Image(data1, affine+1)  # bad affine
+
+                # Loop over every possible axis, including None (explicit and implied)
+                for axis in (list(range(-(dim-2), (dim-1))) + [None, '__default__']):
+
+                    # Allow testing default vs. passing explicit param
+                    if axis == '__default__':
+                        np_concat_kwargs = dict(axis=-1)
+                        concat_imgs_kwargs = dict()
+                        axis = None  # Convert downstream
+                    elif axis is None:
+                        np_concat_kwargs = dict(axis=-1)
+                        concat_imgs_kwargs = dict(axis=axis)
+                    else:
+                        np_concat_kwargs = dict(axis=axis)
+                        concat_imgs_kwargs = dict(axis=axis)
+
+                    # Create expected output
+                    try:
+                        # Error will be thrown if the np.concatenate fails.
+                        #   However, when axis=None, the concatenate is possible
+                        #   but our efficient logic (where all images are
+                        #   3D and the same size) fails, so we also
+                        #   have to expect errors for those.
+                        expect_error = data0.ndim != data1.ndim
+                        if axis is None:  # 3D from here and below
+                            all_data = np.concatenate([data0[..., np.newaxis],
+                                                       data1[..., np.newaxis]],
+                                                      **np_concat_kwargs)
+                        else:  # both 3D, appending on final axis
+                            all_data = np.concatenate([data0, data1],
+                                                      **np_concat_kwargs)
+                    except ValueError:
+                        # Shapes are not combinable
+                        expect_error = True
+
+                    # Check filenames and in-memory images work
+                    with InTemporaryDirectory():
+                        # Try mem-based, file-based, and mixed
+                        imgs = [img0_mem, img1_mem, img2_mem]
+                        img_files = [_as_fname(img) for img in imgs]
+                        imgs_mixed = [imgs[0], img_files[1], imgs[2]]
+                        for img0, img1, img2 in (imgs, img_files, imgs_mixed):
+                            try:
+                                all_imgs = concat_images([img0, img1],
+                                                         **concat_imgs_kwargs)
+                            except ValueError as ve:
+                                assert_true(expect_error, ve.message)
+                            else:
+                                assert_false(expect_error, "Expected a concatenation error, but got none.")
+                                assert_array_equal(all_imgs.get_data(), all_data)
+                                assert_array_equal(all_imgs.affine, affine)
+
+                            # check that not-matching affines raise error
+                            assert_raises(ValueError, concat_images, [img0, img2], **concat_imgs_kwargs)
+
+                            # except if check_affines is False
+                            try:
+                                all_imgs = concat_images([img0, img1], **concat_imgs_kwargs)
+                            except ValueError as ve:
+                                assert_true(expect_error, ve.message)
+                            else:
+                                assert_false(expect_error, "Expected a concatenation error, but got none.")
+                                assert_array_equal(all_imgs.get_data(), all_data)
+                                assert_array_equal(all_imgs.affine, affine)
 
 
 def test_closest_canonical():

From 0567ac3a8d83fc8607d093d08ee4d69bac8fdc7a Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 16:01:16 -0700
Subject: [PATCH 13/16] Convert exceptions to string.

---
 nibabel/tests/test_funcs.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 2f486f4036..18b7360737 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -102,7 +102,7 @@ def test_concat():
                                 all_imgs = concat_images([img0, img1],
                                                          **concat_imgs_kwargs)
                             except ValueError as ve:
-                                assert_true(expect_error, ve.message)
+                                assert_true(expect_error, str(ve))
                             else:
                                 assert_false(expect_error, "Expected a concatenation error, but got none.")
                                 assert_array_equal(all_imgs.get_data(), all_data)
@@ -115,7 +115,7 @@ def test_concat():
                             try:
                                 all_imgs = concat_images([img0, img1], **concat_imgs_kwargs)
                             except ValueError as ve:
-                                assert_true(expect_error, ve.message)
+                                assert_true(expect_error, str(ve))
                             else:
                                 assert_false(expect_error, "Expected a concatenation error, but got none.")
                                 assert_array_equal(all_imgs.get_data(), all_data)

From e9298cfc9aa6edf09dcb66a6fd8fbacff0f3271d Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 12 Mar 2015 16:18:18 -0700
Subject: [PATCH 14/16] Remove default argument.

---
 nibabel/funcs.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index e6f6d4e0b9..5afd5cb03c 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -88,7 +88,7 @@ def squeeze_image(img):
                  img.extra)
 
 
-def _shape_equal_excluding(shape1, shape2, exclude_axes=None):
+def _shape_equal_excluding(shape1, shape2, exclude_axes):
     """ Helper function to compare two array shapes, excluding any
     axis specified."""
 

From 84d990de3800bf7e6865904a0fb66f0d2007b00c Mon Sep 17 00:00:00 2001
From: Ben Cipollini <bcipolli@ucsd.edu>
Date: Thu, 26 Mar 2015 18:19:37 -0700
Subject: [PATCH 15/16] Small code review tweaks.

---
 nibabel/funcs.py            | 3 +--
 nibabel/tests/test_funcs.py | 1 -
 2 files changed, 1 insertion(+), 3 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 5afd5cb03c..4f4446b133 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -112,11 +112,10 @@ def concat_images(images, check_affines=True, axis=None):
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
     axis : None or int, optional
-        If None, concatenates on a new dimension.  This rrequires all images
+        If None, concatenates on a new dimension.  This requires all images
            to be the same shape).
         If not None, concatenates on the specified dimension.  This requires
            all images to be the same shape, except on the specified dimension.
-           For 4D images, axis must be between -2 and 3.
     Returns
     -------
     concat_img : ``SpatialImage``
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 18b7360737..5670451457 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -79,7 +79,6 @@ def test_concat():
                         #   but our efficient logic (where all images are
                         #   3D and the same size) fails, so we also
                         #   have to expect errors for those.
-                        expect_error = data0.ndim != data1.ndim
                         if axis is None:  # 3D from here and below
                             all_data = np.concatenate([data0[..., np.newaxis],
                                                        data1[..., np.newaxis]],

From 188c7ea90919e8e050b32b7c923ae79bf5d7a02a Mon Sep 17 00:00:00 2001
From: Matthew Brett <matthew.brett@gmail.com>
Date: Thu, 26 Mar 2015 22:08:17 -0700
Subject: [PATCH 16/16] RF: try doing image 1 concat stuff outside loop

Try doing the checks on the first image etc outside the body of the
loop.
---
 nibabel/funcs.py            | 84 ++++++++++++++++---------------------
 nibabel/tests/test_funcs.py |  1 +
 2 files changed, 38 insertions(+), 47 deletions(-)

diff --git a/nibabel/funcs.py b/nibabel/funcs.py
index 4f4446b133..fc69b6c780 100644
--- a/nibabel/funcs.py
+++ b/nibabel/funcs.py
@@ -88,19 +88,6 @@ def squeeze_image(img):
                  img.extra)
 
 
-def _shape_equal_excluding(shape1, shape2, exclude_axes):
-    """ Helper function to compare two array shapes, excluding any
-    axis specified."""
-
-    if len(shape1) != len(shape2):
-        return False
-
-    idx_mask = np.ones((len(shape1),), dtype=bool)
-    idx_mask[exclude_axes] = False
-    return np.array_equal(np.asarray(shape1)[idx_mask],
-                          np.asarray(shape2)[idx_mask])
-
-
 def concat_images(images, check_affines=True, axis=None):
     ''' Concatenate images in list to single image, along specified dimension
 
@@ -112,50 +99,53 @@ def concat_images(images, check_affines=True, axis=None):
        If True, then check that all the affines for `images` are nearly
        the same, raising a ``ValueError`` otherwise.  Default is True
     axis : None or int, optional
-        If None, concatenates on a new dimension.  This requires all images
-           to be the same shape).
-        If not None, concatenates on the specified dimension.  This requires
-           all images to be the same shape, except on the specified dimension.
+        If None, concatenates on a new dimension.  This requires all images to
+        be the same shape.  If not None, concatenates on the specified
+        dimension.  This requires all images to be the same shape, except on
+        the specified dimension.
     Returns
     -------
     concat_img : ``SpatialImage``
        New image resulting from concatenating `images` across last
        dimension
     '''
-
+    images = [load(img) if not hasattr(img, 'get_data')
+              else img for img in images]
     n_imgs = len(images)
     if n_imgs == 0:
         raise ValueError("Cannot concatenate an empty list of images.")
-
+    img0 = images[0]
+    affine = img0.affine
+    header = img0.header
+    klass = img0.__class__
+    shape0 = img0.shape
+    n_dim = len(shape0)
+    if axis is None:
+        # collect images in output array for efficiency
+        out_shape = (n_imgs, ) + shape0
+        out_data = np.empty(out_shape)
+    else:
+        # collect images in list for use with np.concatenate
+        out_data = [None] * n_imgs
+    # Get part of shape we need to check inside loop
+    idx_mask = np.ones((n_dim,), dtype=bool)
+    if axis is not None:
+        idx_mask[axis] = False
+    masked_shape = np.array(shape0)[idx_mask]
     for i, img in enumerate(images):
-        if not hasattr(img, 'get_data'):
-            img = load(img)
-
-        if i == 0:  # first image, initialize data from loaded image
-            affine = img.affine
-            header = img.header
-            shape = img.shape
-            klass = img.__class__
-
-            if axis is None:  # collect images in output array for efficiency
-                out_shape = (n_imgs, ) + shape
-                out_data = np.empty(out_shape)
-            else:  # collect images in list for use with np.concatenate
-                out_data = [None] * n_imgs
-
-        elif check_affines and not np.all(img.affine == affine):
-            raise ValueError('Affines do not match')
-
-        elif ((axis is None and not np.array_equal(shape, img.shape)) or
-              (axis is not None and not _shape_equal_excluding(shape, img.shape,
-                                                               exclude_axes=[axis]))):
-            # shape mismatch; numpy broadcast / concatenate can hide these.
-            raise ValueError("Image #%d (shape=%s) does not match the first "
-                             "image shape (%s)." % (i, shape, img.shape))
-
-        out_data[i] = img.get_data()
-
-        del img
+        if len(img.shape) != n_dim:
+            raise ValueError(
+                'Image {0} has {1} dimensions, image 0 has {2}'.format(
+                    i, len(img.shape), n_dim))
+        if not np.all(np.array(img.shape)[idx_mask] == masked_shape):
+            raise ValueError('shape {0} for image {1} not compatible with '
+                             'first image shape {2} with axis == {0}'.format(
+                                 img.shape, i, shape0, axis))
+        if check_affines and not np.all(img.affine == affine):
+            raise ValueError('Affine for image {0} does not match affine '
+                             'for first image'.format(i))
+        # Do not fill cache in image if it is empty
+        out_data[i] = img.get_data(caching='unchanged')
 
     if axis is None:
         out_data = np.rollaxis(out_data, 0, out_data.ndim)
diff --git a/nibabel/tests/test_funcs.py b/nibabel/tests/test_funcs.py
index 5670451457..20d11578b3 100644
--- a/nibabel/tests/test_funcs.py
+++ b/nibabel/tests/test_funcs.py
@@ -86,6 +86,7 @@ def test_concat():
                         else:  # both 3D, appending on final axis
                             all_data = np.concatenate([data0, data1],
                                                       **np_concat_kwargs)
+                        expect_error = False
                     except ValueError:
                         # Shapes are not combinable
                         expect_error = True