code: def _execute_csd_function(
      X,
      times,
      frequencies,
      csd_function,
      params,
      n_fft,
      ch_names=None,
      projs=None,
      n_jobs=None,
      *,
      verbose=None,
  ):
      """Estimate cross-spectral density with a given function.
  
      This function will apply the given CSD function in parallel across epochs.
  
      Parameters
      ----------
      X : array-like, shape (n_epochs, n_channels, n_times)
          The time series data consisting of n_epochs separate observations
          of signals with n_channels time-series of length n_times.
      times : float
          Timestamps for each sample.
      frequencies : list of float
          The frequencies of interest for which the CSD is going to be computed.
      csd_function : function
          Function that performs the actual CSD computation
      params : list
          List of parameters to pass the CSD function.
      n_fft : int
          Number of FFT points. This is stored in the CSD object.
      ch_names : list of str | None
          A name for each time series. If ``None`` (the default), the series will
          be named 'SERIES###'.
      projs : list of Projection | None
          List of projectors to store in the CSD object. Defaults to ``None``,
          which means the projectors defined in the Epochs object will be copied.
      %(n_jobs)s
      %(verbose)s
  
      Returns
      -------
      csd : instance of CrossSpectralDensity
          The computed cross-spectral density.
      """
      n_epochs, n_channels, _ = X.shape
  
      logger.info("Computing cross-spectral density from epochs...")
  
      n_freqs = len(frequencies)
      csds_mean = np.zeros(
          (n_channels * (n_channels + 1) // 2, n_freqs),
          dtype=np.complex128,
      )
  
      # Prepare the function that does the actual CSD computation for parallel
      # execution.
      parallel, my_csd, n_jobs = parallel_func(
          csd_function, n_jobs, verbose=verbose
      )
  
      # Compute CSD for each trial
      n_blocks = int(np.ceil(n_epochs / float(n_jobs)))
      for i in ProgressBar(range(n_blocks), mesg="CSD epoch blocks"):
          epoch_block = X[i * n_jobs : (i + 1) * n_jobs]
          csds = parallel(
              my_csd(this_epoch, *params) for this_epoch in epoch_block
          )
  
          # Add CSD matrices in-place
          csds_mean += np.sum(csds, axis=0)
  
      csds_mean /= n_epochs
      logger.info("[done]")
  
      if ch_names is None:
          ch_names = ["SERIES%03d" % (i + 1) for i in range(n_channels)]
  
      return CrossSpectralDensity(
          csds_mean,
          ch_names=ch_names,
          tmin=times[0],
          tmax=times[-1],
          frequencies=frequencies,
          n_fft=n_fft,
          projs=projs,
      )
  
imports: import copy
  import numbers
  import numpy
  from scipy.fft import rfftfreq
  from .._fiff.pick import _picks_to_idx, pick_channels
  from ..parallel import parallel_func
  from ..time_frequency.multitaper import _compute_mt_params, _csd_from_mt, _mt_spectra, _psd_from_mt_adaptive
  from ..utils import ProgressBar, _check_fname, _import_h5io_funcs, _validate_type, copy_function_doc_to_method_doc, logger, verbose, warn
  from ..viz.misc import plot_csd
  from .tfr import EpochsTFR, _cwt_array, _get_nfft, morlet
label_desc: Use image node_id 14 for calling a global function with example usage: _execute_csd_function(X, times, frequencies, csd_function, params, n_fft, ch_names, projs, n_jobs) and returns: CrossSpectralDensity