docs: def _create_convolution_functions(self, matched_channels, data, 
          free_parameters, fixed_parameters=None):
          """
          Pre-create binning matrix factories. The following options need to be
          followed on a per-matched channel basis.
  
          Options here are:
  
          1) If we have no redshift or resolution parameters to solve for, then
             we can just produce a matrix that will be multiplied in.
             Inputs: none
             Outputs: matrix
  
          If fast_binning is turned on:
  
          2) Provide the wavelengths, redshift, resolution (single value) and
             flux. Returns the convolved, interpolated flux.
             Inputs: flux, wavelengths (e.g., redshift), resolution
             Outputs: normalised flux.
  
          If fast_binning is turned off:
  
          3) If we *just* have redshift parameters to solve for, then we can
             produce a _BoxFactory that will take a redshift z and return
             the matrix
             [This function will have a LRU cacher and should be removed after]
             Inputs: redshift
             Outputs: matrix
  
          4) If we have redshift parameters and resolution parameters, then we
             can produce a _BlurryBoxFactory that will take a redshift z and
             Resolution (or resolution coefficients!) and return a binning matrix
             [This function will have a LRU cacher and should be removed after]
             Inputs: redshift, Resolution parameter(s)
             Outputs: matrix
  
          5) If we *just* have resolution parameters to solve for, then we 
             can produce a _BlurryBoxFactory as well, because by default z=0.
             [This function will have a LRU cacher and should be removed after]
             Inputs: resolution parameter(s)
             Outputs: matrix
  
  
          Because these options require different input/outputs, and the LRU
          cachers need to remain wrapped around simple functions, we may need a
  
          Consistent lambda function:
          lambda(obs_wavelength, obs_flux, z=0, *R)
          """
  
          fast_binning = self._configuration.get("settings", {}).get(
              "fast_binning", 1)
  
          logger.info("Creating convolution functions (fast_binning = {})".format(
              fast_binning))
          logger.info("Free parameters: {}".format(free_parameters))
          logger.info("Fixed parameters: {}".format(fixed_parameters))
  
          if fixed_parameters is None:
              fixed_parameters = {}
  
          convolution_functions = []
          for channel, spectrum in zip(matched_channels, data):
              if channel is None:
                  convolution_functions.append(None)
                  continue
  
              # Any redshift or resolution parameters?
              redshift = "z" in free_parameters \
                  or "z_{}".format(channel) in free_parameters
              resolution = "resolution_{}".format(channel) in free_parameters
  
              # Option 1.
              # Create static binning matrices for each channel.
              if not redshift and not resolution:
  
                  if fast_binning:
                      logger.info("Doing static interpolation for channel {}"\
                          .format(channel))
                      convolution_function = lambda w, f, z, *a: \
                          np.interp(w, generate.wavelengths[-1], f,
                              left=np.nan, right=np.nan)
  
  
                  else:
  
                      logger.info("Creating static matrix for channel {0} because "\
                          "no redshift or resolution parameters were found".format(
                              channel))
  
                      # Is there any z?
                      z = fixed_parameters.get(
                          "z_{}".format(channel), fixed_parameters.get("z", 0))
  
                      # Create the binning matrix based on the globally-scoped array
                      # of wavelengths.
                      matrix = specutils.sample.resample(
                          generate.wavelengths[-1] * (1 + z),
                          spectrum.disp)
  
                      # Wrap in a lambda function to be consistent with other options.
                      convolution_function = lambda w, f, *a: f * matrix
  
              else:
  
                  # If fast_binning is turned on (default):
                  if fast_binning:
                      
                      # Option 2: Convolve with single kernel & interpolate.
                      # [TODO] should w.mean()/R be squared?
                      # px_sigma ~= ((w.mean()/R) / 2.35482)/np.diff(w).mean()
                      # px_scale = ((w.mean()/R) / 2.35)/np.diff(w).mean()
                      # px_scale = R_scale / R
  
                      logger.info("Creating simple convolution & interpolating "\
                          "function for channel {0}".format(channel))
  
                      # [TODO] Account for the existing spectral resolution of the
                      # grid.
                      if resolution:
                          R_scale = spectrum.disp.mean() \
                              / (2.35482 * np.diff(spectrum.disp).mean())
                          """
                          convolution_function = lambda w, f, z, R, *a: \
                              np.interp(w, generate.wavelengths[-1] * (1 + z),
                                  gaussian_filter1d(f, max(0, R_scale/R),
                                      mode="constant", cval=np.nan),
                                  left=np.nan, right=np.nan)
                          """
                          def convolution_function(w, f, z, R, *a):
                              if R > 0:
                                  _ = gaussian_filter1d(f, R_scale/R)
                              else:
                                  _ = f
                              return np.interp(w, generate.wavelengths[-1] * (1 + z),
                                  _, left=np.nan, right=np.nan)
  
                      else:
                          convolution_function = lambda w, f, z, *a: \
                              np.interp(w, generate.wavelengths[-1] * (1 + z), f,
                                  left=np.nan, right=np.nan)
  
                  else:
                      if redshift and not resolution:
                          # Option 3: Produce a _BoxFactory
                          logger.info("Producing a Box Factory for convolution "\
                              "in channel {}".format(channel))
  
                          matrix = specutils.sample._BoxFactory(
                              spectrum.disp, generate.wavelengths[-1])
  
                          # Wrap in a lambda function to be consistent.
                          convolution_function = lambda w, f, z, *a: f * matrix(z)
  
                      else:
                          # Could be redshift and resolution, or just resolution.
                          # Options 4 and 5: Produce a _BlurryBoxFactory
                          logger.info("Producing a Blurry Box Factory for "\
                              "convolution in channel {}".format(channel))
  
                          matrix = specutils.sample._BlurryBoxFactory(
                              spectrum.disp, generate.wavelengths[-1])
  
                          # Wrap in a lambda function to be consistent.
                          convolution_function \
                              = lambda w, f, z, R, *a: f * matrix(R, z)
  
              # Append this channel's convolution function.
              convolution_functions.append(convolution_function)
  
          # Put the convolution functions into the global scope.
          generate.binning_matrices.append(convolution_functions)
          return True