code: def _hrf_kernel(hrf_model, tr, oversampling=50, fir_delays=None):
      """Return the list of matching kernels \
      given the specification of the hemodynamic model and time parameters.
  
      Parameters
      ----------
      hrf_model : string, function, list of functions, or None,
          HRF model to be used.
  
      tr : float
          the repetition time in seconds
  
      oversampling : int, default=50
          Temporal oversampling factor to have a smooth hrf.
  
      fir_delays : 1D-array-like, optional
          List of delays (in scans) for finite impulse response models.
  
      Returns
      -------
      hkernel : list of arrays
          Samples of the hrf (the number depends on the hrf_model used).
  
      """
      acceptable_hrfs = [
          "spm",
          "spm + derivative",
          "spm + derivative + dispersion",
          "fir",
          "glover",
          "glover + derivative",
          "glover + derivative + dispersion",
          None,
      ]
      error_msg = (
          "Could not process custom HRF model provided. "
          "Please refer to the related documentation."
      )
      if hrf_model == "spm":
          hkernel = [spm_hrf(tr, oversampling)]
      elif hrf_model == "spm + derivative":
          hkernel = [
              spm_hrf(tr, oversampling),
              spm_time_derivative(tr, oversampling),
          ]
      elif hrf_model == "spm + derivative + dispersion":
          hkernel = [
              spm_hrf(tr, oversampling),
              spm_time_derivative(tr, oversampling),
              spm_dispersion_derivative(tr, oversampling),
          ]
      elif hrf_model == "glover":
          hkernel = [glover_hrf(tr, oversampling)]
      elif hrf_model == "glover + derivative":
          hkernel = [
              glover_hrf(tr, oversampling),
              glover_time_derivative(tr, oversampling),
          ]
      elif hrf_model == "glover + derivative + dispersion":
          hkernel = [
              glover_hrf(tr, oversampling),
              glover_time_derivative(tr, oversampling),
              glover_dispersion_derivative(tr, oversampling),
          ]
      elif hrf_model == "fir":
          hkernel = [
              np.hstack(
                  (
                      np.zeros((f) * oversampling),
                      np.ones(oversampling) * 1.0 / oversampling,
                  )
              )
              for f in fir_delays
          ]
      elif callable(hrf_model):
          try:
              hkernel = [hrf_model(tr, oversampling)]
          except TypeError:
              raise ValueError(error_msg)
      elif isinstance(hrf_model, Iterable) and all(
          callable(_) for _ in hrf_model
      ):
          try:
              hkernel = [model(tr, oversampling) for model in hrf_model]
          except TypeError:
              raise ValueError(error_msg)
      elif hrf_model is None:
          hkernel = [np.hstack((1, np.zeros(oversampling - 1)))]
      else:
          raise ValueError(
              f'"{hrf_model}" is not a known hrf model. '
              "Use either a custom model or "
              f"one of {acceptable_hrfs}"
          )
      return hkernel
  
imports: import warnings
  from collections.abc import Iterable
  import numpy
  from scipy.stats import gamma
  from nilearn._utils import fill_doc
label_desc: Use image node_id 17 for calling a global function with example usage: _hrf_kernel(hrf_model, tr, oversampling, fir_delays) and returns: hkernel