331 lines
16 KiB
Matlab
331 lines
16 KiB
Matlab
function loaddata_edr(sattelite, year, month, type_sensor) %DATA=loaddata_edr(filename)
|
||
%% By L_DelOff
|
||
% Загрузка данных из файла *.EDR
|
||
% Описание файла приведено в:
|
||
% !docs/!AFRL ASCII and Binary File Format Descriptions.pdf
|
||
%% Поехали
|
||
warning off
|
||
%% Скачать данные в формате архива
|
||
download_fun(sattelite, year, month, type_sensor); %скачивание и обработка нужных данных
|
||
%% Распаковать архив
|
||
extract_data_fun(sattelite, year, month, type_sensor); %скачивание и обработка нужных данных
|
||
%% Преобраховать EDR в mat
|
||
convert_data_fun(sattelite, year, month, type_sensor); %скачивание и обработка нужных данных
|
||
end
|
||
|
||
function convert_data_fun(sattelite, year, month, type_sensor)
|
||
foldername=['DATA/f' num2str(sattelite) '/' type_sensor '/' num2str_new(year,4) '/' num2str_new(month,2)];
|
||
target=ls(foldername);
|
||
[str row]=size(target);
|
||
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)>4
|
||
if filename(end-3:end)=='.EDR'
|
||
message=['Расшифровка файла:',filename,'\n'];
|
||
fprintf(message);
|
||
|
||
fid = fopen([foldername,'/',filename], 'rt'); % открыть для чтения
|
||
counter=1;
|
||
while ~feof(fid) % пока не конец файла
|
||
frame=[]; % здесь будет 1 кадр (1 минута) в текстовом формате
|
||
for i=1:114
|
||
if ~feof(fid)
|
||
frame(i).str = fgets(fid); % Последовательно читаем из файла строки
|
||
end
|
||
end
|
||
|
||
if counter==1
|
||
DATA=getframe_edr(frame); % Считываем кадр
|
||
else
|
||
DATA(counter)=getframe_edr(frame); % Считываем кадр
|
||
end
|
||
counter=counter+1;
|
||
if mod(counter,100)==0
|
||
disp(num2str(counter));
|
||
end
|
||
end
|
||
|
||
save([foldername,'/',filename(1:end-4),'.mat'],'DATA');
|
||
|
||
message=['Расшифровка файла закончена\n'];
|
||
fprintf(message);
|
||
end
|
||
end
|
||
end
|
||
end
|
||
|
||
|
||
|
||
function download_fun(sattelite, year, month, type_sensor)
|
||
% https://satdat.ngdc.noaa.gov/dmsp/data/f18/ssies/2020/05/
|
||
% https://satdat.ngdc.noaa.gov/dmsp/data/f18/ssies/2020/05/PS.CKGWC_SC.U_DI.A_GP.SIES3-F18-R99990-B9999090-APGA_AR.GLOBAL_DD.20200501_TP.000001-235959_DF.EDR.gz
|
||
message=['Ищем данные спутника F' num2str(sattelite) ' по дате:' num2str_new(year,4) '/' num2str_new(month,2) '...'];
|
||
fprintf(message)
|
||
url_site=['https://satdat.ngdc.noaa.gov/dmsp/data/f' num2str_new(sattelite,2) '/' type_sensor '/' num2str_new(year,4) '/' num2str_new(month,2) '/'];
|
||
fullList = webread(url_site);
|
||
message='готово\n';
|
||
fprintf(message)
|
||
url=extract_url(fullList);
|
||
|
||
for i=1:length(url)
|
||
message=['число ' num2str(i) '...'];
|
||
fprintf(message);
|
||
url_file = [url_site url(i).file];
|
||
filename = ['DATA/f' num2str(sattelite) '/' type_sensor '/' num2str_new(year,4) '/' num2str_new(month,2) '/' url(i).file];
|
||
mkdir(['DATA/f' num2str(sattelite) '/' type_sensor '/' num2str_new(year,4) '/' num2str_new(month,2)])
|
||
if ~exist(filename,'file')
|
||
outfilename = websave(filename,url_file);
|
||
message='скачано\n';
|
||
else
|
||
message='уже есть\n';
|
||
end
|
||
fprintf(message);
|
||
end
|
||
|
||
a=1;
|
||
end
|
||
|
||
function url=extract_url(text)
|
||
counter=0;
|
||
i=1;
|
||
while i<length(text)
|
||
if text(i)=='"'
|
||
a=i;
|
||
i=i+1;
|
||
while text(i)~='"'
|
||
i=i+1;
|
||
end
|
||
b=i;
|
||
if b-a>10
|
||
c=text(b-6:b-1);
|
||
if c=='EDR.gz'
|
||
counter=counter+1;
|
||
url(counter).file=text(a+1:b-1);
|
||
end
|
||
end
|
||
end
|
||
i=i+1;
|
||
end
|
||
end
|
||
|
||
function FRAME=getframe_edr(frame)
|
||
HEALTH=0; % "здоровье кадра"
|
||
%% строка 1: Пустая строка (Blank line)
|
||
if isspace(frame(1).str) HEALTH=HEALTH+1;
|
||
end
|
||
%% строка 2: Шапка + данные о ПО(если это первый кадр в файле)
|
||
if frame(2).str~="" HEALTH=HEALTH+1; end
|
||
%% строка 3: 5 значений
|
||
z = textscan(frame(3).str,'%f','EmptyValue',-Inf);
|
||
FRAME.N_record_bin = z{1}(1); %1. Номер записи в двоичном файле, откуда был взят этот EDR.(Number of the record in the binary file from where this EDR was taken.)
|
||
FRAME.N_EDR_bin = z{1}(2); %2. Номер EDR в записи в двоичном файле (1-3).(Number of the EDR within the record in the binary file (1-3).)
|
||
FRAME.SatID = z{1}(3); %3. ID спутника.(Satellite Flight ID (two digit integer))
|
||
FRAME.year = fix(z{1}(4)/1e4); %4. Год, месяц, день, час, минута
|
||
FRAME.month = mod(fix(z{1}(4)/1e2),1e2);
|
||
FRAME.day = mod(z{1}(4),1e2);
|
||
FRAME.hour = fix(z{1}(5)/1e2);
|
||
FRAME.minute = mod(z{1}(5),1e2);
|
||
if length(z{1})==5 HEALTH=HEALTH+1; end
|
||
%% строка 4: ЭФЕМЕРИДА (EPHEMERIS)
|
||
if frame(4).str(1:end-1)=="EPHEMERIS" HEALTH=HEALTH+1; end
|
||
%% строки 5-7: по 6 значений
|
||
for i=0:2
|
||
a = textscan(frame(5+i).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
% данные на секунде 00 20 40 (i=0 1 2)
|
||
FRAME.GeoLat(i+1)=a{1}(1); % 1. Географическая широта (градусы, север)(Geographic latitude (degrees, north))
|
||
FRAME.GeoLong(i+1)=a{1}(2); % 2. Географическая долгота (градусы, восток)(Geographic longitude (degrees, east))
|
||
FRAME.ApLat(i+1)=a{1}(3); % 3. Вершинная широта (градусы, север)(Apex latitude (degrees, north))
|
||
FRAME.ApLong(i+1)=a{1}(4); % 4. Вершинная долгота (градусы, восток)(Apex longitude (degrees, east))
|
||
FRAME.ApLocalTime(i+1)=a{1}(5);% 5. Местное время Apex (часы)(градусы, восток)(Apex local time (hours))
|
||
FRAME.SatAlt(i+1)=a{1}(6); % 6. Высота (км)(Satellite altitude (km))
|
||
end
|
||
%% строка 8: Потенциал спутника (“SATELLITE POTENTIAL, LAST = SOURCE”)
|
||
if frame(8).str(1:end-1)=="SATELLITE POTENTIAL, LAST = SOURCE" HEALTH=HEALTH+1; end
|
||
%% строки 9-10: 16 значений (2х8)
|
||
a = textscan(frame(9).str,'%f','EmptyValue',-Inf);
|
||
b = textscan(frame(10).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==8 HEALTH=HEALTH+1; end
|
||
if length(b{1})==8 HEALTH=HEALTH+1; end
|
||
FRAME.SatPot(1:8)=a{1}(1:8); % значения для потенциала спутника (Vbias + VIP) в вольтах каждые 4 секунды
|
||
FRAME.SatPot(9:15)=b{1}(1:7); % values for the satellite potential (Vbias+VIP) in volts every 4 seconds
|
||
FRAME.SatPotSrc=b{1}(8); % Измеритель потенциала (an integer (1-2) indicating the satellite potential source where:
|
||
% 1 - as set by the on-board microprocessor
|
||
% 2 - as set by the SENPOT sensor)
|
||
%% строка 11: Потенциал спутника (“PRIMARY PLASMA DENSITY, THEN SOURCE”)
|
||
if frame(11).str(1:end-1)=="PRIMARY PLASMA DENSITY, THEN SOURCE" HEALTH=HEALTH+1; end
|
||
%% строки 12-21: 60 значений (6х10)
|
||
% Односекундные средние значения плотности первичной плазмы (# / см3)
|
||
% One-second averages of the primary plasma density (#/cm3)
|
||
z=[];
|
||
for i=12:21
|
||
a = textscan(frame(i).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
end
|
||
FRAME.PrimPlasmaDen=z;
|
||
%% строка 22: Источник данных по плазме
|
||
% An integer (1-3) indicating the plasma density source where:
|
||
% 1 - Ion density from SM sensor
|
||
% 2 - Ion density from DM sensor
|
||
% 3 - Electron density from EP sensor (DC Mode)
|
||
a = textscan(frame(22).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==1 HEALTH=HEALTH+1; end
|
||
FRAME.PrimPlasmaDenSrc=a{1};
|
||
%% строка 23: Горизонтальный ионный дрейф (“HORIZONTAL ION DRIFT VELOCS”)
|
||
if frame(23).str(1:end-1)=="HORIZONTAL ION DRIFT VELOCS" HEALTH=HEALTH+1; end
|
||
%% строки 24-33: 60 значений (6х10)
|
||
% Односекундные значения скорости горизонтального дрейфа (м/с) для времен от HHMM00 до HHMM59 с 6 значениями в строке и всего 60 значениями.
|
||
% One-second values for the horizontal drift speed (m/s) for times HHMM00 through HHMM59 with 6 values per line and 60 values total.
|
||
z=[];
|
||
for i=24:33
|
||
a = textscan(frame(i).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
end
|
||
FRAME.HorIonDrift=z;
|
||
%% строка 34: Вертикальный ионный дрейф (“VERTICAL ION DRIFT VELOCS”)
|
||
if frame(34).str(1:end-1)=="VERTICAL ION DRIFT VELOCS" HEALTH=HEALTH+1; end
|
||
%% строки 35-44: 60 значений (6х10)
|
||
% Односекундные значения скорости вертикального дрейфа (м/с) для времен от HHMM00 до HHMM59 с 6 значениями в строке и всего 60 значениями.
|
||
% One-second values for the vertical drift speed (m/s) for times HHMM00 through HHMM59 with 6 values per line and 60 values total.
|
||
z=[];
|
||
for i=35:44
|
||
a = textscan(frame(i).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
end
|
||
FRAME.VerIonDrift=z;
|
||
%% строка 45: (“CKL ANALYSES, THEN SOURCE”)
|
||
if frame(45).str(1:end-1)=="CKL ANALYSES, THEN SOURCE" HEALTH=HEALTH+1; end
|
||
%% строка 46-63:
|
||
for i=1:6 % для секунд: 5 15 25 35 45 55
|
||
z=[];
|
||
a = textscan(frame(46+(3*(i-1))).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==4 HEALTH=HEALTH+1; end
|
||
FRAME.CKL_RMS(i)=a{1}(1); % CKL Analysis: (RMS ΔN)/N (%)
|
||
FRAME.CKL_T1(i)=a{1}(2); % CKL Analysis: T1
|
||
FRAME.CKL_p1(i)=a{1}(3); % CKL Analysis: p1
|
||
FRAME.CKL_CKL(i)=a{1}(4); % CKL Analysis: CKL
|
||
a = textscan(frame(47+(3*(i-1))).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==8 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
a = textscan(frame(48+(3*(i-1))).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==8 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
FRAME.CKL_PDS(i,:)=z(1:end-1); %CKL Analysis: Decimated power density spectrum (PDS) for time period.
|
||
end
|
||
%% строка 64: Источник данных по CLK
|
||
% An integer (1-3) indicating data source used for CKL calculation where:
|
||
% 1 - SM density data only
|
||
% 2 - SM density and filter data
|
||
% 3 - EP DC mode density data
|
||
a = textscan(frame(64).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==1 HEALTH=HEALTH+1; end
|
||
FRAME.CLK_Src=a{1};
|
||
%% строка 112 (потом понадобится)
|
||
a = textscan(frame(112).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==7 HEALTH=HEALTH+1; end
|
||
FRAME.ADC_temp=a{1}(2); %2: ADC temperature (degrees C)
|
||
FRAME.SEP_temp=a{1}(3); %3: SEP temperature
|
||
FRAME.DM_offsetV=a{1}(4); %4: DM offset voltage (volts) if IES-2 or RPA plasma plate potential (volts) if SSIES-3
|
||
FRAME.DM_mode=a{1}(5); %5: DM mode (integer)
|
||
% SSIES-2:
|
||
% -1: Undefined
|
||
% 0: Normal
|
||
% 1-8: H+
|
||
% 9: FIBA
|
||
% SSIES-3:
|
||
% 0: Slow
|
||
% 1: Normal
|
||
FRAME.EP_mode=a{1}(6); % 6: EP mode (0-6 : A,B,BS,C,D,DS,E) (integer)
|
||
FRAME.EP_VIP=a{1}(7); % 7: VIP at EDR start (volts)
|
||
%% строки 65-80
|
||
if (FRAME.EP_mode==0)||(FRAME.EP_mode==1)||(FRAME.EP_mode==2)||(FRAME.EP_mode==6)
|
||
% строка 65 (“EP SWEEP ANALYSES SETS”)
|
||
if frame(65).str(1:end-1)=="EP SWEEP ANALYSES SETS" HEALTH=HEALTH+1; end
|
||
% EP Sweep analyses (every 4 seconds). There are 15 EP sweep analysis
|
||
% sets. Each is valid for either 4 (modes A, B and BS) or 2 (mode E)
|
||
% seconds centered on the time specified in the set. Each analysis set
|
||
% contains the following parameters:
|
||
for i=1:15
|
||
a = textscan(frame(66+(i-1)).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
FRAME.EP_SCT(i)=a{1}(1); % 1: Sweep center time (UT, seconds) (integer)
|
||
FRAME.EP_Eden(i)=a{1}(2); % 2: Electron density (el/cm3)
|
||
FRAME.EP_Etemp(i)=a{1}(3); % 3: Electron temperature (degrees K)
|
||
FRAME.EP_SatPot(i)=a{1}(4); % 4: Satellite potential (volts)
|
||
FRAME.EP_Qualif(i)=a{1}(5); % 5: Analysis qualifier (integer). Set to zero if the on-board microprocessor did not perform the analysis, per the flag in element 415. If the on-board microprocessor was in use, then it is set to the MP EP flags word from Word 11 of Cycle 1 in the telemetry. This word can also be zero.
|
||
FRAME.EP_PES(i)=a{1}(6); % 6: EP photo-electron surrogate value.
|
||
end
|
||
end
|
||
%% строка 81 (“EP ANALYSES SOURCE”)
|
||
if frame(81).str(1:end-1)=="EP ANALYSES SOURCE" HEALTH=HEALTH+1; end
|
||
%% строка 82
|
||
% An Integer (1-2) indicating EP analysis source where
|
||
% 1 – Ground processing analysis
|
||
% 2 – On-board microprocessor analysis
|
||
a = textscan(frame(82).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==1 HEALTH=HEALTH+1; end
|
||
FRAME.EP_src=a{1};
|
||
%% строка 83 (“RPA SWEEP ANALYSES SETS, THEN SOURCE”)
|
||
if frame(83).str(1:end-1)=="RPA SWEEP ANALYSES SETS, THEN SOURCE" HEALTH=HEALTH+1; end
|
||
%% строка 84-98 (“RPA SWEEP ANALYSES SETS, THEN SOURCE”)
|
||
% RPA Sweep analyses (every 4 seconds). There are 15 RPA sweep analysis
|
||
% sets. Each is valid for the 4 seconds centered on the time specified in
|
||
% the set. Each analysis set contains the following parameters:
|
||
for i=1:15
|
||
a = textscan(frame(84+(i-1)).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==8 HEALTH=HEALTH+1; end
|
||
FRAME.RPA_SCT(i)=a{1}(1); % 1: Sweep center time (UT, seconds) (integer)
|
||
FRAME.RPA_O2den(i)=a{1}(2); % 2: O+ density (ion/cm3)
|
||
FRAME.RPA_HHeden(i)=a{1}(3); % 3: Total (H+ + He+) density (ion/cm3)
|
||
FRAME.RPA_LIF(i)=a{1}(4); % 4: Light ion flag (integer)
|
||
% 0 - No light ion
|
||
% 1 - Light ion is H+
|
||
% 2 - Light ion is He+
|
||
% 3± = 3 + 10000 x (H+ fraction)
|
||
FRAME.RPA_IonNemp(i)=a{1}(5); % 5: Ion temperature (degrees K)
|
||
FRAME.RPA_RamIonDrift(i)=a{1}(6);% 6: Ram ion drift velocity (m/s)
|
||
FRAME.RPA_qualif(i)=a{1}(7); % 7: Analysis qualifier (integer)
|
||
% 0 - Analysis terminated unsuccessfully
|
||
% 1 - Successful analysis
|
||
FRAME.RPA_totIonDens(i)=a{1}(8);% 8: RPA-derived total ion density
|
||
end
|
||
%% строка 99
|
||
% An integer (1-2) indicating RPA analysis source where:
|
||
% 1 - Ground processing analysis
|
||
% 2 - On-board microprocessor analysis
|
||
a = textscan(frame(99).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==1 HEALTH=HEALTH+1; end
|
||
FRAME.RPA_src=a{1};
|
||
%% строка 100 (“DM ION DENSITY”)
|
||
if frame(100).str(1:end-1)=="DM ION DENSITY" HEALTH=HEALTH+1; end
|
||
%% строки 101-110: 60 значений (6х10)
|
||
%DM one-second average ion density (ion/cm3) with 6 values per line and 60
|
||
% values total
|
||
z=[];
|
||
for i=101:110
|
||
a = textscan(frame(i).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})==6 HEALTH=HEALTH+1; end
|
||
z=[z a{1}.'];
|
||
end
|
||
FRAME.DM_IonDen=z;
|
||
%% строка 111 (“ENGINEERING DATA”)
|
||
if frame(111).str(1:end-1)=="ENGINEERING DATA" HEALTH=HEALTH+1; end
|
||
%% строка 113 (“FILLER”)
|
||
if frame(113).str(1:end-1)=="FILLER" HEALTH=HEALTH+1; end
|
||
%% строка 114 (“FILLER”)
|
||
a = textscan(frame(114).str,'%f','EmptyValue',-Inf);
|
||
if length(a{1})>=1 HEALTH=HEALTH+1; end
|
||
FRAME.Filler=a{1};
|
||
|
||
FRAME.HEALTH=HEALTH;
|
||
end |