code: def _add_exg(raw, kind, head_pos, interp, n_jobs, random_state):
      assert isinstance(kind, str) and kind in ("ecg", "blink")
      _validate_type(raw, BaseRaw, "raw")
      _check_preload(raw, "Adding %s noise " % (kind,))
      rng = check_random_state(random_state)
      info, times, first_samp = raw.info, raw.times, raw.first_samp
      data = raw._data
      meg_picks = pick_types(info, meg=True, eeg=False, exclude=())
      meeg_picks = pick_types(info, meg=True, eeg=True, exclude=())
      R, r0 = fit_sphere_to_headshape(info, units="m", verbose=False)[
          :2
      ]
      bem = make_sphere_model(
          r0,
          head_radius=R,
          relative_radii=(0.97, 0.98, 0.99, 1.0),
          sigmas=(0.33, 1.0, 0.004, 0.33),
          verbose=False,
      )
      trans = None
      dev_head_ts, offsets = _check_head_pos(
          head_pos, info, first_samp, times
      )
      if kind == "blink":
          # place dipoles at 45 degree angles in z=0 plane
          exg_rr = np.array(
              [
                  [np.cos(np.pi / 3.0), np.sin(np.pi / 3.0), 0.0],
                  [-np.cos(np.pi / 3.0), np.sin(np.pi / 3), 0.0],
              ]
          )
          exg_rr /= np.sqrt(
              np.sum(exg_rr * exg_rr, axis=1, keepdims=True)
          )
          exg_rr *= 0.96 * R
          exg_rr += r0
          # oriented upward
          nn = np.array([[0.0, 0.0, 1.0], [0.0, 0.0, 1.0]])
          # Blink times drawn from an inhomogeneous poisson process
          # by 1) creating the rate and 2) pulling random numbers
          blink_rate = (
              1 + np.cos(2 * np.pi * 1.0 / 60.0 * times)
          ) / 2.0
          blink_rate *= 12.5 / 60.0
          blink_rate += 4.5 / 60.0
          blink_data = (
              rng.uniform(size=len(times)) < blink_rate / info["sfreq"]
          )
          blink_data = blink_data * (
              rng.uniform(size=len(times)) + 0.5
          )  # amps
          # Activation kernel is a simple hanning window
          blink_kernel = np.hanning(int(0.25 * info["sfreq"]))
          exg_data = (
              np.convolve(blink_data, blink_kernel, "same")[
                  np.newaxis, :
              ]
              * 1e-7
          )
          # Add rescaled noisy data to EOG ch
          ch = pick_types(info, meg=False, eeg=False, eog=True)
          picks = meeg_picks
          del blink_kernel, blink_rate, blink_data
      else:
          if len(meg_picks) == 0:
              raise RuntimeError(
                  "Can only add ECG artifacts if MEG data "
                  "channels are present"
              )
          exg_rr = np.array([[-R, 0, -3 * R]])
          max_beats = int(np.ceil(times[-1] * 80.0 / 60.0))
          # activation times with intervals drawn from a uniform distribution
          # based on activation rates between 40 and 80 beats per minute
          cardiac_idx = np.cumsum(
              rng.uniform(60.0 / 80.0, 60.0 / 40.0, max_beats)
              * info["sfreq"]
          ).astype(int)
          cardiac_idx = cardiac_idx[cardiac_idx < len(times)]
          cardiac_data = np.zeros(len(times))
          cardiac_data[cardiac_idx] = 1
          # kernel is the sum of three hanning windows
          cardiac_kernel = np.concatenate(
              [
                  2 * np.hanning(int(0.04 * info["sfreq"])),
                  -0.3 * np.hanning(int(0.05 * info["sfreq"])),
                  0.2 * np.hanning(int(0.26 * info["sfreq"])),
              ],
              axis=-1,
          )
          exg_data = (
              np.convolve(cardiac_data, cardiac_kernel, "same")[
                  np.newaxis, :
              ]
              * 15e-8
          )
          # Add rescaled noisy data to ECG ch
          ch = pick_types(info, meg=False, eeg=False, ecg=True)
          picks = meg_picks
          del cardiac_data, cardiac_kernel, max_beats, cardiac_idx
          nn = np.zeros_like(exg_rr)
          nn[:, 0] = 1  # arbitrarily rightward
      del meg_picks, meeg_picks
      noise = rng.standard_normal(exg_data.shape[1]) * 5e-6
      if len(ch) >= 1:
          ch = ch[-1]
          data[ch, :] = exg_data * 1e3 + noise
      else:
          ch = None
      src = setup_volume_source_space(
          pos=dict(rr=exg_rr, nn=nn), sphere_units="mm"
      )
      _log_ch("%s simulated and trace" % kind, info, ch)
      del ch, nn, noise
  
      used = np.zeros(len(raw.times), bool)
      get_fwd = _SimForwards(
          dev_head_ts,
          offsets,
          info,
          trans,
          src,
          bem,
          0.005,
          n_jobs,
          picks,
      )
      interper = _Interp2(offsets, get_fwd, interp)
      proc_lims = np.concatenate(
          [np.arange(0, len(used), 10000), [len(used)]]
      )
      for start, stop in zip(proc_lims[:-1], proc_lims[1:]):
          fwd, _ = interper.feed(stop - start)
          data[picks, start:stop] += np.einsum(
              "svt,vt->st", fwd, exg_data[:, start:stop]
          )
          assert not used[start:stop].any()
          used[start:stop] = True
      assert used.all()
  
imports: import os
  from collections.abc import Iterable
  from pathlib import Path
  import numpy
  from .._fiff.constants import FIFF
  from .._fiff.meas_info import Info
  from .._fiff.pick import pick_channels, pick_channels_forward, pick_info, pick_types
  from .._ola import _Interp2
  from ..bem import fit_sphere_to_headshape, make_sphere_model, read_bem_solution
  from ..chpi import _get_hpi_initial_fit, get_chpi_info, head_pos_to_trans_rot_t, read_head_pos
  from ..cov import Covariance, make_ad_hoc_cov, read_cov
  from ..event import _get_stim_channel
  from ..forward import _compute_forwards, _magnetic_dipole_field_vec, _merge_fwds, _prep_meg_channels, _prepare_for_forward, _stc_src_sel, _to_forward_dict, _transform_orig_meg_coils, convert_forward_solution, restrict_forward_to_stc
  from ..io import BaseRaw, RawArray
  from ..source_estimate import _BaseSourceEstimate
  from ..source_space._source_space import _ensure_src, _set_source_space_vertices, setup_volume_source_space
  from ..surface import _CheckInside
  from ..transforms import _get_trans, transform_surface_to
  from ..utils import _check_preload, _pl, _validate_type, check_random_state, logger, verbose
  from .source import SourceSimulator
label_desc: Use image node_id 8 for calling a global function with example usage: _add_exg(raw, kind, head_pos, interp, n_jobs, random_state) without return types