% readaburst
rot = 19

bn = input('Enter a burst number: ');
if(1),
  if(isunix~=1),
    dn = 'd:\data\GH2001\d1\MIB\';
    hdrfn = [dn,'MIB1001.HDR'];
  else
        dn = '/data/GH2001/d1/MIB/';
	hdrfn = [dn,'MIB1ADV.hdr'];
  end
  adrfn = [dn,'MIB1001.adr'];
  gb = 3:771';
  adv = read_hdradv( hdrfn, gb) ;
  burst = adr2mat( adrfn, [bn bn]);
end

adv_elev = 0.01*ecorr(adv.probe_ht);
zparos = adv_elev+(1.19-0.807);
zadv = adv_elev - 0.16;
az = median(burst.heading{1})

u = 0.01*burst.Vx{1};
w = 0.01*burst.Vz{1};
v = 0.01*burst.Vy{1};
p = 0.001*burst.extpress{1}-10.14;

%figure
%MIBtripod_pca(az, 100*u, 100*v );

% rotate into geographic coords
[u,v]=rotadv(u,v,az+rot);
burst.samprate
ws = wave_stats(u,v,p,burst.samprate)
S = ubstats( u, v, burst.samprate )

data = [p u v ];
fs = burst.metadata.samprate(1);
zdep = mean(p)+zparos(bn);
layout = [0. 0. 0.; 0. 0. 0.; zparos(bn) zadv(bn) zadv(bn) ];
datatypes = {'pres','velx','vely'};
fbins = [0.02:.01:.2];
dirs = [-180:1:180];
% S = dirspec(data,layout,datatypes,zdep,fs,fbins,dirs, ...
%	    'EMEP',...
%	    'NFFT',256,'XAXISDIR',262.5,'PLOTTYPE',1);
	  
% plotspec(fbins,dirs,S,4,265.5);
% figure
% pcolor(fbins,dirs,S)


[Gpp, Guu, Gvv, Hrms, Hrmsp] = puv(p, u, v, zdep, zparos(bn), zadv(bn), fs)