diff --git a/tests/models/test_phase_thick_3d.py b/tests/models/test_phase_thick_3d.py
index f2baadc5..52c42f63 100644
--- a/tests/models/test_phase_thick_3d.py
+++ b/tests/models/test_phase_thick_3d.py
@@ -1,7 +1,9 @@
+import pytest
 from waveorder.models import phase_thick_3d
 
 
-def test_calculate_transfer_function():
+@pytest.mark.parametrize("axial_flip", (True, False))
+def test_calculate_transfer_function(axial_flip):
     z_padding = 5
     H_re, H_im = phase_thick_3d.calculate_transfer_function(
         zyx_shape=(20, 100, 101),
@@ -12,6 +14,7 @@ def test_calculate_transfer_function():
         index_of_refraction_media=1.0,
         numerical_aperture_illumination=0.45,
         numerical_aperture_detection=0.55,
+        axial_flip=axial_flip
     )
 
     assert H_re.shape == (20 + 2 * z_padding, 100, 101)
diff --git a/waveorder/models/inplane_anisotropic_thin_pol3d.py b/waveorder/models/inplane_anisotropic_thin_pol3d.py
index b828ca6e..6ded20cf 100644
--- a/waveorder/models/inplane_anisotropic_thin_pol3d.py
+++ b/waveorder/models/inplane_anisotropic_thin_pol3d.py
@@ -50,7 +50,7 @@ def apply_transfer_function(
 def apply_inverse_transfer_function(
     czyx_data,
     intensity_to_stokes_matrix,
-    illumination_wavelength,  # TOOD: MOVE THIS PARAM TO OTF? (leaky param)
+    wavelength_illumination,  # TOOD: MOVE THIS PARAM TO OTF? (leaky param)
     cyx_no_sample_data=None,  # if not None, use this data for background correction
     project_stokes_to_2d=False,
     remove_estimated_background=False,  # if True estimate background from czyx_data and remove it
@@ -97,7 +97,7 @@ def apply_inverse_transfer_function(
         *background_corrected_stokes
     )
 
-    # Return retardance in distance units (matching illumination_wavelength)
-    retardance = adr_parameters[0] * illumination_wavelength / (2 * np.pi)
+    # Return retardance in distance units (matching wavelength_illumination)
+    retardance = adr_parameters[0] * wavelength_illumination / (2 * np.pi)
 
     return retardance, adr_parameters[1], adr_parameters[2], adr_parameters[3]
diff --git a/waveorder/models/phase_thick_3d.py b/waveorder/models/phase_thick_3d.py
index 39f29f66..eb186acb 100644
--- a/waveorder/models/phase_thick_3d.py
+++ b/waveorder/models/phase_thick_3d.py
@@ -38,6 +38,7 @@ def calculate_transfer_function(
     index_of_refraction_media,
     numerical_aperture_illumination,
     numerical_aperture_detection,
+    axial_flip=False,
 ):
     radial_frequencies = util.generate_radial_frequencies(
         zyx_shape[1:], yx_pixel_size
@@ -46,6 +47,8 @@ def calculate_transfer_function(
     z_position_list = torch.fft.ifftshift(
         (torch.arange(z_total) - z_total // 2) * z_pixel_size
     )
+    if axial_flip:
+        z_position_list = torch.flip(z_position_list, dims=(0,))
 
     ill_pupil = optics.generate_pupil(
         radial_frequencies,
@@ -143,7 +146,7 @@ def apply_inverse_transfer_function(
     imaginary_potential_transfer_function,
     z_padding,
     z_pixel_size,  # TODO: MOVE THIS PARAM TO OTF? (leaky param)
-    illumination_wavelength,  # TOOD: MOVE THIS PARAM TO OTF? (leaky param)
+    wavelength_illumination,  # TOOD: MOVE THIS PARAM TO OTF? (leaky param)
     absorption_ratio=0.0,
     method="Tikhonov",
     reg_re=1e-3,
@@ -194,4 +197,4 @@ def apply_inverse_transfer_function(
     if z_padding != 0:
         f_real = f_real[z_padding:-z_padding]
 
-    return f_real * z_pixel_size / 4 / np.pi * illumination_wavelength
+    return f_real * z_pixel_size / 4 / np.pi * wavelength_illumination
diff --git a/waveorder/optics.py b/waveorder/optics.py
index dc692164..43260605 100644
--- a/waveorder/optics.py
+++ b/waveorder/optics.py
@@ -326,7 +326,7 @@ def generate_propagation_kernel(
 
 
 def generate_greens_function_z(
-    radial_frequencies, pupil_support, illumination_wavelength, z_position_list
+    radial_frequencies, pupil_support, wavelength_illumination, z_position_list
 ):
     """
 
@@ -340,7 +340,7 @@ def generate_greens_function_z(
         pupil_support   : torch.tensor
                         the array that defines the support of the pupil function with the size of (Y, X)
 
-        illumination_wavelength       : float
+        wavelength_illumination       : float
                         wavelength of the light in the immersion media
 
         z_position_list      : torch.tensor or list
@@ -354,9 +354,9 @@ def generate_greens_function_z(
     """
 
     oblique_factor = (
-        (1 - illumination_wavelength**2 * radial_frequencies**2)
+        (1 - wavelength_illumination**2 * radial_frequencies**2)
         * pupil_support
-    ) ** (1 / 2) / illumination_wavelength
+    ) ** (1 / 2) / wavelength_illumination
 
     greens_function_z = (
         -1j