Earthquake/analysis.m

137 lines
4.7 KiB
Matlab
Raw Blame History

This file contains ambiguous Unicode characters!

This file contains ambiguous Unicode characters that may be confused with others in your current locale. If your use case is intentional and legitimate, you can safely ignore this warning. Use the Escape button to highlight these characters.

function param=analysis(param)
%% By L_DelOff
% собираю данные по землетрясению
%% Поехали
type_sensor=param.type_sensor;
DATA.DATA=[];
%% Объединяю данные в диапазоне
if 1
for i=1:length(param.dates)
year=param.dates(i).year;
month=param.dates(i).month;
day=param.dates(i).day;
DATA_temp=load_file(param.sattelite_range, year, month, day,...
type_sensor); %загрузка данных
DATA.DATA=[DATA.DATA DATA_temp.DATA];
end
end
%% Сортировка данных
R_quake=param.R_quake; % радиус действия относительно магнитуды,км (еще тыщу ставил)
% расположение эпицентра
fiA=param.fiA; % широта
LA=param.LA; % долгота
fname_rep=param.filename_report;
if exist(fname_rep,'file')
load(fname_rep,'report'); % тут будут все данные для графика
else
%report=[];
%report.data=[];
report.data.N_sat=[];
report.data.R=[];
end
for i=1:length(param.dates)
HHeden=[];
O2den=[];
Eden=[];
Etemp=[];
for j=1:length(DATA.DATA)
date_true=1;
if DATA.DATA(j).year~=param.dates(i).year % сортировка по году
date_true=0;
end
if DATA.DATA(j).month~=param.dates(i).month % по месяцу
date_true=0;
end
if (DATA.DATA(j).day~=param.dates(i).day)% по числам
date_true=0;
end
if date_true
for k=1:3 % координаты спутника только 3 раза в минуту(кадр) записываются
fiB=DATA.DATA(j).GeoLat(k); % 1. Географическая широта (градусы, север)(Geographic latitude (degrees, north))
LB=DATA.DATA(j).GeoLong(k); % 2. Географическая долгота (градусы, восток)(Geographic longitude (degrees, east))
R=mydistance(fiA,LA,fiB,LB);
if R<=R_quake
HHeden(end+1)= mymean(DATA.DATA(j).RPA_HHeden((k-1)*5+1:(k-1)*5+5),1e10);
O2den(end+1) = mymean(DATA.DATA(j).RPA_O2den((k-1)*5+1:(k-1)*5+5),1e10);
Eden(end+1) = mymean(DATA.DATA(j).EP_Eden((k-1)*5+1:(k-1)*5+5),1e10);
Etemp(end+1) = mymean(DATA.DATA(j).EP_Etemp((k-1)*5+1:(k-1)*5+5),1e10);
end
end
end
end
num=give_number(report,param.sattelite_range,R_quake);
report.data(num).N_sat=param.sattelite_range;
report.data(num).R=R_quake;
report.data(num).date{i}=[myformat(param.dates(i).day,2),'.',...
myformat(param.dates(i).month,2),'.',...
myformat(param.dates(i).year,4)];
report.data(num).RPA_HHeden(i)=mymean(HHeden,1e10);
report.data(num).RPA_O2den(i)=mymean(O2den,1e10);
report.data(num).EP_Eden(i)=mymean(Eden,1e10);
report.data(num).EP_Etemp(i)=mymean(Etemp,1e10);
end
param.report=report;
save(fname_rep,'report');
end
function num=give_number(report,f,R)
flag=0; % Проверяет, есть ли в отчете данные со спутника f
[~,b]=size(report.data);
for i=1:b
if report.data(i).N_sat==f
if report.data(i).R==R
flag=1;
num=i;
end
end
end
if flag==0
num=b+1;
end
end
function DATA=load_file(sattelite, year, month, day, type_sensor)
success=0;
foldername=['DATA/f' num2str(sattelite) '/' type_sensor '/' num2str_new(year,4) '/' num2str_new(month,2)];
target=ls(foldername);
[str row]=size(target);
DATA=[];
for i=1:str
filename=target(i,:);%end-3:end);
counter=0;
for j=1:length(filename)
if filename(j)==' '
counter=counter+1;
end
end
filename=filename(1:end-counter);
if length(filename)>7
keyword=[num2str_new(year,4) num2str_new(month,2) num2str_new(day,2)];
flag=0;
for j=1:length(filename)-length(keyword)
if filename(j:j+length(keyword)-1)==keyword
flag=1;
end
end
if (prod(filename(end-3:end)=='.mat')==1)&&(flag==1)
disp(['Найден подходящий файл: ',filename])
DATA = load([foldername,'/',filename], 'DATA'); % открыть для чтения
success=1;
end
end
end
if ~success
disp('Нужных файлов не найдено, перехожу к загрузке')
loaddata_edr(sattelite, year, month, type_sensor);
DATA=load_file(sattelite, year, month, day, type_sensor);
end
end