%file name "OFEWWV.M"
% A "m" file for Matlab. After making a "x_dat" file from "grfftm.f", you
% can plot waveforms using Matlab. A post script data is stored in the 
% "x.ps" file.

clear
close
factor=2.5

fig=figure('Units','centimeters','Position',[0.0 0.0 21.0 29.5]...
,'PaperUnits','centimeters','PaperPosition',[0.0 0.0 21.0 29.5]);
acc_inp=uigetfile('');
%iwanuma_v=load (iwanuma2);
%size(iwanuma_v)
NameHeader= '20110311144600'
%nt=size(iwanuma_v,1);
acc=load (acc_inp);
%load acc;
dat=acc;
n=size(dat,1)

%ACCwave in EW Direction of Office Towe
ch1=9
ch2=35
ch3=38
ch4=31
ch5=19

c1='ch09: O28F EW';
c2='ch35: O22F EW';
c3='ch38: O15F EW';
c4='ch31: O08F EW';
c5='ch19:-100m EW';
amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*4;
  amp2(j)=amp2(j)+amx*3;
  amp3(j)=amp3(j)+amx*2;  
  amp4(j)=amp4(j)+amx*1;
end
ymin=min(amp5)-amx/4;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-')
text(xmax*0.8,ymin+amx*0.1,c5)
text(xmax*0.8,ymin+amx*1.1,c4)
text(xmax*0.8,ymin+amx*2.1,c3)
text(xmax*0.8,ymin+amx*3.1,c2)
text(xmax*0.8,ymin+amx*4.1,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
%axis([140 200 ymin*1.1 ymax*1.01])
title([NameHeader, ' Quake (O28F EW: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print Office_EW_Acc.ps

%ACCwave in NS1 Direction of Office Towe
%clear
ch1=10
ch2=36
ch3=39
ch4=32
ch5=20

c1='ch10: O28F NS1';
c2='ch36: O22F NS1';
c3='ch39: O15F NS1';
c4='ch32: O08F NS1';
c5='ch20:-100m NS';
amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*4;
  amp2(j)=amp2(j)+amx*3;
  amp3(j)=amp3(j)+amx*2;  
  amp4(j)=amp4(j)+amx*1;
end
ymin=min(amp5)-amx/4;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-')
text(xmax*0.8,ymin+amx*0.1,c5)
text(xmax*0.8,ymin+amx*1.1,c4)
text(xmax*0.8,ymin+amx*2.1,c3)
text(xmax*0.8,ymin+amx*3.1,c2)
text(xmax*0.8,ymin+amx*4.1,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
%axis([140 200 ymin*1.1 ymax*1.01])
title([NameHeader, ' Quake (O28F NS1: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print Office_NS1_Acc.ps

%ACCwave in NS2 Direction of Office Towe
%clear
ch1=34
ch2=37
ch3=40
ch4=33
ch5=20

c1='ch34: O28F NS2';
c2='ch37: O22F NS2';
c3='ch40: O15F NS2';
c4='ch33: O08F NS2';
c5='ch20:-100m NS';
amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*4;
  amp2(j)=amp2(j)+amx*3;
  amp3(j)=amp3(j)+amx*2;  
  amp4(j)=amp4(j)+amx*1;
end
ymin=min(amp5)-amx/4;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-')
text(xmax*0.8,ymin+amx*0.1,c5)
text(xmax*0.8,ymin+amx*1.1,c4)
text(xmax*0.8,ymin+amx*2.1,c3)
text(xmax*0.8,ymin+amx*3.1,c2)
text(xmax*0.8,ymin+amx*4.1,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
%axis([140 200 ymin*1.1 ymax*1.01])
title([NameHeader, ' Quake (O28F NS2: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print Office_NS2_Acc.ps

%ACCwave in EW Direction of University Tower
ch1=1
ch2=4
ch3=25
ch4=26
ch5=29
ch6=13
ch7=19

c1='ch01: U29F EW';
c2='ch04: U22F EW';
c3='ch25: U16F EW';
c4='ch26: U08F EW';
c5='ch29: U01F EW';
c6='ch13: UB6F EW';
c7='ch19:-100m EW';

amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  amp6(j)=dat(j,ch6+1);
  amp7(j)=dat(j,ch7+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*6;
  amp2(j)=amp2(j)+amx*5;
  amp3(j)=amp3(j)+amx*4;
  amp4(j)=amp4(j)+amx*3;
  amp5(j)=amp5(j)+amx*2;
  amp6(j)=amp6(j)+amx*1;
end
ymin=min(amp7)-amx/5;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-',time,amp6,'-',time,amp7,'-')
text(xmax*0.8,ymin+amx*0.2,c7)
text(xmax*0.8,ymin+amx*1.2,c6)
text(xmax*0.8,ymin+amx*2.2,c5)
text(xmax*0.8,ymin+amx*3.2,c4)
text(xmax*0.8,ymin+amx*4.2,c3)
text(xmax*0.8,ymin+amx*5.2,c2)
text(xmax*0.8,ymin+amx*6.2,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
title([NameHeader, ' Quake (U29F EW: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print University_EW_Acc.ps

%ACCwave in NS1 Direction of University Tower

ch1=2
ch2=5
ch3=7
ch4=27
ch5=30
ch6=14
ch7=20

c1='ch02: U29F NS1';
c2='ch05: U22F NS1';
c3='ch07: U16F NS1';
c4='ch27: U08F NS1';
c5='ch30: U01F NS1';
c6='ch14: UB6F NS1';
c7='ch20:-100m NS';

amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  amp6(j)=dat(j,ch6+1);
  amp7(j)=dat(j,ch7+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*6;
  amp2(j)=amp2(j)+amx*5;
  amp3(j)=amp3(j)+amx*4;
  amp4(j)=amp4(j)+amx*3;
  amp5(j)=amp5(j)+amx*2;
  amp6(j)=amp6(j)+amx*1;
end
ymin=min(amp7)-amx/5;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-',time,amp6,'-',time,amp7,'-')
text(xmax*0.8,ymin+amx*0.2,c7)
text(xmax*0.8,ymin+amx*1.2,c6)
text(xmax*0.8,ymin+amx*2.2,c5)
text(xmax*0.8,ymin+amx*3.2,c4)
text(xmax*0.8,ymin+amx*4.2,c3)
text(xmax*0.8,ymin+amx*5.2,c2)
text(xmax*0.8,ymin+amx*6.2,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
%title([980820 'Quake (U29F NS1: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
title([NameHeader, ' Quake (U29F NS2: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print University_NS1_Acc.ps

%ACCwave in NS2 Direction of University Tower

ch1=3
ch2=6
ch3=8
ch4=28
ch5=12
ch6=15
ch7=20

c1='ch03: U29F NS2';
c2='ch06: U22F NS2';
c3='ch08: U16F NS2';
c4='ch28: U08F NS2';
c5='ch12: U01F NS2';
c6='ch15: U06B NS2';
c7='ch20:-100m NS';

amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  amp6(j)=dat(j,ch6+1);
  amp7(j)=dat(j,ch7+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*6;
  amp2(j)=amp2(j)+amx*5;
  amp3(j)=amp3(j)+amx*4;
  amp4(j)=amp4(j)+amx*3;
  amp5(j)=amp5(j)+amx*2;
  amp6(j)=amp6(j)+amx*1;
end
ymin=min(amp7)-amx/5;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-',time,amp6,'-',time,amp7,'-')
text(xmax*0.8,ymin+amx*0.2,c7)
text(xmax*0.8,ymin+amx*1.2,c6)
text(xmax*0.8,ymin+amx*2.2,c5)
text(xmax*0.8,ymin+amx*3.2,c4)
text(xmax*0.8,ymin+amx*4.2,c3)
text(xmax*0.8,ymin+amx*5.2,c2)
text(xmax*0.8,ymin+amx*6.2,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
title([NameHeader, ' Quake (U29F NS2: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print University_NS2_Acc.ps

%ACCwave in UD Direction of building
ch1=11
ch2=16
ch3=17
ch4=18
ch5=21

c1='ch11: U29F UD';
c2='ch16: UB6F UD';
c3='ch17: O28F UD';
c4='ch18: OB6F UD';
c5='ch21:-100m UD';
amax=0;
for j=1:n,
  time(j)=dat(j,1);
  amp1(j)=dat(j,ch1+1);
  amp2(j)=dat(j,ch2+1);
  amp3(j)=dat(j,ch3+1);
  amp4(j)=dat(j,ch4+1);
  amp5(j)=dat(j,ch5+1);
  if amp1(j) > amax
     amax=amp1(j);
     permax=time(j);
  end
end
amax
permax

amx=amax*factor;
for j=1:n,
  amp1(j)=amp1(j)+amx*4;
  amp2(j)=amp2(j)+amx*3;
  amp3(j)=amp3(j)+amx*2;  
  amp4(j)=amp4(j)+amx*1;
end
ymin=min(amp5)-amx/4;
ymax=max(amp1)
xmin=min(time)
xmax=max(time)

plot(time,amp1,'-',time,amp2,'-',time,amp3,'-',time,amp4,'-',time,amp5,'-')
text(xmax*0.8,ymin+amx*0.1,c5)
text(xmax*0.8,ymin+amx*1.1,c4)
text(xmax*0.8,ymin+amx*2.1,c3)
text(xmax*0.8,ymin+amx*3.1,c2)
text(xmax*0.8,ymin+amx*4.1,c1)
axis([xmin xmax*1.01 ymin ymax*1.02])
%axis([140 200 ymin*1.1 ymax*1.01])
title([NameHeader, '  Quake (University 29F UD: Max =',num2str(amax),' at ',num2str(permax),' sec )'])
%ylabel('Velocity (cm/sec)')
%ylabel('Displacement (cm)')
ylabel('Acceleration (gal)')
xlabel('time (sec)')
grid

print Building_UD_Acc.ps