diff --git a/iodata/test/test_molden.py b/iodata/test/test_molden.py
index 774da705..f7e3f371 100644
--- a/iodata/test/test_molden.py
+++ b/iodata/test/test_molden.py
@@ -260,30 +260,34 @@ def test_load_molden_nh3_molden_cart():
     assert_allclose(charges, molden_charges, atol=1.0e-3)
 
 
-def test_load_molden_cfour():
-    # The file tested here is created with CFOUR 2.1.
-    file_list = [
-        "h_sonly_sph_cfour.molden",
-        "h_ponly_sph_cfour.molden",
-        "h_donly_sph_cfour.molden",
-        "h_fonly_sph_cfour.molden",
-        "h_gonly_sph_cfour.molden",
-        "h_sonly_cart_cfour.molden",
-        "h_ponly_cart_cfour.molden",
-        "h_donly_cart_cfour.molden",
-        "h_fonly_cart_cfour.molden",
-        "h_gonly_cart_cfour.molden",
-        "h2o_ccpvdz_cfour.molden",
-    ]
-
-    for i in file_list:
-        with as_file(files("iodata.test.data").joinpath(i)) as fn_molden:
-            print(str(fn_molden))
+# The file tested here is created with CFOUR 2.1.
+@pytest.mark.parametrize(
+    ("path", "should_warn"),
+    [
+        ("h_sonly_sph_cfour.molden", False),
+        ("h_ponly_sph_cfour.molden", False),
+        ("h_donly_sph_cfour.molden", True),
+        ("h_fonly_sph_cfour.molden", True),
+        ("h_gonly_sph_cfour.molden", True),
+        ("h_sonly_cart_cfour.molden", False),
+        ("h_ponly_cart_cfour.molden", False),
+        ("h_donly_cart_cfour.molden", True),
+        ("h_fonly_cart_cfour.molden", True),
+        ("h_gonly_cart_cfour.molden", True),
+        ("h2o_ccpvdz_cfour.molden", True),
+    ],
+)
+def test_load_molden_cfour(path, should_warn):
+    with as_file(files("iodata.test.data").joinpath(path)) as fn_molden:
+        if should_warn:
+            with pytest.warns(FileFormatWarning):
+                mol = load_one(str(fn_molden))
+        else:
             mol = load_one(str(fn_molden))
-            # Check normalization
-            olp = compute_overlap(mol.obasis, mol.atcoords)
-            check_orthonormal(mol.mo.coeffsa, olp)
-            check_orthonormal(mol.mo.coeffsb, olp)
+    # Check normalization
+    olp = compute_overlap(mol.obasis, mol.atcoords)
+    check_orthonormal(mol.mo.coeffsa, olp)
+    check_orthonormal(mol.mo.coeffsb, olp)
 
 
 def test_load_molden_nh3_orca():