diff --git a/CHANGELOG.rst b/CHANGELOG.rst
index 6e1d981ce2f7ea6176dffd2e29a9dc70b7d0495a..ae4a1b6fe971ff6aa14b7a7edfa57440d5b328a3 100644
--- a/CHANGELOG.rst
+++ b/CHANGELOG.rst
@@ -31,6 +31,13 @@ Deprecated
   library on [PyPI](https://pypi.org/project/file-tree/).
 
 
+Fixed
+^^^^^
+
+* Fixed an edge-case in the :mod:`.gifti` module, where a surface with a
+  single triangle was being loaded incorrectly.
+
+
 
 3.5.3 (Tuesday 9th February 2021)
 ---------------------------------
diff --git a/fsl/data/gifti.py b/fsl/data/gifti.py
index d75c73a79a1f66d2fb8617848d48d4e11f7e9a36..077b7c4e14eea3c84064b45d9d69f9253a5ed742 100644
--- a/fsl/data/gifti.py
+++ b/fsl/data/gifti.py
@@ -229,7 +229,7 @@ def loadGiftiMesh(filename):
                          'at least one pointset array'.format(filename))
 
     vertices = [ps.data for ps in pointsets]
-    indices  = triangles[0].data
+    indices  = np.atleast_2d(triangles[0].data)
 
     if len(vdata) == 0: vdata = None
     else:               vdata = prepareGiftiVertexData(vdata, filename)
diff --git a/fsl/data/mesh.py b/fsl/data/mesh.py
index b03a628e2eec368f592a1d84a913b1302d9e2d00..41aee35425b16753f1a3e078b20d78d5c958894a 100644
--- a/fsl/data/mesh.py
+++ b/fsl/data/mesh.py
@@ -573,7 +573,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]
 
@@ -687,7 +687,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):
diff --git a/fsl/transform/nonlinear.py b/fsl/transform/nonlinear.py
index 74788a3dc0ebf8b2bcd98a053837858975b027f9..7a9b8e786bdf35a6cfb7d1eeca15d554477a1842 100644
--- a/fsl/transform/nonlinear.py
+++ b/fsl/transform/nonlinear.py
@@ -280,7 +280,7 @@ class DeformationField(NonLinearTransform):
         # Mask out the coordinates
         # that are out of bounds of
         # the deformation field
-        voxels  = np.round(voxels).astype(np.int)
+        voxels  = np.round(voxels).astype(np.int32)
         voxmask = (voxels >= [0, 0, 0]) & (voxels < self.shape[:3])
         voxmask = voxmask.all(axis=1)
         voxels  = voxels[voxmask]
@@ -508,9 +508,9 @@ class CoefficientField(NonLinearTransform):
         u = np.remainder(i, 1)
         v = np.remainder(j, 1)
         w = np.remainder(k, 1)
-        i = np.floor(i).astype(np.int)
-        j = np.floor(j).astype(np.int)
-        k = np.floor(k).astype(np.int)
+        i = np.floor(i).astype(np.int32)
+        j = np.floor(j).astype(np.int32)
+        k = np.floor(k).astype(np.int32)
 
         disps = np.zeros(coords.shape)
 
diff --git a/tests/test_atlases_query.py b/tests/test_atlases_query.py
index fedb0a00d9fa30405c50d2a1a880d78288e20f4e..c90a122e7e3a5b4f12e51e0e5df36428f89f2cb4 100644
--- a/tests/test_atlases_query.py
+++ b/tests/test_atlases_query.py
@@ -85,7 +85,7 @@ def _get_zero_mask(aimg):
         elif isinstance(aimg, fslatlases.ProbabilisticAtlas):
 
             # Keep memory usage down
-            zmask = np.ones(aimg.shape[:3], dtype=np.bool)
+            zmask = np.ones(aimg.shape[:3], dtype=bool)
             for vol in range(aimg.shape[-1]):
                 zmask = np.logical_and(zmask, aimg[..., vol] == 0)
 
