%% plot 24 hour data time series and helicorder record
%%%%%%%%%%%%%%%%%%%%%
% read in 24 hrs data
%%%%%%%%%%%%%%%%%%%%%

%cd '/Users/Jeff/Documents/DATA/KILAUEA2008/SEGY' % move to data directory
cd '/Users/andrea/Desktop/SEGY_tungu09'

for day=0:5
    start_time = [166+day 0 0 0]; % [jday hour min sec] of start
    window = 3600*24; % 24 hours or 86400 s of data
    das = 'A860'; % name of seismic station to read in

    [seisV,hdr] = segycut(start_time,window,das,'1');
    ts = (1:length(seisV))/hdr.sps; % create time series vector
    % filtered seisV
    sps = 100; % get sample rate from header
    npoles = 2; % specify order of filter
    low = 1/600; % corner frequency at 10 Hz

    [B2,A2] = butter(npoles,[low]/(sps/2),'high')
    seisVfiltered = filter(B2,A2,seisV);


% using reshape

    matrixSeis= reshape(seisVfiltered,10*60*100,groups);
  
    meanRSM = mean(abs(matrixSeis)); % array

    time=linspace(0,24,length(meanRSM));
% now plot data filtered traces in several bands
    subplot(6,1,day+1)
    plot(time,meanRSM,'k-')
    grid on
    title( strcat('TGRS station RASM julian day  ',num2str(166+day)))
    xlabel('time (hours)'),ylabel('meanRSM')
    xlim([0 24])
    set(gca,'xtick',[0 5 10 15 20 24])
end