%clear
close

amax=0.00003
load waveb_d;
size(waveb_d)
nob=waveb_d(1,1)
nt=waveb_d(1,2)
dt=waveb_d(1,3)

for j=1:nt,
  t(j)=dt*(j-1);
end
f=max(t)


amaxx=0.0;
amaxy=0.0;
amaxz=0.0;
for i=1:nob,
  for j=1:nt,
    axx=abs(waveb_d(nt*(i-1)+1+j,1));
    if axx > amaxx amaxx=axx; end
    ayy=abs(waveb_d(nt*(i-1)+1+j,2));
    if ayy > amaxy amaxy=ayy; end
    azz=abs(waveb_d(nt*(i-1)+1+j,2));
    if azz > amaxz amaxz=azz; end
  end
end
amaxx=amaxx*0.5
amaxy=amaxy*0.5
amaxz=amaxz*0.5
fac=0.5
ax=amaxx/fac
bx=amaxx/fac
ay=amaxy/fac
by=amaxy/fac
az=amaxz/fac
bz=amaxz/fac

subplot(1,3,1)
hold on
for i=1:nob,
  for j=1:nt,
    wx(j)=waveb_d(nt*(i-1)+j+1,1)+ax+amaxx*(i-1)/2;
  end
  plot(t,wx)
end
hold off
title('NS')
set(gca, 'box','on')
axis([0 f 0 ax+bx+amaxx*(nob-1)/2])
ylabel('Velocity (cm/s)')
xlabel('time (sec)')

subplot(1,3,2)
hold on
for i=1:nob,
  for j=1:nt,
    wy(j)=waveb_d(nt*(i-1)+j+1,2)+ay+amaxy*(i-1)/2;
  end
  plot(t,wy)
end
hold off
title('EW')
set(gca, 'box','on')
axis([0 f 0 ay+by+amaxy*(nob-1)/2])
xlabel('time (sec)')

subplot(1,3,3)
hold on
for i=1:nob,
  for j=1:nt,
    wz(j)=waveb_d(nt*(i-1)+j+1,3)+az+amaxz*(i-1)/2;
  end
  plot(t,wz)
end
hold off
title('UD')
set(gca, 'box','on')
axis([0 f 0 az+bz+amaxz*(nob-1)/2])
xlabel('time (sec)')

print grfftbx1.ps


