publicImage/gprMax
1
0
派生 0
镜像自地址 https://gitee.com/sunhf/gprMax.git 已同步 2025-08-07 04:56:51 +08:00
代码 工单 软件包 项目 版本发布 百科 动态

Added initial MATLAB script to convert gprMax output to some common GPR formats (RD3, DZT, DT1).

浏览代码
这个提交包含在：
Craig Warren
2018-04-04 16:31:01 +01:00
父节点 1acd8e461a
当前提交 57520e0865
共有 1 个文件被更改，包括 204 次插入 和 189 次删除
下载 Patch 文件 下载 Diff 文件 展开所有文件 折叠所有文件

101
tools/MATLAB_scripts/gprMax_output_converter.m → tools/MATLAB_scripts/outputfile_converter.m
查看文件

@@ -1,14 +1,15 @@
%gprMax_output_converter, Angelis Dimitrios 2017, 2018
%Converts the gprMax merged output HDF5 file to RD3, DZT, DT1 files.
% outputfile_converter.m - converts gprMax merged output HDF5 file to RD3, DZT, DT1 files.
%
% Author: Dimitrios Angelis
% Copyright: 2017-2018
% Last modified: 25/03/2018
%Last Modified 25/03/2018
clear, clc, close all;
clear; close all; clc;
%Select file. Keep file name, path name and file extension ================
[infile, path] = uigetfile('*.out', 'Select gprMax File', 'Multiselect', 'Off');
% Select file. Keep file name, path name and file extension ==============
[infile, path] = uigetfile('*.out', 'Select gprMax output file', 'Multiselect', 'Off');
if isequal(infile, 0)
erd = errordlg('No Input File');
erd = errordlg('No output file');
but = findobj(erd, 'Style', 'Pushbutton');
delete(but);
pause(1); close(erd);
@@ -20,13 +21,13 @@ HDR.pname = path;
HDR.fext = 'out';
%Read data from HDF5 file =================================================
% Read data from HDF5 file ===============================================
dataex = h5read(infile, '/rxs/rx1/Ex');
dataey = h5read(infile, '/rxs/rx1/Ey');
dataez = h5read(infile, '/rxs/rx1/Ez');
%Field check ==============================================================
% Field check ============================================================
if dataey == 0 & dataez == 0
data = dataex';
elseif dataex == 0 & dataez == 0
@@ -47,15 +48,15 @@ else
end
%Single to double precision ===============================================
% Convert to double precision ============================================
data = double(data);
%The HDF5 file doesn't contain information about the Tx Rx distance and the
%trace step. We need to provide this information ==========================
% The HDF5 file does not contain information about the Tx-Rx separation ==
% distance and the trace step. The user needs to provide this information
while 1
prompt = {'Antenna Seperation (m)', 'Trace Step (m)'};
dlg_title = 'gprMax Information';
prompt = {'Antenna Separation (m)', 'Trace Step (m)'};
dlg_title = 'gprMax information';
answer = inputdlg(prompt, dlg_title, [1 45]);
answer = str2double(answer);
if isempty(answer)
@@ -69,17 +70,18 @@ while 1
end
%Create header ===================================================== Basic information
% Create header with basic informatio =====================================
info = h5info(infile);
attributes = info.Attributes;
% Antenna name
if length(attributes) > 4
gMaxV = cell2mat(h5readatt(infile, '/', 'gprMax'));
HDR.antenna = (['gprMax ', gMaxV]); %Antenna name
HDR.antenna = (['gprMax ', gMaxV]);
else
HDR.antenna = ('gprMax'); %Antenna name
HDR.antenna = ('gprMax');
end
HDR.ant_sep = answer(1); %Tx Rx distance
HDR.ant_sep = answer(1); % Tx-Rx separation distance
[HDR.num_samp, HDR.num_trac] = size(data); % Samples & traces
HDR.trac_int = answer(2); % Trace interval
HDR.samp_int = h5readatt(infile, '/', 'dt') * 10^9; % Sampling interval
@@ -93,9 +95,9 @@ t = HDR.samp_int : HDR.samp_int : HDR.num_samp * HDR.samp_int;
%**************************************************************************
%******************************** Optional ********************************
%Resample to 1024 samples =================================================
% Resample to 1024 samples ================================================
% I usually perform this step for either 512 or 1024 samples (line 98)
%because many programs cant load files with many samples.
% because many programs cannot load files with more samples.
tx1 = 1 : HDR.num_samp;
fs1 = 1024 / HDR.num_samp;
data = resample(data, tx1, fs1);
@@ -111,8 +113,9 @@ HDR.samp_freq = (1 / HDR.samp_int) * 10^3; %New sampli
%**************************************************************************
%******************************** Optional ********************************
%Data scale ===============================================================
data = data * 32767.5 ./ max(max(abs(data))); %signal * ((1 - 1 / 2^bitrate) * 32768) / max(signal)
% Data scale ==============================================================
data = data * 32767.5 ./ max(max(abs(data)));
%signal * ((1 - 1 / 2^bitrate) * 32768) / max(signal)
data(data >= 32767.5) = 32767;
data = round(data);
@@ -120,7 +123,7 @@ data = round(data);
%**************************************************************************
%Plots ====================================================================
% Plots ===================================================================
pmin = min(data(:));
pmax = max(data(:));
@@ -134,9 +137,14 @@ clims = [pmin pmax];
colormap (flipud(bone(256)));
im1 = imagesc(x, t, data, clims);
set(im1, 'cdatamapping', 'scaled');
title(HDR.fname); xlabel('Distance (m)'); ylabel('Time (ns)');
ax1 = gca; ax1.XAxisLocation = 'Top'; ax1.FontSize = 12;
box off; movegui(f1, 'northeast');
title(HDR.fname);
xlabel('Distance (m)');
ylabel('Time (ns)');
ax1 = gca;
ax1.XAxisLocation = 'Top';
ax1.FontSize = 12;
box off;
movegui(f1, 'northeast');
% Frequency Spectrum
@@ -153,12 +161,16 @@ f2 = figure('Name', 'Frequency Spectrum', ...
'Toolbar', 'Figure');
area(freq, amp, 'FaceColor', 'black');
title(HDR.fname); xlabel('Frequency (MHz)'); ylabel('Amplitude');
ax2 = gca; ax2.FontSize = 12;
box off; movegui(f2, 'southeast');
title(HDR.fname);
xlabel('Frequency (MHz)');
ylabel('Amplitude');
ax2 = gca;
ax2.FontSize = 12;
box off;
movegui(f2, 'southeast');
%Export option rd3 or dzt or dt1 ==========================================
% Export option rd3 or dzt or dt1 =========================================
while 1
choice = questdlg('File Type', 'Export', ...
'RD3', 'DZT', 'DT1', 'Default');
@@ -178,10 +190,10 @@ end
wb = waitbar(0, 'Exporting...', 'Name', 'Exporting File');
%RAD / RD3, Mala GeoScience ===============================================
%Rad is the header file. In this file are contained all the important
% RAD / RD3, Mala GeoScience ==============================================
% Rad is the header file. In this file is all the important
% information such as the number of traces, samples, measurement intervals...
%Rd3 is the data file. This file contain only the data (amplitudes) in a
% Rd3 is the data file. This file contains only the data (amplitudes) in a
% binary form.
if flg == 1
% Header structure
@@ -273,14 +285,14 @@ if flg == 1
fclose(fid);
%DZT, Geophysical Survey Systems Inc. (GSSI) ==============================
% DZT, Geophysical Survey Systems Inc. (GSSI) =============================
% Dzt is a binary file that consists of the file header with all the
% important information (number of traces, samples, channels, etc.) followed
% by the data section.
%Every information that is needed is contained in this file except the
%TxRx distance (antenna seperation). There is a possibility that the
% All information that is needed is contained in this file except the
% TxRx distance (antenna separation). There is a possibility that the
% official GSSI software has stored this information and by using the
%antenna name presents the correct one. All the other software doesn't
% antenna name presents the correct one. All the other software does not
% detect the TxRx distance.
elseif flg == 2
@@ -294,7 +306,7 @@ elseif flg == 2
HDR.scans_per_second = 0; % Scans per second
HDR.scans_per_meter = 1 / HDR.trac_int; % Scans per meter
HDR.meters_per_mark = 0; % Meters per mark
HDR.zero_time_adjustment = 0; %Zero time adjustment (nanoseconds)
HDR.zero_time_adjustment = 0; % Time zero adjustment (nanoseconds)
HDR.time_window = HDR.time_window; % Time window (with no corrections i.e zero time)
HDR.scans_per_pass = 0; % Scan per pass for 2D files
@@ -404,8 +416,8 @@ elseif flg == 2
fclose(fid);
%HD / DT1, Sensors & Software Inc. ========================================
%Hd is the header file. In this file can be found all the important
% HD / DT1, Sensors & Software Inc. =======================================
% Hd is the header file. In this file is all the important
% information such as the number of traces, samples, stacks, etc.
% Dt1 is the data file written in binary form. This file contains as many
% records as there are traces. Each record consists of a header section and
@@ -431,7 +443,7 @@ elseif flg == 3
HDR.trac_int = HDR.trac_int; % Trace interval
HDR.pos_units = 'm'; % Position units
HDR.nominal_freq = 'Unknown'; % Nominal frequency
HDR.ant_sep = HDR.ant_sep; %Antenna seperation
HDR.ant_sep = HDR.ant_sep; % Antenna separation
HDR.pulser_voltage = 'Unknown'; % Pulser voltage (V)
HDR.stacks = 0; % Number of stacks
HDR.survey_mode = 'Reflection'; % Survey mode
@@ -530,4 +542,7 @@ elseif flg == 3
end
fclose(fid);
end
waitbar(1, wb, 'Done!!!'); pause(1); close(wb);
waitbar(1, wb, 'Done!!!');
pause(1);
close(wb);