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)

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