TFit3Redraw(x,y,InstFunction,linewidth,absorbance,InstWidth,StrayLight,noise,NoiseArray,separation,width,IzeroShift)

function TFit3Redraw(x,y,InstFunction,linewidth,absorbance,InstWidth,StrayLight,noise,NoiseArray,separation,width,IzeroShift)
% Performs calculation and plots graph for TFit3
% Example:
% TFit3Redraw(x,y,InstFunction,linewidth,absorbance,InstWidth,StrayLight,noise,NoiseArray,separation,width,IzeroShift)

% Define frequency and wavelength coordinate vectors
global z
ArrayLength=length(x);
j=[-ArrayLength/2:(ArrayLength/2)-1]';
f=(ArrayLength/2)-abs(j);

% InstWidth is the instrumental broadening function
InstFunction=gaussian(f,0,InstWidth);  % InstFunction = Gaussian instrument function centered on zero
fa=fft(InstFunction);  % Fourier transform of same, for later convolution purposes
% Simulate noisy instrumentally-broadened transmission profile, yobsd,
% by convoluting true transmission spectrum y with instrument function InstFunction
% and adding noise.  
position1=ArrayLength/2-separation;
position2=ArrayLength/2;
position3=ArrayLength/2+separation;
TSpectrum1 = gaussian(x,position1,width)+0.5.*gaussian(x,position1...
    +ArrayLength/10,width)+0.3.*gaussian(x,position1-ArrayLength/8,width);
TSpectrum2 = gaussian(x,position2,width)+0.2.*gaussian(x,position2...
    +ArrayLength/10,width)+0.5.*gaussian(x,position2-ArrayLength/10,width);
TSpectrum3 = gaussian(x,position3,width)+0.6.*gaussian(x,position3...
    +ArrayLength/10,width)+0.2.*gaussian(x,position3-ArrayLength/8,width);
y1=10 .^ (-TSpectrum1); yb1=broaden(y1,InstFunction,fa);
y2=10 .^ (-TSpectrum2); yb2=broaden(y2,InstFunction,fa);
y3=10 .^ (-TSpectrum3); yb3=broaden(y3,InstFunction,fa);
ReferenceSpectrum1 = -log10(yb1); ReferenceSpectrum2 = -log10(yb2);...
    ReferenceSpectrum3 = -log10(yb3);
y1=10 .^ (-absorbance(1) .* TSpectrum1); yb1=broaden(y1,InstFunction,fa);
y2=10 .^ (-absorbance(2) .* TSpectrum2); yb2=broaden(y2,InstFunction,fa);
y3=10 .^ (-absorbance(3) .* TSpectrum3); yb3=broaden(y3,InstFunction,fa);
Spectrum1 = -log10(yb1); Spectrum2 = -log10(yb2);Spectrum3 = -log10(yb3);

y=10 .^ (-absorbance(1) .* TSpectrum1 - absorbance(2) .* TSpectrum2 ...
    - absorbance(3) .* TSpectrum3);
yo=broaden(y,InstFunction,fa);
yobsd=StrayLight+yo+((noise/InstWidth).*randn(size(yo))).*sqrt(yo);   % Add simulated photon noise
yobsd=yobsd.*(1-StrayLight);
yobsd=yobsd.*(1+IzeroShift.*randn); % Random shifts in Izero 

% Three-component weighted regression
Background=ones(size(ReferenceSpectrum1));
w=y;  % Weight
WeightedMatrix=([w w w w] .* [Background ReferenceSpectrum1 ReferenceSpectrum2 ...
        ReferenceSpectrum3])\(-log10(yobsd) .* w); 
WeightedMatrix=abs(WeightedMatrix);

% TFit method
options = optimset('TolX',0.01);
start=WeightedMatrix(2:4)';
Spectra=[TSpectrum1,TSpectrum2,TSpectrum3];
lam=FMINSEARCH('fitM',start,options,yobsd,Spectra,InstFunction,StrayLight);
TFit=abs(lam(1:3));

% Plot individual spectra on Figure 2
figure(2);plot(x,yb1,'g',x,yb2,'b',x,yb3,'r',x,yobsd,'k:',x, gaussian(x,ArrayLength/2,InstWidth),'m--')
title('Transmission Spectra')
text(0,1.15,'     Green: Component 1   Blue: Component 2   Red: Component 3')  
text(0,1.1,'             Dotted: Mixture       Dashed: Inst Function')
axis([0 ArrayLength 0 1.2]); % Update plot

% Compare mathods on Figure 1
figure(1);loglog([.001 10],[.001 10],'k',absorbance(1),WeightedMatrix(2),'gx',absorbance(2),WeightedMatrix(3),'bx',...
absorbance(3),WeightedMatrix(4),'rx',absorbance(1),TFit(1),'go',absorbance(2),TFit(2),'bo',absorbance(3),TFit(3),'ro')
Rdisplay=round(1000*WeightedMatrix)/1000;
Tdisplay=round(1000*TFit)/1000;
xlabel(['R1=' num2str(Rdisplay(2)) '  R2=' num2str(Rdisplay(3)) '  R3=' num2str(Rdisplay(4)) '     T1=' num2str(Tdisplay(1)) '  T2= ' num2str(Tdisplay(2)) '  T3=' num2str(Tdisplay(3)) ]);
Adisplay=round(1000*absorbance)/1000;
title(['A1=' num2str(Adisplay(1)) '   A2=' num2str(Adisplay(2)) '   A3=' num2str(Adisplay(3)) '  Sepn=' num2str(round(separation)) '  InstWidth=' num2str(round(InstWidth)) '  Noise = ' num2str(round(10000*noise)/100) '%' ]); 
text(.002,7,'x = Weighted Regression method')
text(.002,4,'o = Transmission Fitting method')
axis([.001 10 .001 10]);