clear
clf
load x.out
load t.out
a=size(t);
n=round(a/2);
load uhov.out
load vhov.out
load whov.out
t2=t(n:a);
u2hov=uhov(n:a,:);
v2hov=vhov(n:a,:);
w2hov=whov(n:a,:);
t=t2;
uhov=u2hov;
vhov=v2hov;
whov=w2hov;
xs=size(x);
x=x(2:xs(2)-1);
xs=size(x);
xg=round(x*xs(2)/x(xs(2)));
im=round(xs(2));
im2=round(im/2);
xs2=xs/2;
ts=size(t);
jm=round(ts(1));
ts2=ts(1)/2;
j=round(ts2(1));
ftmax=ts/t(ts(1));
%ftmax=t*ts(1)./t(ts(1));
tg=4.*t*ftmax(1)/t(jm);
i=input('Enter u,v, or w','s');
if i=='u' 
   ut=uhov(:,3);
   utf=fft(ut);
   a=abs(utf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
%   c=2.*log(b);
   c=b;
   c2=c(1:j);
   c2(1)=0.0;
   t2=tg(1:j);
   t2=t2-t2(1);
   h=plot(t2,c2);
   set(h,'linewidth',2)
   xlabel('Frequency (cycles/day)','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Zonal velocity','fontweight','bold')
   pause
   uk=uhov(jm,:);
   ukf=fft(uk);
   a=abs(ukf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
%   c=2.*log(b);
   c=b;
   c2=c(1:im2);
   c2(1)=0.0;
   x2=xg(1:im2);
   x2=x2-1;
   h=plot(x2,c2);
   set(h,'linewidth',2)
   xlabel('Wavenumber','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Zonal velocity','fontweight','bold')
elseif i == 'v'
   ut=vhov(:,3);
   utf=fft(ut);
   a=abs(utf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
%   c=2.*log(b);
   c=b;
   c2=c(1:j);
   c2(1)=0.0;
   t2=tg(1:j);
   t2=t2-t2(1);
   h=plot(t2,c2);
   set(h,'linewidth',2)
   xlabel('Frequency (cycles/day)','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Meridional velocity','fontweight','bold')
   pause
   uk=vhov(jm,:);
   ukf=fft(uk);
   a=abs(ukf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
%   c=2.*log(b);
   c=b;
   c2=c(1:im2);
   c2(1)=0.0;
   x2=xg(1:im2);
   x2=x2-1;
   h=plot(x2,c2);
   set(h,'linewidth',2)
   xlabel('Wavenumber','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Meridional velocity','fontweight','bold')
elseif i == 'w'
   ut=whov(:,3);
   utf=fft(ut);
   a=abs(utf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
%   c=2.*log(b);
   c=b;
   c2=c(1:j);
   c2(1)=0.0;
   t2=tg(1:j);
   t2=t2-t2(1);
   h=plot(t2,c2);
   set(h,'linewidth',2)
   xlabel('Frequency (cycles/day)','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Vertical velocity','fontweight','bold')
   pause
   uk=whov(jm,:);
   ukf=fft(uk);
   a=abs(ukf);
   amin=min(a);
   a(1)=amin;
   b=a+1-amin;
   c=b;
   c2=c(1:im2);
   c2(1)=0.0;
   x2=xg(1:im2);
   x2=x2-1;
   h=plot(x2,c2);
   set(h,'linewidth',2)
   xlabel('Wavenumber','fontweight','bold')
   ylabel('Power','fontweight','bold')
   title('Vertical velocity','fontweight','bold')
  end   
   
