% m-code for picking potential signals.
% Different picking algorithm will be tested.


clear
YEAR = '2000';
jday1 = 1;
jday2 = 1;
ARRAY1 = 'FNET';
ARRAY2 = 'BDSN';
ARRAY3 = 'GRSN';
BAND = 'BAND1';
EXT = 'AMPEQ';
Erad = 6371;
refgvel = 3.8;
%amp_cre = 1.3e-8;
amp_fac = [3 2 2];
idur4lta = 30;

% Outfile name
OUTNAME = '2000.JAN.EQ';
fnameout1 = sprintf('%s.%s.pick',OUTNAME,ARRAY1);
fnameout2 = sprintf('%s.%s.pick',OUTNAME,ARRAY2);
fnameout3 = sprintf('%s.%s.pick',OUTNAME,ARRAY3);


% Reading amp data.


jjday = jday1:jday2;
mmjday = length(jjday);

fid1 = fopen(fnameout1,'w');

for i = 1:mmjday
jday = jjday(i);

if (jday < 10)
julday = sprintf('00%i',jday);
elseif (jday < 100)
julday = sprintf('0%i',jday);
else
julday = sprintf('%i',jday);
end

% Picking for first ARRAY
tempname01 = sprintf('%s.%s.%s.%s.%s.amp',YEAR,julday,ARRAY1,BAND,EXT);
[dumm,dumm1,tmstla,tmstlo,stime,sdtime,sang,sdang,tempout01] = ...
readamp(tempname01);

mm01 = size(tempout01);
time01 = 0 : sdtime : sdtime*(mm01(1)-1);
time01 = (time01 + stime)/(3600*24) + jday;

baz01 = 0 : sdang : sdang*(mm01(2)-1);
baz01 = baz01 + sang;
[tmaxamp01,imaxamp01]=max(tempout01');
maxamp01 = tmaxamp01';
maxbaz01 = baz01(imaxamp01);
 
[time11,maxamp11,maxbaz11] = ...   
    cal_pick_lta(time01,maxamp01,maxbaz01,amp_fac,idur4lta);
%amp_cre = 1.0e-8;
%[time11,maxamp11,maxbaz11] = ...
%    cal_pick(time01,maxamp01,maxbaz01,amp_fac1,amp_cre);


% Select picks for given julian day.
mm11 = size(time11);
jj = 1;
for j = 1 : mm11(2)
    if (time11(j) < jday+1) & (time11(j) >= jday)
        time21(jj)=time11(j);
        maxamp21(jj)=maxamp11(j);
        maxbaz21(jj)=maxbaz11(j);
        jj = jj + 1;
    end
end

% Picking for second ARRAY
tempname02 = sprintf('%s.%s.%s.%s.%s.amp',YEAR,julday,ARRAY2,BAND,EXT);
[dumm,dumm1,tmstla,tmstlo,stime,sdtime,sang,sdang,tempout02] = ...
readamp(tempname02);

mm02 = size(tempout02);
time02 = 0 : sdtime : sdtime*(mm02(1)-1);
time02 = (time02 + stime)/(3600*24) + jday;

baz02 = 0 : sdang : sdang*(mm02(2)-1);
baz02 = baz02 + sang;
[tmaxamp02,imaxamp02]=max(tempout02');
maxamp02 = tmaxamp02';
maxbaz02 = baz02(imaxamp02);
%amp_fac2 = [2 0 1.5];

 
[time12,maxamp12,maxbaz12] = ...
    cal_pick_lta(time02,maxamp02,maxbaz02,amp_fac,idur4lta);

% Select picks for given julian day.
mm12 = size(time12);
jj = 1;
for j = 1 : mm12(2)
    if (time12(j) < jday+1) & (time12(j) >= jday)
        time22(jj)=time12(j);
        maxamp22(jj)=maxamp12(j);
        maxbaz22(jj)=maxbaz12(j);
        jj = jj + 1;
    end
end


% Picking for third ARRAY
tempname03 = sprintf('%s.%s.%s.%s.%s.amp',YEAR,julday,ARRAY3,BAND,EXT);
[dumm,dumm1,tmstla,tmstlo,stime,sdtime,sang,sdang,tempout03] = ...
readamp(tempname03);

mm03 = size(tempout03);
time03 = 0 : sdtime : sdtime*(mm03(1)-1);
time03 = (time03 + stime)/(3600*24) + jday;

baz03 = 0 : sdang : sdang*(mm03(2)-1);
baz03 = baz03 + sang;
[tmaxamp03,imaxamp03]=max(tempout03');
maxamp03 = tmaxamp03';
maxbaz03 = baz03(imaxamp03);
 
[time13,maxamp13,maxbaz13] = ...
    cal_pick_lta(time03,maxamp03,maxbaz03,amp_fac,idur4lta);

% Select picks for given julian day.
mm13 = size(time13);
jj = 1;
for j = 1 : mm13(2)
    if (time13(j) < jday+1) & (time13(j) >= jday)
        time23(jj)=time13(j);
        maxamp23(jj)=maxamp13(j);
        maxbaz23(jj)=maxbaz13(j);
        jj = jj + 1;
    end
end

end

clf
subplot(3,1,1),plot(time01,maxamp01)
hold on
plot(time21,maxamp21,'r.')
xlim([jday1 jday2+1])
subplot(3,1,2),plot(time02,maxamp02)
hold on
plot(time22,maxamp22,'r.')
xlim([jday1 jday2+1])
subplot(3,1,3),plot(time03,maxamp03)
hold on
plot(time23,maxamp23,'r.')
xlim([jday1 jday2+1])
