% Test read and plot WILBER seismograms  Z component only, test filter
%  SAGE 2000, L. Braile
% Mag 6.5 event off coast of central Chile
% Origin Time = 1998 sept 3 (day 246) 17:37:59.0
% -29.2934, -71.5471
% Station NNA (Nana, Peru) -11.9875, -76.8422
% Dist = 17.93 deg (1993 km), Az = 343 deg
% Start time = 17:40:54.0 (00:2:55 after event)
% dt = 0.05s, Nyquist = 10 HZ, npts = 10520
% Z component
load NNAz1.txt -ascii
%filter with butterworth bandpass filter
fcutlo = 0.001;  % low cutoff frequency expressed as a fraction of the Nyquist
fcuthi = 0.01;  % high cutoff frequency expressed as a fraction of the Nyquist
%  fcutlo < fcuthi
w = [fcutlo fcuthi];
[b,a]=butter(2,w);
x=NNAz1;
mx=mean(x);
x=x-mx;
z=filtfilt(b,a,x); % use forward and backward filtering to avoid phase shift
% Apply scal factor to filtered trace
sf = 1;
z = z*sf;
% create time vector
t=0:0.05:525.95;
t=t/60 + 2.9167; % 2 min 55 sec after event
% Plot a trace to mark the minute marks
dm=120:1200:9720;
xm=zeros(10520,1);
for i=1:9
   j=dm(i);
   xm(j)=20000;
end
xm=xm-180000;
% plot unfiltered and filtered seismograms
% shift filtered trace up
z = z + 150000;
plot(t,xm,'k',t,x,'k',t,z,'r')
hold on
set(gca,'fontsize',13)
str=['                                         '];
set(gca,'yticklabel',str')
xlabel('Travel time (minutes)','fontsize',14)
ylabel('Amplitude','fontsize',14)
title('Magnitude 6.5 earthquake, near coast of central Chile, 29.2934\rm^{o\rm} S, 71.5471\rm^{o\rm} W','fontsize',13)
text(3.1,395000,'Station = NNA (Nana, Peru, 11.9875\rm^{o\rm} S, 76.8422\rm^{o\rm} W)','fontsize',13)
text(3.1,425000,'Origin time = 17:37:59.0 GMT 1998/09/03, Depth = 27 km','fontsize',13)
text(3.1,365000,'Distance = 17.93\rm^{o\rm} (1993 km), Azimuth = 343\rm^{o\rm}','fontsize',13)
text(3.57,-145000,'17:42:00','fontsize',13)
text(3.47,-115000,'(hr:min:sec, GMT, 3 September 1998)','fontsize',13)
text(5.57,-145000,'17:44:00','fontsize',13)
text(7.57,-145000,'17:46:00','fontsize',13)
text(9.57,-145000,'17:48:00','fontsize',13)
text(4.07,-35000,'\uparrow{\itP}','fontsize',16,'color','c')
text(7.3,-35000,'{\itS}\uparrow','fontsize',16,'color','c')
text(2.85,15000,'Z','fontsize',14,'color','k')
text(3.3,15000,'\uparrow{\itup}','fontsize',14,'color','k')
text(3.3,-25000,'\downarrow{\itdown}','fontsize',14,'color','k')
x1=[8.3 8.3 10.95 10.95];
y1=[-100000 -110000 -110000 -100000];
plot(x1,y1,'-c','linewidth',1)
text(9.7,-90000,'Rayleigh','fontsize',14,'fontangle','italic','color','c')
axis([2.8 11 -200000 450000])
set(gca,'ytick',[-200000 450000])
hold off