code: def _compute_corrs(
      X,
      y,
      smin,
      smax,
      n_jobs=None,
      fit_intercept=False,
      edge_correction=True,
  ):
      """Compute auto- and cross-correlations."""
      if fit_intercept:
          # We could do this in the Fourier domain, too, but it should
          # be a bit cleaner numerically to do it here.
          X_offset = np.mean(X, axis=0)
          y_offset = np.mean(y, axis=0)
          if X.ndim == 3:
              X_offset = X_offset.mean(axis=0)
              y_offset = np.mean(y_offset, axis=0)
          X = X - X_offset
          y = y - y_offset
      else:
          X_offset = y_offset = 0.0
      if X.ndim == 2:
          assert y.ndim == 2
          X = X[:, np.newaxis, :]
          y = y[:, np.newaxis, :]
      assert X.shape[:2] == y.shape[:2]
      len_trf = smax - smin
      len_x, n_epochs, n_ch_x = X.shape
      len_y, n_epcohs, n_ch_y = y.shape
      assert len_x == len_y
  
      n_fft = next_fast_len(2 * X.shape[0] - 1)
  
      _, cuda_dict = _setup_cuda_fft_multiply_repeated(
          n_jobs, [1.0], n_fft, "correlation calculations"
      )
      del n_jobs  # only used to set as CUDA
  
      # create our Toeplitz indexer
      ij = np.empty((len_trf, len_trf), int)
      for ii in range(len_trf):
          ij[ii, ii:] = np.arange(len_trf - ii)
          x = np.arange(n_fft - 1, n_fft - len_trf + ii, -1)
          ij[ii + 1 :, ii] = x
  
      x_xt = np.zeros([n_ch_x * len_trf] * 2)
      x_y = np.zeros((len_trf, n_ch_x, n_ch_y), order="F")
      n = n_epochs * (n_ch_x * (n_ch_x + 1) // 2 + n_ch_x)
      logger.info("Fitting %d epochs, %d channels" % (n_epochs, n_ch_x))
      pb = ProgressBar(n, mesg="Sample")
      count = 0
      pb.update(count)
      for ei in range(n_epochs):
          this_X = X[:, ei, :]
          # XXX maybe this is what we should parallelize over CPUs at some point
          X_fft = cuda_dict["rfft"](this_X, n=n_fft, axis=0)
          X_fft_conj = X_fft.conj()
          y_fft = cuda_dict["rfft"](y[:, ei, :], n=n_fft, axis=0)
  
          for ch0 in range(n_ch_x):
              for oi, ch1 in enumerate(range(ch0, n_ch_x)):
                  this_result = cuda_dict["irfft"](
                      X_fft[:, ch0] * X_fft_conj[:, ch1],
                      n=n_fft,
                      axis=0,
                  )
                  # Our autocorrelation structure is a Toeplitz matrix, but
                  # it's faster to create the Toeplitz ourselves than use
                  # linalg.toeplitz.
                  this_result = this_result[ij]
                  # However, we need to adjust for coeffs that are cut off,
                  # i.e. the non-zero delays should not have the same AC value
                  # as the zero-delay ones (because they actually have fewer
                  # coefficients).
                  #
                  # These adjustments also follow a Toeplitz structure, so we
                  # construct a matrix of what has been left off, compute their
                  # inner products, and remove them.
                  if edge_correction:
                      _edge_correct(
                          this_result, this_X, smax, smin, ch0, ch1
                      )
  
                  # Store the results in our output matrix
                  x_xt[
                      ch0 * len_trf : (ch0 + 1) * len_trf,
                      ch1 * len_trf : (ch1 + 1) * len_trf,
                  ] += this_result
                  if ch0 != ch1:
                      x_xt[
                          ch1 * len_trf : (ch1 + 1) * len_trf,
                          ch0 * len_trf : (ch0 + 1) * len_trf,
                      ] += this_result.T
                  count += 1
                  pb.update(count)
  
              # compute the crosscorrelations
              cc_temp = cuda_dict["irfft"](
                  y_fft * X_fft_conj[:, slice(ch0, ch0 + 1)],
                  n=n_fft,
                  axis=0,
              )
              if smin < 0 and smax >= 0:
                  x_y[:-smin, ch0] += cc_temp[smin:]
                  x_y[len_trf - smax :, ch0] += cc_temp[:smax]
              else:
                  x_y[:, ch0] += cc_temp[smin:smax]
              count += 1
              pb.update(count)
  
      x_y = np.reshape(x_y, (n_ch_x * len_trf, n_ch_y), order="F")
      return x_xt, x_y, n_ch_x, X_offset, y_offset
  
imports: import numpy
  from scipy import linalg
  from scipy.signal import fftconvolve
  from scipy.sparse.csgraph import laplacian
  from ..cuda import _setup_cuda_fft_multiply_repeated
  from ..filter import next_fast_len
  from ..fixes import jit
  from ..utils import ProgressBar, _check_option, _validate_type, logger, warn
  from .base import BaseEstimator
label_desc: Use image node_id 1 for calling a global function with example usage: _compute_corrs(X, y, smin, smax, n_jobs, fit_intercept, edge_correction) and returns: x_xt,  x_y,  n_ch_x,  X_offset,  y_offset