Программа по построению волновых форм спектрограмм землетрясений и
микросейсм до и после землетрясений
p='f:\ATR\';='f:\';
%nomer='13';'01'='2014204214130320_POLE__1_2.atr';='2014204214130320_POLE__1_3.atr';='2014204224130320_POLE__1_2.atr';='2014204224130320_POLE__1_3.atr';=datenum(2014,07,23,21,41,30.320);=datenum(2014,07,23,22,48,17);='M=1.1,Az=000';=50; timePosle=100;=-125; z2=-90;'02'='2014205135635320_POLE__1_2.atr';='2014205135635320_POLE__1_3.atr';='2014205145635320_POLE__1_2.atr';='2014205145635320_POLE__1_3.atr';=datenum(2014,07,24,13,56,35.320);=datenum(2014,07,24,14,48,9); %!!!!!='M=2.4,Az=097';=50; timePosle=100;1=-120; z2=-80;
case '04'
filenameX1='2014205165635320_POLE__1_2.atr';
filenameY1='2014205165635320_POLE__1_3.atr';='2014205175635320_POLE__1_2.atr';='2014205175635320_POLE__1_3.atr';=datenum(2014,07,24,16,56,35.320);=datenum(2014,07,24,17,12,3);='M=3.8,Az=027';=50; timePosle=250;=-120; z2=-50;'05'='2014205175635320_POLE__1_2.atr';='2014205175635320_POLE__1_3.atr';='2014205185635320_POLE__1_2.atr';='2014205185635320_POLE__1_3.atr';=datenum(2014,07,24,17,56,35.320);=datenum(2014,07,24,18,40,42);='M=1.8,Az=107';=50; timePosle=100;=-125; z2=-80;'06'='2014206055635320_POLE__1_2.atr';='2014206055635320_POLE__1_3.atr';='2014206065635320_POLE__1_2.atr';='2014206065635320_POLE__1_3.atr';=datenum(2014,07,25,05,56,35.320);=datenum(2014,07,25,06,07,2);='M=0.9,Az=239';=50; timePosle=70;=-115; z2=-95;
Продолжение приложения Б'07'='2014206115635320_POLE__1_2.atr';='2014206115635320_POLE__1_3.atr';='2014206125635320_POLE__1_2.atr';='2014206125635320_POLE__1_3.atr';=datenum(2014,07,25,11,56,35.320);=datenum(2014,07,25,12,43,6);='M=0.7,Az=229';=50; timePosle=40;=-125; z2=-80;'08'='2014206175635320_POLE__1_2.atr';='2014206175635320_POLE__1_3.atr';='2014206185635320_POLE__1_2.atr';='2014206185635320_POLE__1_3.atr';=datenum(2014,07,25,17,56,35.320);=datenum(2014,07,25,18,04,27);='M=2.0,Az=117';=50; timePosle=150;=-125; z2=-85;'09'='2014206205635320_POLE__1_2.atr';='2014206205635320_POLE__1_3.atr';='2014206215635320_POLE__1_2.atr';='2014206215635320_POLE__1_3.atr';=datenum(2014,07,25,20,56,35.320);=datenum(2014,07,25,21,24,7);='M=1.6,Az=249';=50; timePosle=100;
z1=-125; z2=-95;
case '10'='2014207175635320_POLE__1_2.atr';='2014207175635320_POLE__1_3.atr';='2014207185635320_POLE__1_2.atr';='2014207185635320_POLE__1_3.atr';=datenum(2014,07,26,17,56,35.320);=datenum(2014,07,26,18,06,31);='M=1.7,Az=027';=50; timePosle=100;=-125; z2=-85;'11'='2014208055635320_POLE__1_2.atr';='2014208055635320_POLE__1_3.atr';='2014208065635320_POLE__1_2.atr';='2014208065635320_POLE__1_3.atr';=datenum(2014,07,27,05,56,35.320);=datenum(2014,07,27,06,32,58);='M=2.1,Az=237';=50; timePosle=130;=-120; z2=-90;'12'='2014208075635320_POLE__1_2.atr';='2014208075635320_POLE__1_3.atr';='2014208085635320_POLE__1_2.atr';='2014208085635320_POLE__1_3.atr';=datenum(2014,07,27,07,56,35.320);=datenum(2014,07,27,08,24,25);='M=1.7,Az=217';=50; timePosle=115;
z1=-125; z2=-95;'13'='2014208125635320_POLE__1_2.atr';='2014208125635320_POLE__1_3.atr';='2014208135635320_POLE__1_2.atr';='2014208135635320_POLE__1_3.atr';=datenum(2014,07,27,12,56,35.320);=datenum(2014,07,27,13,12,9);='M=4.0,Az=329';=50; timePosle=300;=-125; z2=-55;'14'='2014210225635320_POLE__1_2.atr';='2014210225635320_POLE__1_3.atr';='2014210235635320_POLE__1_2.atr';='2014210235635320_POLE__1_3.atr';=datenum(2014,07,29,22,56,35.320);=datenum(2014,07,29,23,32,35);='M=2.1,Az=269';=50; timePosle=150;=-125; z2=-90;'15'='2014210225635320_POLE__1_2.atr';='2014210225635320_POLE__1_3.atr';='2014210235635320_POLE__1_2.atr';='2014210235635320_POLE__1_3.atr';=datenum(2014,07,29,22,56,35.320);=datenum(2014,07,29,23,35,26);
Продолжение приложения Б='M=1.9,Az=247';=50; timePosle=130;=-125; z2=-95;'16'='2014211005635320_POLE__1_2.atr';='2014211005635320_POLE__1_3.atr';='2014211015635320_POLE__1_2.atr';='2014211015635320_POLE__1_3.atr';=datenum(2014,07,30,00,56,35.320);=datenum(2014,07,30,01,47,28);='M=1.8,Az=217';=50; timePosle=120;=-125; z2=-85;'17'='2014211075635320_POLE__1_2.atr';='2014211075635320_POLE__1_3.atr';='2014211085635320_POLE__1_2.atr';='2014211085635320_POLE__1_3.atr';=datenum(2014,07,30,07,56,35.320);=datenum(2014,07,30,08,20,48);='M=1.8,Az=099';=50; timePosle=130;=-125; z2=-80;'18'='2014212005635320_POLE__1_2.atr';='2014212005635320_POLE__1_3.atr';='2014212015635320_POLE__1_2.atr';='2014212015635320_POLE__1_3.atr';=datenum(2014,07,31,00,56,35.320);
Продолжение приложения Б=datenum(2014,07,31,01,54,13);='M=1.7,Az=237';=50; timePosle=110;=-125; z2=-90;'19'='2014212025635320_POLE__1_2.atr';='2014212025635320_POLE__1_3.atr';='2014212035635320_POLE__1_2.atr';='2014212035635320_POLE__1_3.atr';=datenum(2014,07,31,02,56,35.320);=datenum(2014,07,31,03,43,5);='M=2.7,Az=217';=50; timePosle=120;=-125; z2=-90;'20'='2014213035956320_POLE__1_2.atr';='2014213035956320_POLE__1_3.atr';='2014213045956320_POLE__1_2.atr';='2014213045956320_POLE__1_3.atr';=datenum(2014,08,01,03,59,56.320);=datenum(2014,08,01,04,46,0);='M=2.1,Az=262';=50; timePosle=150;=-125; z2=-75;'21'='2014213035956320_POLE__1_2.atr';='2014213035956320_POLE__1_3.atr';='2014213045956320_POLE__1_2.atr';='2014213045956320_POLE__1_3.atr';=datenum(2014,08,01,03,59,56.320);
begtime=datenum(2014,08,01,04,49,29);='M=1.4,Az=266';=30; timePosle=70;=-125; z2=-90;'22'='2014213035956320_POLE__1_2.atr';='2014213035956320_POLE__1_3.atr';='2014213045956320_POLE__1_2.atr';='2014213045956320_POLE__1_3.atr';=datenum(2014,08,01,03,59,56.320);=datenum(2014,08,01,04,51,32);='M=0.8,Az=262';=30; timePosle=70;=-125; z2=-100;=.3; dt=3;=' ';=9; L=86400; fs=200;=1.583e-6/2001.67*1000;
% load= importdata([p filenameX1], delimeter, nStrok);=data666.data;= importdata([p filenameX2], delimeter, nStrok);=[dataX; data666.data]*m;= importdata([p filenameY1], delimeter, nStrok);=data666.data;= importdata([p filenameY2], delimeter, nStrok);=[dataY; data666.data]*m;
Продолжение приложения Б;
% filter
[b,a]=butter(3,fmin/fs*2,'high');=filter(b,a,dataX-mean(dataX));=filter(b,a,dataY-mean(dataY));
% cut=round((begtime-timeDo/L-filetime)*L*fs);=round((begtime+timePosle/L-filetime)*L*fs);=dataXf(t1:t2); Y=dataYf(t1:t2);
% spectrogram(X,fs,1000,900,0,z1,z2), ylim([1,30]);(gcf,'-dpng',[p2 nomer '_' ZTname '_X.png'],'-r300')(gcf);(Y,fs,1000,900,0,z1,z2), ylim([1,30]);(gcf,'-dpng',[p2 nomer '_' ZTname '_Y.png'],'-r300')(gcf);sp01(fname,fs,wl,overlap, meaning,zmin,zmax)
%load file(fname)=' ';=9;= importdata(fname, delimeter, nStrok);=data666.data;data666;=fname;
Продолжение приложения Б
% computingmeaning>1=AntiTrendFast(data1,meaning);
[~,F,T,P]=spectrogram(data1,hann(wl),overlap,wl,fs);(1:3,:)=[];(1:3,:)=[];
%plot 1st;('position',[0.04 0.75 0.94
0.22]);(gca,'fontSize',9)((1:length(data1))./fs,data1); axis tight; grid
on;=get(gca,'ylim');(1)=aaa(1)-0.02*(aaa(2)-aaa(1));(2)=aaa(2)+0.02*(aaa(2)-aaa(1));(aaa);
%plot 2nd('position',[0.04 0.05 0.94 0.62]);(gca,'fontSize',9)
(T,F,10*log10(P),'edgecolor','none');=0; aa2=length(data1)/fs; aa3=F(1,1); aa22=F(size(F)); aa4=aa22(1);([aa1 aa2 aa3 aa4]);(gca,'yscale','log');('east');(jet(4096));(gca,'clim',[zminzmax]);;
a=[.0001 .0001 .0002 .0002 .0003 .0003 .0004 .0004 .0005 .0005 .0006 .0006 .0007 .0007 .0008 .0008 .0009 .0009 ...
.001 .001 .002 .002 .003 .003 .004 .004 .005 .005 .006 .006 .007 .007 .008 .008 .009 .009 ...
.01 .01 .02 .02 .03 .03 .04 .04 .05 .05 .06 .06 .07 .07 .08 .08 .09 .09 ...%9*2*6=108
.1 .1 .2 .2 .3 .3 .4 .4 .5 .5 .6 .6 .7 .7 .8 .8 .9 .9 ...
1 2 2 3 3 4 4 5 5 6 6 7 7 8 8 9 9 ...
10 20 20 30 30 40 40 50 50 60 60 70 70 80 80 90 90];
b=[aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 ...
aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 ...
aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1];
c=3000*ones(1,108);
ПРИЛОЖЕНИЕ В
Программа в система matlab
для построения спектрограмм по трем каналам
%% параметры
% пути к файлам
%p='D:\_____\';% для ноута 11''
p='s:\ATR\';%дляПК
='2014206055635320_POLE__1_2.atr';='2014206055635320_POLE__1_2.atr';
% время=datenum(2014,07,25,05,56,35.320);=datenum(2014,07,25,06,07,02.125);
dt=4;
% фильтр=3;=20;
% прочее=' ';=9;=86400;=200;
%% загружаем
% 1-й= importdata([p filenameX], delimeter, nStrok);=data666.data;
% 2-й
data666 = importdata([p filenameY], delimeter, nStrok);=data666.data;data666;
%% фильтруем
[b,a]=butter(3,[fminfmax]/fs*2);=filter(b,a,dataX-mean(dataX));=filter(b,a,dataY-mean(dataY));
%% вырезаем
% EQtime - времяземлетрясения
% filetime - времяначалафайла
% они даны в сутках, например 735812.1666240741
% L - кол-во секунд в сутках, =86400
% dt - длина интересующей нас записи, с
dataXfC=dataXf(((EQtime-filetime)*L*fs-dt*fs/2):((EQtime-filetime)*L*fs+dt*fs/2));=dataYf(((EQtime-filetime)*L*fs-dt*fs/2):((EQtime-filetime)*L*fs+dt*fs/2));
%% ищемразброс=max((dataXfC.^2+dataYfC.^2).^.5);
%% вращаем, рисуем=99999999999999999999999;=0;(1);=1:18=(i-1)*10;
[ X1,Y1 ] = Func_rotate(dataXfC,dataYfC,alf );('position',[.05 .05*(i+.5) .2 .049]);
plot(1:length(X1),X1)on, axis tight, ylim([-r r]);(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.3 .05*(i+.5) .2 .049]);(1:length(Y1),Y1)on, axis tight, ylim([-r r]);(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.55 .05*(i+.5) .2 .049]);(1:length(X1),X1.^3)on, axis tight, ylim([-r.^3 r.^3]);(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.8 .05*(i+.5) .2 .049]);(1:length(Y1),Y1.^3)on, axis tight, ylim([-r.^3 r.^3]);(num2str(alf));>1, set(gca,'xticklabel',''), end
(Y1)<am=std(Y1);=alf;(Y1)>amax
Продолжение приложения В
amax=std(Y1);
iii=alf;(gcf,'units','normalized','position',[0.15 0.15 .8 .75]);(gcf,'-dpng',[p filenameX '.png'],'-r300')
%close(gcf)
% prompt = {'alfa'};
% dlg_title = 'Чемуравноalfa?';
% num_lines = 1;
% def = {num2str(iii)};
% answer = inputdlg(prompt,dlg_title,num_lines,def);
%% вращаем, рисуем ещё раз=99999999999999999999999;
amax=0;(2);=1:18=iii-18+(i-1)*2;
[ X1,Y1 ] = Func_rotate(dataXfC,dataYfC,alf );
('position',[.05 .05*(i+.5) .2 .049]);(1:length(X1),X1)on, axis tight, ylim([-r r]);(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.3 .05*(i+.5) .2 .049]);(1:length(Y1),Y1)
Продолжение приложения В
gridon, axistight, ylim([-rr]);
ylabel(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.55 .05*(i+.5) .2 .049]);(1:length(X1),X1.^3)on, axis tight, ylim([-r.^3 r.^3]);(num2str(alf));>1, set(gca,'xticklabel',''), end
('position',[.8 .05*(i+.5) .2 .049]);(1:length(Y1),Y1.^3)on, axis tight, ylim([-r.^3 r.^3]);(num2str(alf));>1, set(gca,'xticklabel',''), end
(Y1)<am=std(Y1);=alf;(Y1)>amax=std(Y1);=alf;(gcf,'units','normalized','position',[0.15 0.15 .8 .75]);(num2str(iii))(gcf,'-dpng',[p filenameX '__.png'],'-r300')
%close(gcf)
Продолжение приложения В[ X1,Y1 ] = Func_rotate( X,Y,alfa )=alfa*pi/180;=X*cos(a)-Y*sin(a);=X*sin(a)+Y*cos(a);
ПРИЛОЖЕНИЕ Г
Программа для расчета азимута землетрясений, разработанная в системе matlab
function sp01(fname,fs,wl,overlap, meaning,zmin,zmax)
%load file=' ';=9;= importdata(fname, delimeter,
nStrok);=data666.data;data666;
% computingmeaning>1=AntiTrendFast(data1,meaning);
[~,F,T,P]=spectrogram(data1,wl,overlap,wl,fs);(1:3,:)=[];(1:3,:)=[];
%plot 1st;('position',[0.04 0.75 0.94 0.22]);(gca,'fontSize',9)((1:length(data1))./fs,data1); axis tight; grid on;=get(gca,'ylim');(1)=aaa(1)-0.02*(aaa(2)-aaa(1));(2)=aaa(2)+0.02*(aaa(2)-aaa(1));(aaa);
%plot 2nd('position',[0.04 0.05 0.94 0.62]);(gca,'fontSize',9)
(T,F,10*log10(P),'edgecolor','none');=0; aa2=length(data1)/fs; aa3=F(1,1); aa22=F(size(F)); aa4=aa22(1);([aa1 aa2 aa3 aa4]);(gca,'yscale','log');('east');(jet(4096));
(gca,'clim',[zminzmax]);
;
a=[.0001 .0001 .0002 .0002 .0003 .0003 .0004 .0004 .0005 .0005 .0006 .0006 .0007 .0007 .0008 .0008 .0009 .0009 ...
.001 .001 .002 .002 .003 .003 .004 .004 .005 .005 .006 .006 .007 .007 .008 .008 .009 .009 ...
.01 .01 .02 .02 .03 .03 .04 .04 .05 .05 .06 .06 .07 .07 .08 .08 .09 .09 ...%9*2*6=108
.1 .1 .2 .2 .3 .3 .4 .4 .5 .5 .6 .6 .7 .7 .8 .8 .9 .9 ...
1 2 2 3 3 4 4 5 5 6 6 7 7 8 8 9 9 ...
10 20 20 30 30 40 40 50 50 60 60 70 70 80 80 90 90];
b=[aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 ...
aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 ...
aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 ...
aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1 aa1-1 aa2+1 aa2+1 aa1-1];
c=3000*ones(1,108);(b,a,c,':b'); hold off;