diff --git a/tests/test_gifti.py b/tests/test_gifti.py
index 6d23e665e105f24e4756148ec896ea6c86436ff4..9266e86d98ee3804cdbaf3f7e3da0229270bad4f 100644
--- a/tests/test_gifti.py
+++ b/tests/test_gifti.py
@@ -381,3 +381,21 @@ def test_GiftiMesh_needsFixing():
         surf = gifti.GiftiMesh(fname, fixWinding=True)
 
         assert np.all(np.isclose(surf.indices, idxs_fixed))
+
+
+def test_loadGiftiMesh_onetriangle():
+    verts = np.array([[0, 0, 0], [1, 1, 1], [0, 1, 0]])
+    idxs  = np.array([[0, 1, 2]])
+    verts = nib.gifti.GiftiDataArray(verts, intent='NIFTI_INTENT_POINTSET')
+    idxs  = nib.gifti.GiftiDataArray(idxs,  intent='NIFTI_INTENT_TRIANGLE')
+    gimg  = nib.gifti.GiftiImage(darrays=[verts, idxs])
+
+    with tempdir():
+        fname = op.abspath('test.gii')
+        gimg.to_filename(fname)
+
+        gimg, tris, verts, _ = gifti.loadGiftiMesh('test.gii')
+        verts = verts[0]
+
+        assert verts.shape == (3, 3)
+        assert tris.shape  == (1, 3)
diff --git a/tests/test_image_resample.py b/tests/test_image_resample.py
index 6961e7aa71c63a09d21f93ca5e128d8622b678a9..ccc0b8935780c2793b060a5065133e792fb823b1 100644
--- a/tests/test_image_resample.py
+++ b/tests/test_image_resample.py
@@ -75,7 +75,7 @@ def test_resample(seed):
                    (np.isclose(np.modf(origtestcoords)[0], 0.5)))
             out = np.any(out, axis=1) | (resvals == 0)
 
-            origtestcoords = np.array(origtestcoords.round(), dtype=np.int)
+            origtestcoords = np.array(origtestcoords.round(), dtype=int)
 
             origtestcoords = origtestcoords[~out, :]
             restestcoords  = restestcoords[ ~out, :]
@@ -166,8 +166,8 @@ def test_resampleToPixdims():
 
     img          = fslimage.Image(make_random_image(dims=(10, 10, 10)))
     imglo, imghi = affine.axisBounds(img.shape, img.voxToWorldMat)
-    oldpix       = np.array(img.pixdim, dtype=np.float)
-    oldshape     = np.array(img.shape,  dtype=np.float)
+    oldpix       = np.array(img.pixdim, dtype=float)
+    oldshape     = np.array(img.shape,  dtype=float)
 
     for i, origin in it.product(range(25), ('centre', 'corner')):
 
@@ -215,16 +215,16 @@ def test_resampleToReference2():
     # More specific test - output
     # data gets transformed correctly
     # into reference space
-    img          = np.zeros((5, 5, 5), dtype=np.float)
+    img          = np.zeros((5, 5, 5), dtype=float)
     img[1, 1, 1] = 1
     img          = fslimage.Image(img)
 
     refv2w = affine.scaleOffsetXform([1, 1, 1], [-1, -1, -1])
-    ref    = np.zeros((5, 5, 5), dtype=np.float)
+    ref    = np.zeros((5, 5, 5), dtype=float)
     ref    = fslimage.Image(ref, xform=refv2w)
     res    = resample.resampleToReference(img, ref, order=0)
 
-    exp          = np.zeros((5, 5, 5), dtype=np.float)
+    exp          = np.zeros((5, 5, 5), dtype=float)
     exp[2, 2, 2] = 1
 
     assert np.all(np.isclose(res[0], exp))
diff --git a/tests/test_imagewrapper.py b/tests/test_imagewrapper.py
index 868245fefd24738197b47c2eab4e3214ffc28837..9de727fc98040f2ab7dcb21c17ff198bb65faf0a 100644
--- a/tests/test_imagewrapper.py
+++ b/tests/test_imagewrapper.py
@@ -1195,7 +1195,7 @@ def test_3D_indexing(shape=None, img=None):
 
     assert type(img[0, 0, 0]) == np.float64
 
-    mask1 = np.zeros(shape, dtype=np.bool)
+    mask1 = np.zeros(shape, dtype=bool)
 
     mask1[0, 0, 0] = True
     mask1[1, 0, 0] = True
