code: def test_simulate_raw_bem(raw_data):
      """Test simulation of raw data with BEM."""
      pytest.importorskip("nibabel")
      raw, src_ss, stc, trans, sphere = raw_data
      src = setup_source_space(
          "sample", "oct1", subjects_dir=subjects_dir
      )
      for s in src:
          s["nuse"] = 3
          s["vertno"] = src[1]["vertno"][:3]
          s["inuse"].fill(0)
          s["inuse"][s["vertno"]] = 1
      # use different / more complete STC here
      vertices = [s["vertno"] for s in src]
      stc = SourceEstimate(
          np.eye(sum(len(v) for v in vertices)),
          vertices,
          0,
          1.0 / raw.info["sfreq"],
      )
      stcs = [stc] * 15
      raw_sim_sph = simulate_raw(raw.info, stcs, trans, src, sphere)
      raw_sim_bem = simulate_raw(raw.info, stcs, trans, src, bem_fname)
      # some components (especially radial) might not match that well,
      # so just make sure that most components have high correlation
      assert_array_equal(raw_sim_sph.ch_names, raw_sim_bem.ch_names)
      picks = pick_types(raw.info, meg=True, eeg=True)
      n_ch = len(picks)
      corr = np.corrcoef(raw_sim_sph[picks][0], raw_sim_bem[picks][0])
      assert_array_equal(corr.shape, (2 * n_ch, 2 * n_ch))
      med_corr = np.median(np.diag(corr[:n_ch, -n_ch:]))
      assert med_corr > 0.65
      # do some round-trip localization
      for s in src:
          transform_surface_to(s, "head", trans)
      locs = np.concatenate([s["rr"][s["vertno"]] for s in src])
      tmax = (len(locs) - 1) / raw.info["sfreq"]
      cov = make_ad_hoc_cov(raw.info)
      # The tolerance for the BEM is surprisingly high (28) but I get the same
      # result when using MNE-C and Xfit, even when using a proper 5120 BEM :(
      for use_raw, bem, tol in (
          (raw_sim_sph, sphere, 2),
          (raw_sim_bem, bem_fname, 31),
      ):
          events = find_events(use_raw, "STI 014")
          assert len(locs) == 6
          evoked = Epochs(
              use_raw, events, 1, 0, tmax, baseline=None
          ).average()
          assert len(evoked.times) == len(locs)
          fits = fit_dipole(evoked, cov, bem, trans, min_dist=1.0)[
              0
          ].pos
          diffs = np.sqrt(np.sum((locs - fits) ** 2, axis=-1)) * 1000
          med_diff = np.median(diffs)
          assert med_diff < tol, "%s: %s" % (bem, med_diff)
      # also test event timings with SourceSimulator
      first_samp = raw.first_samp
      events = find_events(raw, initial_event=True, verbose=False)
      evt_times = events[:, 0]
      assert len(events) == 3
      labels_sim = [[], [], []]  # random l+r hemisphere points
      labels_sim[0] = Label([src_ss[0]["vertno"][1]], hemi="lh")
      labels_sim[1] = Label([src_ss[0]["vertno"][4]], hemi="lh")
      labels_sim[2] = Label([src_ss[1]["vertno"][2]], hemi="rh")
      wf_sim = np.array([2, 1, 0])
      for this_fs in (0, first_samp):
          ss = SourceSimulator(
              src_ss, 1.0 / raw.info["sfreq"], first_samp=this_fs
          )
          for i in range(3):
              ss.add_data(labels_sim[i], wf_sim, events[np.newaxis, i])
          assert ss.n_times == evt_times[-1] + len(wf_sim) - this_fs
      raw_sim = simulate_raw(
          raw.info, ss, src=src_ss, bem=bem_fname, first_samp=first_samp
      )
      data = raw_sim.get_data()
      amp0 = data[:, evt_times - first_samp].max()
      amp1 = data[:, evt_times + 1 - first_samp].max()
      amp2 = data[:, evt_times + 2 - first_samp].max()
      assert_allclose(amp0 / amp1, wf_sim[0] / wf_sim[1], rtol=1e-5)
      assert amp2 == 0
      assert raw_sim.n_times == ss.n_times
  
imports: from copy import deepcopy
  from pathlib import Path
  import numpy
  import pytest
  from numpy.testing import assert_allclose, assert_array_equal
  from mne import Epochs, SourceEstimate, VolSourceEstimate, convert_forward_solution, create_info, find_events, fit_dipole, make_ad_hoc_cov, make_bem_solution, make_forward_solution, make_sphere_model, pick_types, read_bem_solution, read_cov, read_source_spaces, read_trans, setup_source_space, setup_volume_source_space, transform_surface_to
  from mne._fiff.constants import FIFF
  from mne.bem import _surfaces_to_bem
  from mne.chpi import compute_chpi_amplitudes, compute_chpi_locs, compute_head_pos, get_chpi_info, read_head_pos
  from mne.datasets import testing
  from mne.io import RawArray, read_raw_fif
  from mne.label import Label
  from mne.simulation import add_chpi, add_ecg, add_eog, add_noise, simulate_raw, simulate_sparse_stc
  from mne.simulation.source import SourceSimulator
  from mne.source_space._source_space import _compare_source_spaces
  from mne.surface import _get_ico_surface
  from mne.tests.test_chpi import _assert_quats
  from mne.utils import catch_logging
label_desc: Use image node_id 8 for calling a global function with example usage: test_simulate_raw_bem(raw_data) without return types