diff --git a/mne/bem.py b/mne/bem.py index d361272fd49..2b0d9524c2f 100644 --- a/mne/bem.py +++ b/mne/bem.py @@ -2393,56 +2393,75 @@ def make_scalp_surfaces( if mri == "T1.mgz": mri = mri if (subj_path / "mri" / mri).exists() else "T1" - logger.info("1. Creating a dense scalp tessellation with mkheadsurf...") + threshold = _ensure_int(threshold, "threshold") - def check_seghead(surf_path=subj_path / "surf"): - surf = None - for k in ["lh.seghead", "lh.smseghead"]: - this_surf = surf_path / k - if this_surf.exists(): - surf = this_surf - break - return surf + # Check for existing files + _check_file(subj_path / "mri" / "seghead.mgz", overwrite) - my_seghead = check_seghead() - threshold = _ensure_int(threshold, "threshold") - if my_seghead is None: - this_env = deepcopy(os.environ) - this_env["SUBJECTS_DIR"] = str(subjects_dir) - this_env["SUBJECT"] = subject - this_env["subjdir"] = str(subj_path) - if "FREESURFER_HOME" not in this_env: - raise RuntimeError( - "The FreeSurfer environment needs to be set up to use " - "make_scalp_surfaces to create the outer skin surface " - "lh.seghead" - ) - run_subprocess( - [ - "mkheadsurf", - "-subjid", - subject, - "-srcvol", - mri, - "-thresh1", - str(threshold), - "-thresh2", - str(threshold), - ], - env=this_env, - ) + seghead_path = subj_path / "surf" / "lh.seghead" + _check_file(seghead_path, overwrite) - surf = check_seghead() - if surf is None: - raise RuntimeError("mkheadsurf did not produce the standard output file.") + smseghead_path = subj_path / "surf" / "lh.smseghead" + _check_file(smseghead_path, overwrite) bem_dir = subjects_dir / subject / "bem" - if not bem_dir.is_dir(): - os.mkdir(bem_dir) fname_template = bem_dir / (f"{subject}-head-{{}}.fif") dense_fname = str(fname_template).format("dense") - logger.info(f"2. Creating {dense_fname} ...") _check_file(dense_fname, overwrite) + + for level in _tri_levels.items(): + dec_fname = str(fname_template).format(level) + if overwrite: + os.remove(dec_fname) + else: + if no_decimate: + if os.path.exists(dec_fname): + raise OSError( + f"Trying to generate new scalp surfaces" + f"but {dec_fname} already exists." + f"To avoid mixing different scalp surface solutions," + f"delete this file or use overwrite to automatically delete it." + ) + else: + _check_file(dec_fname, overwrite) + + logger.info("1. Creating a dense scalp tessellation with mkheadsurf...") + + this_env = deepcopy(os.environ) + this_env["SUBJECTS_DIR"] = str(subjects_dir) + this_env["SUBJECT"] = subject + this_env["subjdir"] = str(subj_path) + if "FREESURFER_HOME" not in this_env: + raise RuntimeError( + "The FreeSurfer environment needs to be set up to use " + "make_scalp_surfaces to create the outer skin surface " + "lh.seghead" + ) + run_subprocess( + [ + "mkheadsurf", + "-subjid", + subject, + "-srcvol", + mri, + "-thresh1", + str(threshold), + "-thresh2", + str(threshold), + ], + env=this_env, + ) + if os.path.exists(seghead_path): + surf = seghead_path + elif os.path.exists(smseghead_path): + surf = smseghead_path + else: + raise ValueError("mkheadsurf did not produce the standard output file.") + + logger.info(f"2. Creating {dense_fname} ...") + if not bem_dir.is_dir(): + os.mkdir(bem_dir) + # Helpful message if we get a topology error msg = ( "\n\nConsider using pymeshfix directly to fix the mesh, or --force " @@ -2452,6 +2471,7 @@ def check_seghead(surf_path=subj_path / "surf"): [surf], [FIFF.FIFFV_BEM_SURF_ID_HEAD], [1], incomplete=incomplete, extra=msg )[0] write_bem_surfaces(dense_fname, surf, overwrite=overwrite) + if os.getenv("_MNE_TESTING_SCALP", "false") == "true": tris = [len(surf["tris"])] # don't actually decimate for ii, (level, n_tri) in enumerate(_tri_levels.items(), 3): @@ -2467,7 +2487,6 @@ def check_seghead(surf_path=subj_path / "surf"): ) dec_fname = str(fname_template).format(level) logger.info(f"{ii}.2 Creating {dec_fname}") - _check_file(dec_fname, overwrite) dec_surf = _surfaces_to_bem( [dict(rr=points, tris=tris)], [FIFF.FIFFV_BEM_SURF_ID_HEAD], diff --git a/mne/commands/tests/test_commands.py b/mne/commands/tests/test_commands.py index ae885d3c21d..dc53e045dd9 100644 --- a/mne/commands/tests/test_commands.py +++ b/mne/commands/tests/test_commands.py @@ -12,7 +12,6 @@ import pytest from numpy.testing import assert_allclose, assert_equal -import mne from mne import ( concatenate_raws, read_bem_solution, @@ -167,53 +166,74 @@ def test_kit2fiff(): check_usage(mne_kit2fiff, force_help=True) -@pytest.mark.slowtest @pytest.mark.ultraslowtest @testing.requires_testing_data +@requires_freesurfer("mkheadsurf") def test_make_scalp_surfaces(tmp_path, monkeypatch): """Test mne make_scalp_surfaces.""" pytest.importorskip("nibabel") pytest.importorskip("pyvista") check_usage(mne_make_scalp_surfaces) - has = "SUBJECTS_DIR" in os.environ - # Copy necessary files to avoid FreeSurfer call + tempdir = str(tmp_path) - surf_path = op.join(subjects_dir, "sample", "surf") - surf_path_new = op.join(tempdir, "sample", "surf") - os.mkdir(op.join(tempdir, "sample")) - os.mkdir(surf_path_new) - subj_dir = op.join(tempdir, "sample", "bem") - os.mkdir(subj_dir) + t1_path = op.join(subjects_dir, "sample", "mri", "T1.mgz") + t1_path_new = op.join(tempdir, "sample", "mri", "T1.mgz") - cmd = ("-s", "sample", "--subjects-dir", tempdir) - monkeypatch.setattr( - mne.bem, - "decimate_surface", - lambda points, triangles, n_triangles: (points, triangles), - ) - dense_fname = op.join(subj_dir, "sample-head-dense.fif") - medium_fname = op.join(subj_dir, "sample-head-medium.fif") + headseg_path = op.join(tempdir, "sample", "mri", "seghead.mgz") + surf_path = op.join(tempdir, "sample", "surf", "lh.seghead") + dense_fname = op.join(tempdir, "sample", "bem", "sample-head-dense.fif") + medium_fname = op.join(tempdir, "sample", "bem", "sample-head-medium.fif") + sparse_fname = op.join(tempdir, "sample", "bem", "sample-head-sparse.fif") + + os.makedirs(op.join(tempdir, "sample", "mri"), exist_ok=True) + os.makedirs(op.join(tempdir, "sample", "surf"), exist_ok=True) + + shutil.copy(t1_path, t1_path_new) + cmd = ("-s", "sample", "--subjects-dir", tempdir, "--no-decimate") with ArgvSetter(cmd, disable_stdout=False, disable_stderr=False): monkeypatch.delenv("FREESURFER_HOME", raising=False) with pytest.raises(RuntimeError, match="The FreeSurfer environ"): mne_make_scalp_surfaces.run() - shutil.copy(op.join(surf_path, "lh.seghead"), surf_path_new) + monkeypatch.setenv("FREESURFER_HOME", tempdir) mne_make_scalp_surfaces.run() + + assert op.isfile(headseg_path) + assert op.isfile(surf_path) + assert op.isfile(dense_fname) + assert not op.isfile(medium_fname) + assert not op.isfile(sparse_fname) + + cmd = ("-s", "sample", "--subjects-dir", tempdir) + with ArgvSetter(cmd, disable_stdout=False, disable_stderr=False): + with pytest.raises(RuntimeError, match="use --overwrite to overwrite it"): + mne_make_scalp_surfaces.run() + + cmd = ("-s", "sample", "--subjects-dir", tempdir, "--overwrite") + with ArgvSetter(cmd, disable_stdout=False, disable_stderr=False): + mne_make_scalp_surfaces.run() + assert op.isfile(headseg_path) + assert op.isfile(surf_path) assert op.isfile(dense_fname) assert op.isfile(medium_fname) - with pytest.raises(OSError, match="overwrite"): + assert op.isfile(sparse_fname) + + cmd = ("-s", "sample", "--subjects-dir", tempdir, "--no-decimate", "--overwrite") + with ArgvSetter(cmd, disable_stdout=False, disable_stderr=False): + mne_make_scalp_surfaces.run() + assert op.isfile(headseg_path) + assert op.isfile(surf_path) + assert op.isfile(dense_fname) + assert not op.isfile(medium_fname) + assert not op.isfile(sparse_fname) + + os.remove(headseg_path) + os.remove(surf_path) + os.remove(dense_fname) + cmd = ("-s", "sample", "--subjects-dir", tempdir, "--no-decimate") + with ArgvSetter(cmd, disable_stdout=False, disable_stderr=False): + with pytest.raises(RuntimeError, match="Trying to generate new scalp surfaces"): mne_make_scalp_surfaces.run() - # actually check the outputs - head_py = read_bem_surfaces(dense_fname) - assert_equal(len(head_py), 1) - head_py = head_py[0] - head_c = read_bem_surfaces( - op.join(subjects_dir, "sample", "bem", "sample-head-dense.fif") - )[0] - assert_allclose(head_py["rr"], head_c["rr"]) - if not has: - assert "SUBJECTS_DIR" not in os.environ @pytest.mark.slowtest diff --git a/mne/tests/test_bem.py b/mne/tests/test_bem.py index 9aa8a85fd4d..55a8ab07cb3 100644 --- a/mne/tests/test_bem.py +++ b/mne/tests/test_bem.py @@ -43,7 +43,12 @@ from mne.io import read_info from mne.surface import _get_ico_surface, read_surface from mne.transforms import translation -from mne.utils import _record_warnings, catch_logging, check_version +from mne.utils import ( + _record_warnings, + catch_logging, + check_version, + requires_freesurfer, +) fname_raw = Path(__file__).parents[1] / "io" / "tests" / "data" / "test_raw.fif" subjects_dir = testing.data_path(download=False) / "subjects" @@ -519,11 +524,27 @@ def test_io_head_bem(tmp_path): assert np.allclose(head["tris"], head_defect["tris"]) -@pytest.mark.slowtest # ~4 s locally +@pytest.mark.slowtest +@requires_freesurfer("mkheadsurf") +@testing.requires_testing_data def test_make_scalp_surfaces_topology(tmp_path, monkeypatch): - """Test topology checks for make_scalp_surfaces.""" + """Test make_scalp_surfaces and topology checks.""" pytest.importorskip("pyvista") pytest.importorskip("nibabel") + + # tests on 'sample' + subject = "sample" + subjects_dir = testing.data_path(download=False) + with pytest.raises(OSError, match="use --overwrite to overwrite it"): + make_scalp_surfaces( + subject, subjects_dir, force=False, verbose=True, overwrite=False + ) + + make_scalp_surfaces( + subject, subjects_dir, force=False, verbose=True, overwrite=True + ) + + # tests on custom surface subjects_dir = tmp_path subject = "test" surf_dir = subjects_dir / subject / "surf" @@ -551,6 +572,7 @@ def _decimate_surface(points, triangles, n_triangles): make_scalp_surfaces( subject, subjects_dir, force=False, verbose=True, overwrite=True ) + bem_dir = subjects_dir / subject / "bem" sparse_path = bem_dir / f"{subject}-head-sparse.fif" assert not sparse_path.is_file()