function [xmax,d1max,d2max,d3max,d4max]=quad_peak(xax,data1,data2,data3,data4);
%
% Call [xmax,d1max,d2max,d3max,d4max]=quad_peak(xax,data1,data2,data3,data4);
%
% Finds max of each row of data1, fits quadratic in region of maximum,
% and returns xmax as value of xax corresponding to maximum.  Also finds peak values in data1, data2,
% and data3 at xmax.
%
% Inputs:
% xax       = x-axis for each row of data1;
% data1     = data matrix, each row of which is searched for the peak;
% Outputs:
% xmax      = output xax value at which peak occurred;
% d1max     = value of data matrix corresponding to max;
% data2 ... = additional data matrices that are evaluated at xmax;
% d2max ... = values of data2 ... at xmax
%
[nrow,ncol]=size(data1);
xmax=zeros(nrow,1); d1max=zeros(nrow,1);
if nargin<3, data2=[]; d2max=[]; else d2max=zeros(nrow,1); end
if nargin<4, data3=[]; d3max=[]; else d3max=zeros(nrow,1); end
if nargin<5, data4=[]; d4max=[]; else d4max=zeros(nrow,1); end
for jrow=1:nrow,
	[dummy,jcolmax]=max(data1(jrow,:));
	if jcolmax==1 | jcolmax==ncol,
		xmax(jrow)=xax(jcolmax);
		d1max(jrow)=data1(jrow,jcolmax);
		if ~isempty(data2), d2max(jrow)=data2(jrow,jcolmax); end
		if ~isempty(data3), d3max(jrow)=data3(jrow,jcolmax); end
		if ~isempty(data4), d4max(jrow)=data4(jrow,jcolmax); end
	else
		iicubfit=[jcolmax-1:jcolmax+1];
		abc=polyfit(xax(iicubfit),data1(jrow,iicubfit),2);
		xmax(jrow)=-abc(2)/(2*abc(1));
		d1max(jrow)=abc(1)*xmax(jrow)^2 + abc(2)*xmax(jrow) + abc(3);
		if ~isempty(data2),
% Fit data2 to a quadratic over max region and evaluate at xmax(jrow):
			abc=polyfit(xax(iicubfit),data2(jrow,iicubfit),2);
			d2max(jrow)=abc(1)*xmax(jrow)^2 + abc(2)*xmax(jrow) + abc(3);
		end
		if ~isempty(data3),
% Fit data3 to a quadratic over max region and evaluate at xmax(jrow):
			abc=polyfit(xax(iicubfit),data3(jrow,iicubfit),2);
			d3max(jrow)=abc(1)*xmax(jrow)^2 + abc(2)*xmax(jrow) + abc(3);
		end
		if ~isempty(data4),
% Fit data4 to a quadratic over max region and evaluate at xmax(jrow):
			abc=polyfit(xax(iicubfit),data4(jrow,iicubfit),2);
			d4max(jrow)=abc(1)*xmax(jrow)^2 + abc(2)*xmax(jrow) + abc(3);
		end
	end
end
%
% COPYRIGHT (c) 2008
% 
% Applied Research Laboratories, The University of Texas at Austin
% (ARL:UT)
% All rights reserved.
% 
% This material may be reproduced by or for the U.S. Government
% only. This notice must appear in all copies of this file and its
% derivatives. Further distribution outside the U.S. Government is
% prohibited without the explicit written authorization from the
% Director of the Environmental Sciences Laboratory, ARL:UT.
% 
% THIS SOFTWARE IS PROVIDED "AS IS" WITHOUT WARRANTY OF ANY KIND,
% EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO THE
% IMPLIED WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
% PURPOSE, OR NONINFRINGEMENT.
% 
% IN NO EVENT WILL ARL:UT BE LIABLE FOR ANY DAMAGES WHATSOEVER ARISING
% OUT OF THE USE OR INABILITY TO USE THIS SOFTWARE, EVEN IF ARL:UT HAS
% BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGES.
