%tddop.m

fc=75;                        % carrier frequency
cyc_dig=2;                    % cycles per digit ( or Q)
Pn=4;                         % periods per transmission
law = [ 1 0 0 0 0 1 1];          % 103 oct
sw_iap=0;                     % iap align circums
Wn=40;
%vsearch=-5:0.2:5;
vsearch=-10:2:10;
vg=-5.0;
vsearch=vsearch+vg;

nr=length(law)-1;             % length of sequence register
L = 2.^nr - 1;                % digits per period
sam_cyc=4;                    % samples per cycle
sam_dig=sam_cyc*cyc_dig;      % samples per digit
fs=sam_cyc*fc;                % sampling frequency
N=sam_cyc*cyc_dig*L;          % samples per period
T=N/fs;                       % period length
dt=1/fs;
taxis=0:dt:T-dt;

snrdb=15.0;
%sg( 0.0, 0.0, fc, cyc_dig, Pn, 15.0, law);
sg( 1.301, -vg, fc, cyc_dig, Pn, 15.0, law);

sam=readadf('proto.sam');

fprintf(1, 'rms input power: %6.2f (bits)\n', std(sam));

dm=msq( sam, fc, cyc_dig, Pn, Wn, sw_iap, vsearch, law);

[ dmmax, dind]=max( abs(dm)');
[ dmpk, vind]=max( dmmax); 
dmang=180*angle(dm(vind, dind(vind)))/pi;
vmax=vsearch(vind);
tmax=taxis(dind(vind));

y1=abs(dm(vind, dind(vind)-1));
y2=abs(dm(vind, dind(vind)));
y3=abs(dm(vind, dind(vind)+1));

tte=taxis(dind(vind)) + dt*(y1-y3)/(y1-2*y2+y3)/2;

md2rms=sqrt(log(2));

dmmed=median( abs(dm)')/md2rms;
dmmedfl=median(dmmed);
dmmedgr=dmmed(vind);

fprintf(1,'pk: %6.2f ang: %7.2f tmax: %7.4f tte %7.4f vmax: %7.3f\n',...
        20*log10(dmpk), dmang, tmax, tte, vmax);
fprintf(1,'floor: %6.2f groove: %6.2f\n',...
        20*log10(dmmedfl), 20*log10(dmmedgr))
fprintf(1,'snr_floor: %6.2f snr_groove: %6.2f\n',...
        20*log10(dmpk/dmmedfl), 20*log10(dmpk/dmmedgr))

brt=1.0; sat=1.0; cmap=[];
for i=0:63
  hue = 240.0-270.0*i/64;
  [red,green,blue]= hsb( hue, sat, brt);
  cmap=[ cmap; [red, green, blue]];
end
clear brt sat i hue red green blue;

clf reset
set( gcf, 'nextplot', 'replace');
%set( gcf, 'units', 'inches');
%set( gcf, 'position', [ 0 3 7 4]);
%set( gcf, 'units', 'pixels');
colormap( cmap);

dim3=1;
if dim3
%  tind=find( taxis>=0 & taxis<=0.5);
  tind=find( taxis>=max(taxis)-0.5 & taxis<=max(taxis));
  gcs=surf( taxis(tind), vsearch, (abs(dm(:,tind))));
  set( gca, 'fontname', 'times');
  xlabel('Time Delay (seconds)')
  ylabel('Opening Velocity (m/s)')
  colormap( cmap)
  shading interp
  set( gcs, 'edgecolor', 'k')
  set( gcs, 'meshstyle', 'row') 
  view( -15, 70)
  % view(viewmtx( -15, 65, 60))
  pause
end

plot( taxis, 20*log10(abs(dm)), 'b')
hold on
plot( taxis, 20*log10(abs(dm(vind,:))), 'r')
plot( [min(taxis) max(taxis)], 20*log10([dmmedfl dmmedfl]), 'g')
plot( [min(taxis) max(taxis)], 20*log10([dmmedgr dmmedgr]), 'g')
hold off
set( gca, 'fontname', 'times');
xlabel('Time Delay (seconds)')
set( gca, 'xlim', [min(taxis) max(taxis)]);
grid
pause

if dim3
  plot( vsearch, 20*log10(dmmax), 'r');
  hold on
  plot( vsearch, 20*log10(dmmed), 'g');
  hold off
  set( gca, 'fontname', 'times');
  xlabel('Opening velocity (m/s)')
  grid
  pause
end


















