function P=SE_term_maxing(npt_max,ntimes,npt_xax)
%
% Call P=SE_term_maxing(npt_max,ntimes,npt_xax);
%
% Finds the losses due to maximizing over a set of Gaussian-distributed
% random numbers.  The set has maximum size of npt_max.  The calculation
% is performed in a Monte Carlo sense ntimes times, and the size of the set of
% numbers is sampled npt_xax times, from 1 to npt_max.  A cubic fit of the
% function max_loss in dB versus log(# pts maxed over), and the resulting
% coefficients are returned in P.
%
% The result for P=SE_term_maxing(1000,500,100) is:
%       P=[0.1778   -1.7162    7.0875    0.1623], or
%       max_loss (in dB) = 0.1778*x^3 + -1.7162*x^2 + 7.0875*x + 0.1623,
%           where x=log10(# npts maxed over).
%
% To compute maxing loss (in dB) for npts max'ed over:
%   max_loss=polyval(P,log10(npt_maxed_over));
% 
%
% 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.
x=randn(npt_max,ntimes);
x2=x.^2;
npt_log_ax=linspace(log10(1),log10(npt_max),npt_xax);
npt_ax=unique(round(10.^npt_log_ax)); npt_xax=length(npt_ax);
npt_log_ax=log10(npt_ax);
max_loss=zeros(1,npt_xax);
for j=1:npt_xax,
    npt=npt_ax(j);
    max_loss(j)=10*log10(mean(max(x2(1:npt,:),[],1)));
end
figure
plot_cfill(npt_log_ax,max_loss);
% Fit a cubic, which works quite well:
P=polyfit(npt_log_ax,max_loss,3);
hold on;
yfit=polyval(P,npt_log_ax);
plot_cfill(npt_log_ax,yfit,'none','b',5,'-',1.5);
xlabel('log # pts maxed over');
ylabel('Maxing Loss - dB');
plot_set3(gca,1.5,'Arial',14,1);
plot_size(gca,'in',5,2);
grid on;
axis tight;