@@ -1252,7 +1252,7 @@ def test_3D_4D_indexing():
     assert type(img[0, 0, 0, 0]) == np.float64
     assert type(img[0, 0, 0, :]) == np.float64
 
-    mask = np.zeros(padShape, dtype=np.bool)
+    mask = np.zeros(padShape, dtype=bool)
     mask[0, 0, 0, 0] = True
 
     assert type(img[mask])      == np.ndarray
@@ -1282,13 +1282,13 @@ def test_3D_len_one_indexing(shape=None, img=None):
     assert type(img[0, 0])    == np.ndarray
     assert type(img[0, 0, 0]) == np.float64
 
-    mask = np.zeros(shape[:2], dtype=np.bool)
+    mask = np.zeros(shape[:2], dtype=bool)
     mask[0, 0] = True
 
     assert type(img[mask])      == np.ndarray
     assert      img[mask].shape == (1, )
 
-    mask = np.zeros(shape, dtype=np.bool)
+    mask = np.zeros(shape, dtype=bool)
     mask[0, 0, 0] = True
 
     assert type(img[mask])      == np.ndarray
@@ -1352,7 +1352,7 @@ def test_4D_indexing(shape=None, img=None):
 
     assert type(img[0, 0, 0, 0]) == np.float64
 
-    mask1 = np.zeros(shape, dtype=np.bool)
+    mask1 = np.zeros(shape, dtype=bool)
 
     mask1[0, 0, 0, 0] = True
     mask1[1, 0, 0, 0] = True
diff --git a/tests/test_mesh.py b/tests/test_mesh.py
index f5ae5ead39be0dc3621a847c2083636cbf23e6a7..3a3ff94b87dd8ea5ad55423290050f2b304d5f08 100644
--- a/tests/test_mesh.py
+++ b/tests/test_mesh.py
@@ -234,6 +234,13 @@ def test_normals():
         -fnormals, fslmesh.calcFaceNormals(verts, triangles_cw)))
     assert np.all(np.isclose(
         fnormals, fslmesh.calcFaceNormals(verts, triangles_ccw)))
+
+    # Make sure result is (1, 3) for input of (1, 3)
+    onetri = np.atleast_2d(triangles_ccw[0, :])
+    result = fslmesh.calcFaceNormals(verts, onetri)
+    assert result.shape == (1, 3)
+    assert np.all(np.isclose(fnormals[0, :], result))
+
     assert np.all(np.isclose(
         -vnormals, fslmesh.calcVertexNormals(verts, triangles_cw, -fnormals)))
     assert np.all(np.isclose(
diff --git a/tests/test_transform/test_affine.py b/tests/test_transform/test_affine.py
index 469c6c13d98950e379a9a042467b0a4f4d517dcb..f1dceb61312c7f05f97c9d871dc8a4a0ae95f492 100644
--- a/tests/test_transform/test_affine.py
+++ b/tests/test_transform/test_affine.py
@@ -395,7 +395,7 @@ def test_transform():
         mask = np.array([[1, 0, 0, 1],
                          [0, 1, 0, 1],
                          [0, 0, 1, 1],
-                         [0, 0, 0, 1]], dtype=np.bool)
+                         [0, 0, 0, 1]], dtype=bool)
 
         return np.all((xform != 0) == mask)
 
diff --git a/tests/test_transform/test_nonlinear.py b/tests/test_transform/test_nonlinear.py
index 8bfb0c3aba0d7d0746cf48d56b19631dfc1c2536..1cb8cf8f7df8e2fea674961d00b898f72b1858f7 100644
--- a/tests/test_transform/test_nonlinear.py
+++ b/tests/test_transform/test_nonlinear.py
@@ -161,7 +161,7 @@ def test_convertDeformationSpace():
         refcoords = [np.random.randint(0, basefield.shape[0], 5),
                      np.random.randint(0, basefield.shape[1], 5),
                      np.random.randint(0, basefield.shape[2], 5)]
-        refcoords = np.array(refcoords, dtype=np.int).T
+        refcoords = np.array(refcoords, dtype=int).T
         refcoords = affine.transform(refcoords, ref.voxToScaledVoxMat)
         srccoords = basefield.transform(refcoords)