Code covered by the BSD License

# Method for quantitative absorption spectroscopy, version 2.1

### Tom O'Haver (view profile)

14 Aug 2006 (Updated )

A computational method for quantitative analysis by multiwavelength absorption spectroscopy

CalculateTfit(x,y,InstFunction,width,absorbance,InstWidth,straylight,noise,NoiseArray,IzeroShift)
```function CalculateTfit(x,y,InstFunction,width,absorbance,InstWidth,straylight,noise,NoiseArray,IzeroShift)
% Performs calculation and plots graph for DemoTFit
% example:
% CalculateTfit(x,y,InstFunction,width,absorbance,InstWidth,straylight,noise,NoiseArray,IzeroShift)
global z
global c
% Define frequency and wavelength coordinate vectors
ArrayLength=length(x);
j=[-ArrayLength/2:(ArrayLength/2)-1]';
f=(ArrayLength/2)-abs(j);

% Calculate noisy instrumentally-broadened transmission profile, yobsd,
% by convoluting true transmission spectrum y with instrument function InstFunction
% Note:  To model gaussian absorption, change 'lorentzian' to 'gaussian'
TrueSpectrum=lorentzian(x,(ArrayLength/2),width);
% Unbroadened transmission profile, plotted as green line on the graph.
y=10 .^ (absorbance .* (-TrueSpectrum));
fy=fft(y);
% Broadening function , plotted as magenta dashed line
InstFunction=gaussian(f,0,InstWidth);  % defined Gaussian instrument function centered on zero
fa=fft(InstFunction);
fy1=fy.*fa;             % Convolution by multiplying Fourier transforms
yobsd=real(ifft(fy1));  % and inverse transforming the result.
yo=yobsd./sum(InstFunction);
NoiseArray=randn(size(yo));
% ybosd is the broadened transmission profile with noise, plotted as red dots on the graph.
% Use only one of the following two lines, to model either photon or constant noise
yobsd=straylight+yo+(noise.*NoiseArray).*sqrt(yo).*1/InstWidth;  % Adds simulated photon noise
% yobsd=straylight+yo+(noise.*NoiseArray).*1/InstWidth;  % Adds simulated constant noise
yobsd=yobsd.*(1-straylight);
yobsd=yobsd.*(1+IzeroShift.*randn); % Random shifts in Izero
h=figure(1);

SingleWavelengthAbsorbance=-log10(yobsd(ArrayLength./2));  % Single-wavelength Beer's Law absornance
% Weighted regression
Background=ones(size(y));
weight=y;
WeightedRegression=([weight weight] .* [Background TrueSpectrum])\(-log10(yobsd) .* weight);

% TFit method
options = optimset('TolX',0.000001);
start=10; % Because of the very large dynamic range of absorbance, two start values are
if SingleWavelengthAbsorbance<1,start=1;,end  % used to prevent stalling on local optima.
lam=FMINSEARCH('fitM',start,options,yobsd,TrueSpectrum,InstFunction,straylight);
plot(x,real(yobsd),'r.',x,real(y),'g',x,real(c)*z,'b',x,gaussian(x,ArrayLength/2,InstWidth),'m:');
title([ 'Abs.Width = ' num2str(round(10*width)/10)  '   Inst.Width = ' num2str(round(10*InstWidth)/10) '    straylight= ' num2str(round(1000*mean(straylight))/10) '%   Noise = ' num2str(round(1000*noise)/10)  '%' ]);