diff --git a/examples/ilastik3d.py b/examples/ilastik3d.py
index 8733d0e6c687cfcad83a623dc3b27f1db7e5074d..0e3814757becfcec3493069c7c999aada2558888 100644
--- a/examples/ilastik3d.py
+++ b/examples/ilastik3d.py
@@ -65,4 +65,6 @@ if __name__ == '__main__':
 
     resp_models = requests.get(uri + 'models')
 
+    raise Exception('x, y axes are currently flipped in pixelmap')
+
     server_process.terminate()
\ No newline at end of file
diff --git a/model_server/image.py b/model_server/image.py
index 5d16d980f26cf037bcab7e72379a722eba5ec0f0..8cbd37c1ed265260a356db75b7b095d4ea0dbcf1 100644
--- a/model_server/image.py
+++ b/model_server/image.py
@@ -81,14 +81,14 @@ class CziImageFileAccessor(GenericImageFileAccessor):
         if (sd.get('S') and (sd['S'] > 1)) or (sd.get('T') and (sd['T'] > 1)):
             raise DataShapeError(f'Cannot handle image with multiple positions or time points: {sd}')
 
-        idx = {k: sd[k] for k in ['X', 'Y', 'C', 'Z']}
-        xycz = np.moveaxis(
+        idx = {k: sd[k] for k in ['Y', 'X', 'C', 'Z']}
+        yxcz = np.moveaxis(
             cf.asarray(),
             [cf.axes.index(ch) for ch in idx],
             [0, 1, 2, 3]
         )
 
-        self._data = self.conform_data(xycz.reshape(xycz.shape[0:4]))
+        self._data = self.conform_data(yxcz.reshape(yxcz.shape[0:4]))
 
     def __del__(self):
         self.czifile.close()
@@ -98,7 +98,7 @@ def write_accessor_data_to_file(fpath: Path, accessor: GenericImageDataAccessor)
     try:
         zcxy = np.moveaxis(
             accessor.data,
-            [3, 2, 0, 1],
+            [3, 2, 1, 0],
             [0, 1, 2, 3]
         )
         tifffile.imwrite(fpath, zcxy, imagej=True)
diff --git a/tests/test_image.py b/tests/test_image.py
index 83d25ba1488080aa657c7bfe212c84a87387a1f0..00a746fe4d467eef398abcf0f4eebaa6726ef241 100644
--- a/tests/test_image.py
+++ b/tests/test_image.py
@@ -26,9 +26,9 @@ class TestCziImageFileAccess(unittest.TestCase):
         nc = 4
         nz = 11
         c = 3
-        cf = InMemoryDataAccessor(np.random.rand(w, h, nc, nz))
+        cf = InMemoryDataAccessor(np.random.rand(h, w, nc, nz))
         sc = cf.get_one_channel_data(c)
-        self.assertEqual(sc.shape, (w, h, 1, nz))
+        self.assertEqual(sc.shape, (h, w, 1, nz))
 
     def test_write_single_channel_tif(self):
         ch = 4
@@ -58,4 +58,23 @@ class TestCziImageFileAccess(unittest.TestCase):
     def test_conform_data_longer_than_xycz(self):
         data = np.random.rand(256, 512, 12, 8, 3)
         with self.assertRaises(DataShapeError):
-            acc = InMemoryDataAccessor(data)
\ No newline at end of file
+            acc = InMemoryDataAccessor(data)
+
+
+    def test_write_multichannel_image_preserves_axes(self):
+        h = 256
+        w = 512
+        c = 3
+        nz = 10
+
+        yxcz = (2**8 * np.random.rand(h, w, c, nz)).astype('uint8')
+        acc = InMemoryDataAccessor(yxcz)
+        self.assertTrue(
+            write_accessor_data_to_file(
+                output_path / f'rand3d.tif',
+                acc
+            )
+        )
+        # need to sort out x,y flipping since np convention yxcz flips axes in 3d tif
+        self.assertEqual(acc.shape_dict['X'], w, acc.shape_dict)
+        self.assertEqual(acc.shape_dict['Y'], h, acc.shape_dict)
\ No newline at end of file