@@ -0,0 +1 @@
+affine-RAS.afni

@@ -0,0 +1 @@
+affine-RAS.afni

@@ -0,0 +1,2 @@
+# 3dvolreg matrices (DICOM-to-DICOM, row-by-row):
+0.999999	-0.000999999	-0.001	-4	0.00140494	0.621609	0.783327	-2	-0.000161717	-0.783327	0.62161	-1

@@ -57,15 +57,19 @@ def test_affines_save(image_orientation):
 itk = nbl.load(os.path.join(data_path, 'affine-%s-itk.tfm' % image_orientation),
 fmt='itk')
 fsl = np.loadtxt(os.path.join(data_path, 'affine-%s.fsl' % image_orientation))
+afni = np.loadtxt(os.path.join(data_path, 'affine-%s.afni' % image_orientation))
 
 with InTemporaryDirectory():
 xfm.to_filename('M.tfm', fmt='itk')
 xfm.to_filename('M.fsl', fmt='fsl')
+xfm.to_filename('M.afni', fmt='afni')
 
 nb_itk = nbl.load('M.tfm', fmt='itk')
 nb_fsl = np.loadtxt('M.fsl')
+nb_afni = np.loadtxt('M.afni')
 
 assert_equal(itk, nb_itk)
 assert_almost_equal(fsl, nb_fsl)
+assert_almost_equal(afni, nb_afni)
 
 # Create version not aligned to canonical

@@ -204,7 +204,7 @@ def _to_hdf5(self, x5_root):
 
 def to_filename(self, filename, fmt='X5', moving=None):
 """Store the transform in BIDS-Transforms HDF5 file format (.x5)."""
-if fmt.lower() in ['itk', 'ants', 'elastix', 'nifty']:
+if fmt.lower() in ['itk', 'ants', 'elastix']:
 with open(filename, 'w') as f:
 f.write('#Insight Transform File V1.0\n')
 
@@ -229,26 +229,25 @@ def to_filename(self, filename, fmt='X5', moving=None):
 T = self.matrix.copy()
 pre = LPS
 post = LPS
-if any(obliquity(self.reference.affine) * 180 / pi > OBLIQUITY_THRESHOLD_DEG):
+if obliquity(self.reference.affine).min() * 180 / pi > OBLIQUITY_THRESHOLD_DEG:
 print('Reference affine axes are oblique.')
 M = self.reference.affine
 A = shape_zoom_affine(self.reference.shape,
-voxel_sizes(M), x_flip=True)
-pre = M.dot(np.linalg.inv(A))
+voxel_sizes(M), x_flip=False, y_flip=False)
+pre = M.dot(np.linalg.inv(A)).dot(LPS)
 
 if not moving:
 moving = self.reference
 
-if moving and any(obliquity(moving.affine) * 180 / pi > OBLIQUITY_THRESHOLD_DEG):
+if moving and obliquity(moving.affine).min() * 180 / pi > OBLIQUITY_THRESHOLD_DEG:
 print('Moving affine axes are oblique.')
-M = moving.affine
-A = shape_zoom_affine(moving.shape,
-voxel_sizes(M), x_flip=True)
-post = M.dot(np.linalg.inv(A))
+M2 = moving.affine
+A2 = shape_zoom_affine(moving.shape,
+voxel_sizes(M2), x_flip=True, y_flip=True)
+post = A2.dot(np.linalg.inv(M2))
 
 # swapaxes is necessary, as axis 0 encodes series of transforms
-T = np.swapaxes(post.dot(self.matrix.copy().dot(pre)), 0, 1)
-parameters = np.swapaxes(post.dot(T.dot(pre)), 0, 1)
+parameters = np.swapaxes(post.dot(self.matrix.copy().dot(pre)), 0, 1)
 parameters = parameters[:, :3, :].reshape((T.shape[0], -1))
 np.savetxt(filename, parameters, delimiter='\t', header="""\
 3dvolreg matrices (DICOM-to-DICOM, row-by-row):""", fmt='%g')

@@ -1456,7 +1456,7 @@ def finite_range(arr, check_nan=False):
 return np.nanmin(mins), np.nanmax(maxes)
 
 
-def shape_zoom_affine(shape, zooms, x_flip=True):
+def shape_zoom_affine(shape, zooms, x_flip=True, y_flip=False):
 ''' Get affine implied by given shape and zooms
 
 We get the translations from the center of the image (implied by
@@ -1471,6 +1471,9 @@ def shape_zoom_affine(shape, zooms, x_flip=True):
 x_flip : {True, False}
 whether to flip the X row of the affine. Corresponds to
 radiological storage on disk.
+y_flip : {False, True}
+whether to flip the Y row of the affine. Corresponds to
+DICOM storage on disk when x_flip is also True.
 
 Returns
 -------
@@ -1510,6 +1513,9 @@ def shape_zoom_affine(shape, zooms, x_flip=True):
 zooms = full_zooms
 if x_flip:
 zooms[0] *= -1
+
+if y_flip:
+zooms[1] *= -1
 # Get translations from center of image
 origin = (shape - 1) / 2.0
 aff = np.eye(4)