function [tau,D_ax,amb]=ambi_surf(fs,f_bb,cts,NPRI,nD,Dlims,ii_window,nzp_tau)
%
% Call [tau,D_ax,amb]=ambi_surf(fs,f_bb,cts,NPRI,nD,Dlims,ii_window,nzp_tau)
%
% fs            = Sample rate in Hz;
% f_bb          = Basebanding frequency in Hz;
% cts           = Complex time series;
% NPRI          = Number of PRIs to perform correlations with;
% nD            = Number of Dopplers;
% Dlims         = Dmax for [-Dmax,Dmax] or [Dmin,Dmax];
% ii_window     = Time-domain window: 0 for none, 1 for Kaiser, 2 for Hann;
% nzp_tau       = Zero pad factor in time domain: pad by 2^nzp_tau;
%
% Computes and plots D-tau ambiguity surface for a basebanded time series.
%
if NPRI<1 | round(NPRI)~=NPRI, error('Illegal NPRI: must be integer >= 1.'); end
npt_orig=length(cts); PRI=npt_orig/fs,
T_int=PRI*NPRI; nT_int=npt_orig*(NPRI+1);
npt=2^(nextpow2((NPRI+1)*npt_orig) + nzp_tau);
if length(Dlims)==1,
    dmax=Dlims;
    if dmax==0, dmax=.5/T_int, end
    if nD==0, dd=1/T_int; nD=10*ceil(dmax/dd); nD=2*floor(nD/2)+1; end
    D_ax=linspace(-dmax,dmax,nD);
elseif length(Dlims)==2,
    D_ax=linspace(Dlims(1),Dlims(2),nD); dmax=max(abs(Dlims));
    if nD==0, dd=1/T_int; nD=10*ceil(dmax/dd); nD=2*floor(nD/2)+1; end
else
    error('Dlims must be dmax or [dmin,dmax] '); 
end

% Construct C1 as having (NPRI+1) periods of signal, padded out to npt:
C1=fft([repmat(cts(:),[NPRI+1,1]); zeros(npt-nT_int,1)]);
if 0,
	f_window=linspace(.5,4,length(C1)).';
	C1=C1.*f_window;
	exp_window=exp(i*linspace(0,20*pi,length(C1)).');
	C1=C1.*exp_window;
	disp(['Testing linear window effects: ',num2str(10*log10(f_window(end)/f_window(1)))]);
end
C1c=conj(C1);
% Construct C0 as NPRI periods of signal, possibly Hann-shaded, heavily zero 
% padded out to nfft such that interpolation in frequency is very easy:
nfft=2^nextpow2(npt*10), npt_rep=NPRI*npt_orig;
if ii_window==1,
        hann=kaiser(npt_rep,8); hann=hann/norm(hann);
        C0=fft([hann.*repmat(cts(:),[NPRI,1]); zeros(nfft-npt_rep,1)]); 
elseif ii_window==2,
        hann=(8/3)*.5*(1 - cos(2*pi*[0:npt_rep-1]'/(npt_rep-1)));
        C0=fft([hann.*repmat(cts(:),[NPRI,1]); zeros(nfft-npt_rep,1)]); 
else
        C0=fft([repmat(cts(:),[NPRI,1]); zeros(nfft-npt_rep,1)]); 
end
disp(['Finished C0 fine-scale FFT with nfft = ',int2str(nfft)]);
Tw=npt/fs; df=1/Tw; 
Tw0=nfft/fs; df0=1/Tw0;
% Pad in frequency domain to avoid interpolating outside valid f interval:
npad=ceil(1.2*dmax/df0);
f0=f_bb + df0*[-nfft/2-npad:nfft/2-1+npad].';
% Construct C0_wrap, which is zero-padded in frequency and shifted such that
% negative frequencies appear first:
C0_wrap=[zeros(npad,1); fftshift(C0); zeros(npad,1)];
npt2=floor(npt/2); f1=f_bb + df*[-npt2:npt-npt2-1].';
%disp(f0([1,end])); disp(f1([1,end]));

amb=zeros(nD,npt_orig);
iidiag=0;
np2=round(npt_orig/2);
ii=[np2+1:npt_orig,1:np2];
tau=[0:npt_orig-1]/fs - np2/fs;

for jD=1:nD,
	D=D_ax(jD); 
	alpha=(f_bb + D)/f_bb;
	f0_stretch=f0*alpha;    % New f0 axis;
	C2_wrap=interp1(f0_stretch,C0_wrap,f1);
	C2=ifftshift(C2_wrap);
%       disp(f0_stretch([1,end])), disp(f1([1,end])),
%size(C1c), size(C2),
	c_tau=ifft(C1.*conj(C2));
    if mod(jD,10)==1, disp(['Finished jD= ',int2str(jD),' of ',int2str(nD)]); end
%%      size(c_tau), size(amb), keyboard
    if iidiag,
        figure(11);
        plot([0:npt-1]/fs,20*log10(abs(c_tau)),'r'); hold on;
        c1p=ifft(C1); c2p=ifft(C2);
        plot([0:npt-1]/fs,20*log10(abs(c1p)),'b');
        plot([0:npt-1]/fs,20*log10(abs(c2p)),'g');
        pause;
    end
    amb(jD,:)=10*log10(max(1e-20,real(c_tau(ii)).^2 + imag(c_tau(ii)).^2)).';
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.
