Field Spectroscopy Facility Post Processing Toolbox
by Iain Robinson
A toolbox for importing and processing optical spectra acquired with portable spectroradiometers.
%CONVOLVE Convolves a spectra to a set of spectral response functions
% C = CONVOLVE(S, B) convolves the specturm or spectra in the structure
% array B to the set of spectral response functions ("bands") in the
% array B.
% Known issues:
% 1. Cannot be applied where data are overlapping or not sorted.
% FSF Post Processing Toolbox
% Field Spectroscopy Facility, Natural Environment Research Council
% Author: Iain Robinson
% Contact: email@example.com
% Requirements: MATLAB R2009a or later
% the Band class from the FSF Post Processing Toolbox
% Revision: 1.3.3
% Date: 2011-09-04
% Convolve a spectrum to a band set.
function output = convolve(input, bands)
% Copy the input spectrum to the output spectrum. This ensures that
% header data is preserved, as well as other fields in the
% structure such as the spectrum's name and datetime.
output(n) = input(n);
% The wavelength scale of the input spectrum.
inputWavelengths = input(n).wavelength;
inputValues = input(n).data;
% Check that there are no repeated wavelength values.
error('The input spectrum''s wavelength data contains at least one repeated value:\n\t%d\nThis may indicate that there are overlapped data. The convolution cannot be applied to overlapped data.', mode(inputWavelengths))
% Check that input wavelengths are sorted in increasing order.
error('The input spectrum''s wavelength data must be sorted in increasing order of wavelength.');
% The function should also check that the input spectrum is defined at
% equally spaced wavelength invervals. As yet, it does not do this.
% Interpolate bands to the wavelength scale of the spectrum (if
bandWavelengths = bands(1).wavelength;
% Check whether interpolation is required...
if isequal(size(inputWavelengths), size(bandWavelengths)) && isequal(inputWavelengths, bandWavelengths)
interpolatedBands = bands;
% Interpolate bands to the wavelength scale of the input spectra.
% The interpolate function is defined in the Band class.
interpolatedBands = bands.interpolate(inputWavelengths);
% The output wavelength is the band's centre wavelength.
% If the band's centre wavelength is not defined the outputs
% are simply numbered. A better approach may be to estimate the
% band's centre wavelenght from the band data.
outputWavelengths(i) = i;
outputWavelengths(i) = bands(i).centre;
% The output value is the integral of the product
% of the band with the input spectrum.
outputValues(i) = trapz(inputWavelengths, inputValues .* interpolatedBands(i).response);
% Copy the input spectrum to the output spectrum, then overwrite
% wavelength and data with the new values. (Copying the structure
% ensures that metadata fields such as dsatetime are preserved.)
output(n).wavelength = outputWavelengths';
output(n).data = outputValues';