你已经派生过 gprMax
镜像自地址
https://gitee.com/sunhf/gprMax.git
已同步 2025-08-08 07:24:19 +08:00
Update to use pathlib.
这个提交包含在：
Craig Warren
@@ -17,7 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import argparse
|
||||
import os
|
||||
from pathlib import Path
|
||||
import sys
|
||||
|
||||
import h5py
|
||||
@@ -29,14 +29,17 @@ from gprMax.exceptions import CmdInputError
|
||||
"""Plots a comparison of fields between given simulation output and experimental data files."""
|
||||
|
||||
# Parse command line arguments
|
||||
parser = argparse.ArgumentParser(description='Plots a comparison of fields between given simulation output and experimental data files.', usage='cd gprMax; python -m tests.test_compare_experimental modelfile realfile output')
|
||||
parser = argparse.ArgumentParser(description='Plots a comparison of fields between given simulation output and experimental data files.', usage='cd gprMax; python -m tests.test_experimental modelfile realfile output')
|
||||
parser.add_argument('modelfile', help='name of model output file including path')
|
||||
parser.add_argument('realfile', help='name of file containing experimental data including path')
|
||||
parser.add_argument('output', help='output to be plotted, i.e. Ex Ey Ez', nargs='+')
|
||||
args = parser.parse_args()
|
||||
|
||||
modelfile = Path(args.modelfile)
|
||||
realfile = Path(args.realfile)
|
||||
|
||||
# Model results
|
||||
f = h5py.File(args.modelfile, 'r')
|
||||
f = h5py.File(Path(modelfile), 'r')
|
||||
path = '/rxs/rx1/'
|
||||
availablecomponents = list(f[path].keys())
|
||||
|
||||
@@ -48,7 +51,7 @@ else:
|
||||
polarity = 1
|
||||
|
||||
if args.output[0] not in availablecomponents:
|
||||
raise CmdInputError('{} output requested to plot, but the available output for receiver 1 is {}'.format(args.output[0], ', '.join(availablecomponents)))
|
||||
raise CmdInputError(f"{args.output[0]} output requested to plot, but the available output for receiver 1 is {', '.join(availablecomponents)}")
|
||||
|
||||
floattype = f[path + args.output[0]].dtype
|
||||
iterations = f.attrs['Iterations']
|
||||
@@ -64,7 +67,7 @@ f.close()
|
||||
modelmax = np.where(np.abs(model) == 1)[0][0]
|
||||
|
||||
# Real results
|
||||
with open(args.realfile, 'r') as f:
|
||||
with open(realfile, 'r') as f:
|
||||
real = np.loadtxt(f)
|
||||
real[:, 1] = real[:, 1] / np.amax(np.abs(real[:, 1]))
|
||||
realmax = np.where(np.abs(real[:, 1]) == 1)[0][0]
|
||||
@@ -72,19 +75,20 @@ realmax = np.where(np.abs(real[:, 1]) == 1)[0][0]
|
||||
difftime = - (timemodel[modelmax] - real[realmax, 0])
|
||||
|
||||
# Plot modelled and real data
|
||||
fig, ax = plt.subplots(num=args.modelfile + ' versus ' + args.realfile, figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, ax = plt.subplots(num=modelfile.stem + '_vs_' + realfile.stem, figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
ax.plot(timemodel + difftime, model, 'r', lw=2, label='Model')
|
||||
ax.plot(real[:, 0], real[:, 1], 'r', ls='--', lw=2, label='Experiment')
|
||||
ax.set_xlabel('Time [s]')
|
||||
ax.set_ylabel('Amplitude')
|
||||
ax.set_xlim([0, timemodel[-1]])
|
||||
ax.set_ylim([-1, 1])
|
||||
# ax.set_ylim([-1, 1])
|
||||
ax.legend()
|
||||
ax.grid()
|
||||
|
||||
# Save a PDF/PNG of the figure
|
||||
savename = os.path.abspath(os.path.dirname(args.modelfile)) + os.sep + os.path.splitext(os.path.split(args.modelfile)[1])[0] + '_vs_' + os.path.splitext(os.path.split(args.realfile)[1])[0]
|
||||
# fig.savefig(savename + '.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig((savename + '.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
savename = modelfile.stem + '_vs_' + realfile.stem
|
||||
savename = modelfile.parent / savename
|
||||
# fig.savefig(savename.with_suffix('.pdf'), dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(savename.with_suffix('.png'), dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
plt.show()
|
||||
|
||||
@@ -44,12 +44,11 @@ modelset = 'models_basic'
|
||||
# modelset += 'models_advanced'
|
||||
# modelset += 'models_pmls'
|
||||
|
||||
basepath = Path(__file__).parents[0] / 'tests' / modelset
|
||||
basepath = Path(__file__).parents[0] / modelset
|
||||
|
||||
|
||||
# List of available basic test models
|
||||
# testmodels = ['hertzian_dipole_fs_analytical', '2D_ExHyHz', '2D_EyHxHz', '2D_EzHxHy', 'cylinder_Ascan_2D', 'hertzian_dipole_fs', 'hertzian_dipole_hs', 'hertzian_dipole_dispersive', 'magnetic_dipole_fs']
|
||||
testmodels = ['2D_ExHyHz', '2D_EyHxHz', '2D_EzHxHy']
|
||||
testmodels = ['hertzian_dipole_fs_analytical', '2D_ExHyHz', '2D_EyHxHz', '2D_EzHxHy', 'cylinder_Ascan_2D', 'hertzian_dipole_fs', 'hertzian_dipole_hs', 'hertzian_dipole_dispersive', 'magnetic_dipole_fs']
|
||||
|
||||
# List of available advanced test models
|
||||
# testmodels = ['antenna_GSSI_1500_fs', 'antenna_MALA_1200_fs']
|
||||
@@ -72,12 +71,12 @@ for i, model in enumerate(testmodels):
|
||||
|
||||
# Run model
|
||||
file = basepath / model / model
|
||||
gprMax.run(inputfile=file.with_suffix('.in'), gpu=None)
|
||||
gprMax.run(inputfile=file.with_suffix('.in'), gpu=[0])
|
||||
|
||||
# Special case for analytical comparison
|
||||
if model == 'hertzian_dipole_fs_analytical':
|
||||
# Get output for model file
|
||||
filetest = h5py.File(file.with_suffix('.out'), 'r')
|
||||
filetest = h5py.File(file.with_suffix('.h5'), 'r')
|
||||
testresults[model]['Test version'] = filetest.attrs['gprMax']
|
||||
|
||||
# Get available field output component names
|
||||
@@ -109,8 +108,8 @@ for i, model in enumerate(testmodels):
|
||||
# Get output for model and reference files
|
||||
fileref = file.stem + '_ref'
|
||||
fileref = file.parent / Path(fileref)
|
||||
fileref = h5py.File(fileref.with_suffix('.out'), 'r')
|
||||
filetest = h5py.File(file.with_suffix('.out'), 'r')
|
||||
fileref = h5py.File(fileref.with_suffix('.h5'), 'r')
|
||||
filetest = h5py.File(file.with_suffix('.h5'), 'r')
|
||||
testresults[model]['Ref version'] = fileref.attrs['gprMax']
|
||||
testresults[model]['Test version'] = filetest.attrs['gprMax']
|
||||
|
||||
@@ -201,8 +200,8 @@ for i, model in enumerate(testmodels):
|
||||
filediffs = file.parent / Path(filediffs)
|
||||
# fig1.savefig(file.with_suffix('.pdf'), dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(savediffs.with_suffix('.pdf'), dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig1.savefig(file.with_suffix('.png'), dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(filediffs.with_suffix('.png'), dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
fig1.savefig(file.with_suffix('.png'), dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
fig2.savefig(filediffs.with_suffix('.png'), dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
# Summary of results
|
||||
for name, data in sorted(testresults.items()):
|
||||
|
||||
@@ -19,6 +19,7 @@
|
||||
import argparse
|
||||
import glob
|
||||
import os
|
||||
from pathlib import Path
|
||||
|
||||
import h5py
|
||||
import numpy as np
|
||||
@@ -46,14 +47,14 @@ def get_output_data(filename, rxnumber, rxcomponent):
|
||||
|
||||
# Check there are any receivers
|
||||
if nrx == 0:
|
||||
raise CmdInputError('No receivers found in {}'.format(filename))
|
||||
raise CmdInputError(f'No receivers found in {filename}')
|
||||
|
||||
path = '/rxs/rx' + str(rxnumber) + '/'
|
||||
availableoutputs = list(f[path].keys())
|
||||
|
||||
# Check if requested output is in file
|
||||
if rxcomponent not in availableoutputs:
|
||||
raise CmdInputError('{} output requested to plot, but the available output for receiver 1 is {}'.format(rxcomponent, ', '.join(availableoutputs)))
|
||||
raise CmdInputError(f"{rxcomponent} output requested to plot, but the available output for receiver 1 is {', '.join(availableoutputs)}")
|
||||
|
||||
outputdata = f[path + '/' + rxcomponent]
|
||||
outputdata = np.array(outputdata)
|
||||
@@ -71,17 +72,16 @@ def merge_files(basefilename, removefiles=False):
|
||||
outputs (boolean): Flag to remove individual output files after merge.
|
||||
"""
|
||||
|
||||
outputfile = basefilename + '_merged.out'
|
||||
files = glob.glob(basefilename + '*.out')
|
||||
outputfile = basefilename + '_merged.h5'
|
||||
files = glob.glob(basefilename + '*.h5')
|
||||
outputfiles = [filename for filename in files if '_merged' not in filename]
|
||||
modelruns = len(outputfiles)
|
||||
|
||||
# Combined output file
|
||||
fout = h5py.File(outputfile, 'w')
|
||||
|
||||
# Add positional data for rxs
|
||||
for model in range(modelruns):
|
||||
fin = h5py.File(basefilename + str(model + 1) + '.out', 'r')
|
||||
fin = h5py.File(basefilename + str(model + 1) + '.h5', 'r')
|
||||
nrx = fin.attrs['nrx']
|
||||
|
||||
# Write properties for merged file on first iteration
|
||||
@@ -117,7 +117,7 @@ def merge_files(basefilename, removefiles=False):
|
||||
|
||||
if removefiles:
|
||||
for model in range(modelruns):
|
||||
file = basefilename + str(model + 1) + '.out'
|
||||
file = basefilename + str(model + 1) + '.h5'
|
||||
os.remove(file)
|
||||
|
||||
if __name__ == "__main__":
|
||||
|
||||
@@ -17,7 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import argparse
|
||||
import os
|
||||
from pathlib import Path
|
||||
|
||||
import h5py
|
||||
import numpy as np
|
||||
@@ -30,7 +30,9 @@ from gprMax.utilities import fft_power
|
||||
|
||||
|
||||
def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
"""Plots electric and magnetic fields and currents from all receiver points in the given output file. Each receiver point is plotted in a new figure window.
|
||||
"""Plots electric and magnetic fields and currents from all receiver points
|
||||
in the given output file. Each receiver point is plotted in a new figure
|
||||
window.
|
||||
|
||||
Args:
|
||||
filename (string): Filename (including path) of output file.
|
||||
@@ -41,8 +43,10 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
plt (object): matplotlib plot object.
|
||||
"""
|
||||
|
||||
file = Path(filename)
|
||||
|
||||
# Open output file and read some attributes
|
||||
f = h5py.File(filename, 'r')
|
||||
f = h5py.File(file, 'r')
|
||||
nrx = f.attrs['nrx']
|
||||
dt = f.attrs['dt']
|
||||
iterations = f.attrs['Iterations']
|
||||
@@ -50,7 +54,7 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
|
||||
# Check there are any receivers
|
||||
if nrx == 0:
|
||||
raise CmdInputError(f'No receivers found in {filename}')
|
||||
raise CmdInputError(f'No receivers found in {file}')
|
||||
|
||||
# Check for single output component when doing a FFT
|
||||
if fft:
|
||||
@@ -96,7 +100,9 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
pltrange = np.s_[0:pltrange]
|
||||
|
||||
# Plot time history of output component
|
||||
fig, (ax1, ax2) = plt.subplots(nrows=1, ncols=2, num='rx' + str(rx), figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, (ax1, ax2) = plt.subplots(nrows=1, ncols=2, num='rx' + str(rx),
|
||||
figsize=(20, 10), facecolor='w',
|
||||
edgecolor='w')
|
||||
line1 = ax1.plot(time, outputdata, 'r', lw=2, label=outputtext)
|
||||
ax1.set_xlabel('Time [s]')
|
||||
ax1.set_ylabel(outputtext + ' field strength [V/m]')
|
||||
@@ -104,7 +110,9 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
ax1.grid(which='both', axis='both', linestyle='-.')
|
||||
|
||||
# Plot frequency spectra
|
||||
markerline, stemlines, baseline = ax2.stem(freqs[pltrange], power[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax2.stem(freqs[pltrange],
|
||||
power[pltrange], '-.',
|
||||
use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'r')
|
||||
plt.setp(markerline, 'markerfacecolor', 'r', 'markeredgecolor', 'r')
|
||||
@@ -131,7 +139,10 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
|
||||
# Plotting if no FFT required
|
||||
else:
|
||||
fig, ax = plt.subplots(subplot_kw=dict(xlabel='Time [s]', ylabel=outputtext + ' field strength [V/m]'), num='rx' + str(rx), figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, ax = plt.subplots(subplot_kw=dict(xlabel='Time [s]',
|
||||
ylabel=outputtext + ' field strength [V/m]'),
|
||||
num='rx' + str(rx), figsize=(20, 10),
|
||||
facecolor='w', edgecolor='w')
|
||||
line = ax.plot(time, outputdata, 'r', lw=2, label=outputtext)
|
||||
ax.set_xlim([0, np.amax(time)])
|
||||
# ax.set_ylim([-15, 20])
|
||||
@@ -146,7 +157,9 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
|
||||
# If multiple outputs required, create all nine subplots and populate only the specified ones
|
||||
else:
|
||||
fig, ax = plt.subplots(subplot_kw=dict(xlabel='Time [s]'), num='rx' + str(rx), figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, ax = plt.subplots(subplot_kw=dict(xlabel='Time [s]'),
|
||||
num='rx' + str(rx), figsize=(20, 10),
|
||||
facecolor='w', edgecolor='w')
|
||||
if len(outputs) == 9:
|
||||
gs = gridspec.GridSpec(3, 3, hspace=0.3, wspace=0.3)
|
||||
else:
|
||||
@@ -215,8 +228,12 @@ def mpl_plot(filename, outputs=Rx.defaultoutputs, fft=False):
|
||||
ax.grid(which='both', axis='both', linestyle='-.')
|
||||
|
||||
# Save a PDF/PNG of the figure
|
||||
# fig.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_rx' + str(rx) + '.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_rx' + str(rx) + '.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
savename = file.stem + '_rx' + str(rx)
|
||||
savename = file.parent / savename
|
||||
# fig.savefig(savename.with_suffix('.pdf'), dpi=None, format='pdf',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(savename.with_suffix('.png'), dpi=150, format='png',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
f.close()
|
||||
|
||||
@@ -228,8 +245,12 @@ if __name__ == "__main__":
|
||||
# Parse command line arguments
|
||||
parser = argparse.ArgumentParser(description='Plots electric and magnetic fields and currents from all receiver points in the given output file. Each receiver point is plotted in a new figure window.', usage='cd gprMax; python -m tools.plot_Ascan outputfile')
|
||||
parser.add_argument('outputfile', help='name of output file including path')
|
||||
parser.add_argument('--outputs', help='outputs to be plotted', default=Rx.defaultoutputs, choices=['Ex', 'Ey', 'Ez', 'Hx', 'Hy', 'Hz', 'Ix', 'Iy', 'Iz', 'Ex-', 'Ey-', 'Ez-', 'Hx-', 'Hy-', 'Hz-', 'Ix-', 'Iy-', 'Iz-'], nargs='+')
|
||||
parser.add_argument('-fft', action='store_true', help='plot FFT (single output must be specified)', default=False)
|
||||
parser.add_argument('--outputs', help='outputs to be plotted',
|
||||
default=Rx.defaultoutputs,
|
||||
choices=['Ex', 'Ey', 'Ez', 'Hx', 'Hy', 'Hz', 'Ix', 'Iy', 'Iz', 'Ex-', 'Ey-', 'Ez-', 'Hx-', 'Hy-', 'Hz-', 'Ix-', 'Iy-', 'Iz-'],
|
||||
nargs='+')
|
||||
parser.add_argument('-fft', action='store_true', help='plot FFT (single output must be specified)',
|
||||
default=False)
|
||||
args = parser.parse_args()
|
||||
|
||||
plthandle = mpl_plot(args.outputfile, args.outputs, fft=args.fft)
|
||||
|
||||
@@ -17,8 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import argparse
|
||||
import os
|
||||
import sys
|
||||
from pathlib import Path
|
||||
|
||||
import h5py
|
||||
import numpy as np
|
||||
@@ -42,13 +41,15 @@ def mpl_plot(filename, outputdata, dt, rxnumber, rxcomponent):
|
||||
plt (object): matplotlib plot object.
|
||||
"""
|
||||
|
||||
(path, filename) = os.path.split(filename)
|
||||
file = Path(filename)
|
||||
|
||||
fig = plt.figure(num=filename + ' - rx' + str(rxnumber), figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
plt.imshow(outputdata, extent=[0, outputdata.shape[1], outputdata.shape[0] * dt, 0], interpolation='nearest', aspect='auto', cmap='seismic', vmin=-np.amax(np.abs(outputdata)), vmax=np.amax(np.abs(outputdata)))
|
||||
fig = plt.figure(num=file.stem + ' - rx' + str(rxnumber), figsize=(20, 10),
|
||||
facecolor='w', edgecolor='w')
|
||||
plt.imshow(outputdata, extent=[0, outputdata.shape[1], outputdata.shape[0] * dt, 0],
|
||||
interpolation='nearest', aspect='auto', cmap='seismic',
|
||||
vmin=-np.amax(np.abs(outputdata)), vmax=np.amax(np.abs(outputdata)))
|
||||
plt.xlabel('Trace number')
|
||||
plt.ylabel('Time [s]')
|
||||
# plt.title('{}'.format(filename))
|
||||
|
||||
# Grid properties
|
||||
ax = fig.gca()
|
||||
@@ -63,9 +64,10 @@ def mpl_plot(filename, outputdata, dt, rxnumber, rxcomponent):
|
||||
cb.set_label('Current [A]')
|
||||
|
||||
# Save a PDF/PNG of the figure
|
||||
# savefile = os.path.splitext(filename)[0]
|
||||
# fig.savefig(path + os.sep + savefile + '.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(path + os.sep + savefile + '.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(file.with_suffix('.pdf'), dpi=None, format='pdf',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(file.with_suffix('.png'), dpi=150, format='png',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
return plt
|
||||
|
||||
@@ -73,9 +75,11 @@ def mpl_plot(filename, outputdata, dt, rxnumber, rxcomponent):
|
||||
if __name__ == "__main__":
|
||||
|
||||
# Parse command line arguments
|
||||
parser = argparse.ArgumentParser(description='Plots a B-scan image.', usage='cd gprMax; python -m tools.plot_Bscan outputfile output')
|
||||
parser = argparse.ArgumentParser(description='Plots a B-scan image.',
|
||||
usage='cd gprMax; python -m tools.plot_Bscan outputfile output')
|
||||
parser.add_argument('outputfile', help='name of output file including path')
|
||||
parser.add_argument('rx_component', help='name of output component to be plotted', choices=['Ex', 'Ey', 'Ez', 'Hx', 'Hy', 'Hz', 'Ix', 'Iy', 'Iz'])
|
||||
parser.add_argument('rx_component', help='name of output component to be plotted',
|
||||
choices=['Ex', 'Ey', 'Ez', 'Hx', 'Hy', 'Hz', 'Ix', 'Iy', 'Iz'])
|
||||
args = parser.parse_args()
|
||||
|
||||
# Open output file and read number of outputs (receivers)
|
||||
@@ -85,7 +89,7 @@ if __name__ == "__main__":
|
||||
|
||||
# Check there are any receivers
|
||||
if nrx == 0:
|
||||
raise CmdInputError('No receivers found in {}'.format(args.outputfile))
|
||||
raise CmdInputError(f'No receivers found in {args.outputfile}')
|
||||
|
||||
for rx in range(1, nrx + 1):
|
||||
outputdata, dt = get_output_data(args.outputfile, rx, args.rx_component)
|
||||
|
||||
@@ -17,8 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import argparse
|
||||
import os
|
||||
import sys
|
||||
from pathlib import Path
|
||||
|
||||
import h5py
|
||||
import numpy as np
|
||||
@@ -29,7 +28,8 @@ from gprMax.exceptions import CmdInputError
|
||||
|
||||
|
||||
def calculate_antenna_params(filename, tltxnumber=1, tlrxnumber=None, rxnumber=None, rxcomponent=None):
|
||||
"""Calculates antenna parameters - incident, reflected and total volatges and currents; s11, (s21) and input impedance.
|
||||
"""Calculates antenna parameters - incident, reflected and total volatges
|
||||
and currents; s11, (s21) and input impedance.
|
||||
|
||||
Args:
|
||||
filename (string): Filename (including path) of output file.
|
||||
@@ -43,7 +43,8 @@ def calculate_antenna_params(filename, tltxnumber=1, tlrxnumber=None, rxnumber=N
|
||||
"""
|
||||
|
||||
# Open output file and read some attributes
|
||||
f = h5py.File(filename, 'r')
|
||||
file = Path(filename)
|
||||
f = h5py.File(file, 'r')
|
||||
dxdydz = f.attrs['dx_dy_dz']
|
||||
dt = f.attrs['dt']
|
||||
iterations = f.attrs['Iterations']
|
||||
@@ -52,9 +53,9 @@ def calculate_antenna_params(filename, tltxnumber=1, tlrxnumber=None, rxnumber=N
|
||||
time = np.linspace(0, (iterations - 1) * dt, num=iterations)
|
||||
df = 1 / np.amax(time)
|
||||
|
||||
print('Time window: {:g} s ({} iterations)'.format(np.amax(time), iterations))
|
||||
print('Time step: {:g} s'.format(dt))
|
||||
print('Frequency bin spacing: {:g} Hz'.format(df))
|
||||
print(f'Time window: {np.amax(time):g} s ({iterations} iterations)')
|
||||
print(f'Time step: {dt:g} s')
|
||||
print(f'Frequency bin spacing: {df:g} Hz')
|
||||
|
||||
# Read/calculate voltages and currents from transmitter antenna
|
||||
tltxpath = '/tls/tl' + str(tltxnumber) + '/'
|
||||
@@ -81,7 +82,7 @@ def calculate_antenna_params(filename, tltxnumber=1, tlrxnumber=None, rxnumber=N
|
||||
availableoutputs = list(f[rxpath].keys())
|
||||
|
||||
if rxcomponent not in availableoutputs:
|
||||
raise CmdInputError('{} output requested, but the available output for receiver {} is {}'.format(rxcomponent, rxnumber, ', '.join(availableoutputs)))
|
||||
raise CmdInputError(f'{rxcomponent} output requested, but the available output for receiver {rxnumber} is {', '.join(availableoutputs)}')
|
||||
|
||||
rxpath += rxcomponent
|
||||
|
||||
@@ -145,7 +146,8 @@ def calculate_antenna_params(filename, tltxnumber=1, tlrxnumber=None, rxnumber=N
|
||||
|
||||
|
||||
def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref, Irefp, Vtotal, Vtotalp, Itotal, Itotalp, s11, zin, yin, s21=None):
|
||||
"""Plots antenna parameters - incident, reflected and total volatges and currents; s11, (s21) and input impedance.
|
||||
"""Plots antenna parameters - incident, reflected and total volatges and
|
||||
currents; s11, (s21) and input impedance.
|
||||
|
||||
Args:
|
||||
filename (string): Filename (including path) of output file.
|
||||
@@ -166,20 +168,21 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
# To a certain drop from maximum power
|
||||
pltrangemax = np.where((np.amax(Vincp[1::]) - Vincp[1::]) > 60)[0][0] + 1
|
||||
# To a maximum frequency
|
||||
# pltrangemax = np.where(freqs > 6e9)[0][0]
|
||||
pltrangemax = np.where(freqs > 3e9)[0][0]
|
||||
pltrange = np.s_[pltrangemin:pltrangemax]
|
||||
|
||||
# Print some useful values from s11, and input impedance
|
||||
s11minfreq = np.where(s11[pltrange] == np.amin(s11[pltrange]))[0][0]
|
||||
print('s11 minimum: {:g} dB at {:g} Hz'.format(np.amin(s11[pltrange]), freqs[s11minfreq + pltrangemin]))
|
||||
print('At {:g} Hz...'.format(freqs[s11minfreq + pltrangemin]))
|
||||
print('Input impedance: {:.1f}{:+.1f}j Ohms'.format(np.abs(zin[s11minfreq + pltrangemin]), zin[s11minfreq + pltrangemin].imag))
|
||||
# print('Input admittance (mag): {:g} S'.format(np.abs(yin[s11minfreq + pltrangemin])))
|
||||
# print('Input admittance (phase): {:.1f} deg'.format(np.angle(yin[s11minfreq + pltrangemin], deg=True)))
|
||||
print(f's11 minimum: {np.amin(s11[pltrange]):g} dB at {freqs[s11minfreq + pltrangemin]:g} Hz')
|
||||
print(f'At {freqs[s11minfreq + pltrangemin]:g} Hz...')
|
||||
print(f'Input impedance: {np.abs(zin[s11minfreq + pltrangemin]):.1f}{zin[s11minfreq + pltrangemin].imag:+.1f}j Ohms')
|
||||
# print(f'Input admittance (mag): {np.abs(yin[s11minfreq + pltrangemin]):g} S')
|
||||
# print(f'Input admittance (phase): {np.angle(yin[s11minfreq + pltrangemin], deg=True):.1f} deg')
|
||||
|
||||
# Figure 1
|
||||
# Plot incident voltage
|
||||
fig1, ax = plt.subplots(num='Transmitter transmission line parameters', figsize=(20, 12), facecolor='w', edgecolor='w')
|
||||
fig1, ax = plt.subplots(num='Transmitter transmission line parameters',
|
||||
figsize=(20, 12), facecolor='w', edgecolor='w')
|
||||
gs1 = gridspec.GridSpec(4, 2, hspace=0.7)
|
||||
ax = plt.subplot(gs1[0, 0])
|
||||
ax.plot(time, Vinc, 'r', lw=2, label='Vinc')
|
||||
@@ -191,7 +194,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of incident voltage
|
||||
ax = plt.subplot(gs1[0, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vincp[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vincp[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'r')
|
||||
plt.setp(markerline, 'markerfacecolor', 'r', 'markeredgecolor', 'r')
|
||||
@@ -212,7 +216,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of incident current
|
||||
ax = plt.subplot(gs1[1, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Iincp[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Iincp[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'b')
|
||||
plt.setp(markerline, 'markerfacecolor', 'b', 'markeredgecolor', 'b')
|
||||
@@ -233,7 +238,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of total voltage
|
||||
ax = plt.subplot(gs1[2, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vtotalp[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vtotalp[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'r')
|
||||
plt.setp(markerline, 'markerfacecolor', 'r', 'markeredgecolor', 'r')
|
||||
@@ -254,7 +260,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of total current
|
||||
ax = plt.subplot(gs1[3, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Itotalp[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], Itotalp[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'b')
|
||||
plt.setp(markerline, 'markerfacecolor', 'b', 'markeredgecolor', 'b')
|
||||
@@ -275,7 +282,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of reflected voltage
|
||||
# ax = plt.subplot(gs1[4, 1])
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vrefp[pltrange], '-.')
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], Vrefp[pltrange],
|
||||
# '-.', use_line_collection=True)
|
||||
# plt.setp(baseline, 'linewidth', 0)
|
||||
# plt.setp(stemlines, 'color', 'r')
|
||||
# plt.setp(markerline, 'markerfacecolor', 'r', 'markeredgecolor', 'r')
|
||||
@@ -296,7 +304,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot frequency spectra of reflected current
|
||||
# ax = plt.subplot(gs1[5, 1])
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], Irefp[pltrange], '-.')
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], Irefp[pltrange],
|
||||
# '-.', use_line_collection=True)
|
||||
# plt.setp(baseline, 'linewidth', 0)
|
||||
# plt.setp(stemlines, 'color', 'b')
|
||||
# plt.setp(markerline, 'markerfacecolor', 'b', 'markeredgecolor', 'b')
|
||||
@@ -308,10 +317,12 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Figure 2
|
||||
# Plot frequency spectra of s11
|
||||
fig2, ax = plt.subplots(num='Antenna parameters', figsize=(20, 12), facecolor='w', edgecolor='w')
|
||||
fig2, ax = plt.subplots(num='Antenna parameters', figsize=(20, 12),
|
||||
facecolor='w', edgecolor='w')
|
||||
gs2 = gridspec.GridSpec(2, 2, hspace=0.3)
|
||||
ax = plt.subplot(gs2[0, 0])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], s11[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], s11[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'g')
|
||||
plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -326,7 +337,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
# Plot frequency spectra of s21
|
||||
if s21 is not None:
|
||||
ax = plt.subplot(gs2[0, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], s21[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], s21[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'g')
|
||||
plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -340,7 +352,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot input resistance (real part of impedance)
|
||||
ax = plt.subplot(gs2[1, 0])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], zin[pltrange].real, '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], zin[pltrange].real,
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'g')
|
||||
plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -355,7 +368,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot input reactance (imaginery part of impedance)
|
||||
ax = plt.subplot(gs2[1, 1])
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], zin[pltrange].imag, '-.')
|
||||
markerline, stemlines, baseline = ax.stem(freqs[pltrange], zin[pltrange].imag,
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'g')
|
||||
plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -369,7 +383,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot input admittance (magnitude)
|
||||
# ax = plt.subplot(gs2[2, 0])
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], np.abs(yin[pltrange]), '-.')
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], np.abs(yin[pltrange]),
|
||||
# '-.', use_line_collection=True)
|
||||
# plt.setp(baseline, 'linewidth', 0)
|
||||
# plt.setp(stemlines, 'color', 'g')
|
||||
# plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -383,7 +398,8 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
|
||||
# Plot input admittance (phase)
|
||||
# ax = plt.subplot(gs2[2, 1])
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], np.angle(yin[pltrange], deg=True), '-.')
|
||||
# markerline, stemlines, baseline = ax.stem(freqs[pltrange], np.angle(yin[pltrange], deg=True),
|
||||
# '-.', use_line_collection=True)
|
||||
# plt.setp(baseline, 'linewidth', 0)
|
||||
# plt.setp(stemlines, 'color', 'g')
|
||||
# plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g')
|
||||
@@ -396,10 +412,18 @@ def mpl_plot(filename, time, freqs, Vinc, Vincp, Iinc, Iincp, Vref, Vrefp, Iref,
|
||||
# ax.grid(which='both', axis='both', linestyle='-.')
|
||||
|
||||
# Save a PDF/PNG of the figure
|
||||
# fig1.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_tl_params.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_ant_params.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig1.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_tl_params.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(os.path.splitext(os.path.abspath(filename))[0] + '_ant_params.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
savename1 = file.stem + '_tl_params'
|
||||
savename1 = file.parent / savename1
|
||||
savename2 = file.stem + '_ant_params'
|
||||
savename2 = file.parent / savename2
|
||||
# fig1.savefig(savename1.with_suffix('.png'), dpi=150, format='png',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(savename2.with_suffix('.png'), dpi=150, format='png',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig1.savefig(savename1.with_suffix('.pdf'), dpi=None, format='pdf',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig2.savefig(savename2.with_suffix('.pdf'), dpi=None, format='pdf',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
return plt
|
||||
|
||||
|
||||
@@ -17,8 +17,7 @@
|
||||
# along with gprMax. If not, see <http://www.gnu.org/licenses/>.
|
||||
|
||||
import argparse
|
||||
import os
|
||||
import sys
|
||||
from pathlib import Path
|
||||
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
@@ -82,21 +81,22 @@ def mpl_plot(w, timewindow, dt, iterations, fft=False):
|
||||
timeiter.iternext()
|
||||
|
||||
print('Waveform characteristics...')
|
||||
print('Type: {}'.format(w.type))
|
||||
print('Maximum (absolute) amplitude: {:g}'.format(np.max(np.abs(waveform))))
|
||||
print(f'Type: {w.type}')
|
||||
print(f'Maximum (absolute) amplitude: {np.max(np.abs(waveform)):g}')
|
||||
|
||||
if w.freq and not w.type == 'gaussian':
|
||||
print('Centre frequency: {:g} Hz'.format(w.freq))
|
||||
print(f'Centre frequency: {w.freq:g} Hz')
|
||||
|
||||
if w.type == 'gaussian' or w.type == 'gaussiandot' or w.type == 'gaussiandotnorm' or w.type == 'gaussianprime' or w.type == 'gaussiandoubleprime':
|
||||
if (w.type == 'gaussian' or w.type == 'gaussiandot' or w.type == 'gaussiandotnorm'
|
||||
or w.type == 'gaussianprime' or w.type == 'gaussiandoubleprime'):
|
||||
delay = 1 / w.freq
|
||||
print('Time to centre of pulse: {:g} s'.format(delay))
|
||||
print(f'Time to centre of pulse: {delay:g} s')
|
||||
elif w.type == 'gaussiandotdot' or w.type == 'gaussiandotdotnorm' or w.type == 'ricker':
|
||||
delay = np.sqrt(2) / w.freq
|
||||
print('Time to centre of pulse: {:g} s'.format(delay))
|
||||
print(f'Time to centre of pulse: {delay:g} s')
|
||||
|
||||
print('Time window: {:g} s ({} iterations)'.format(timewindow, iterations))
|
||||
print('Time step: {:g} s'.format(dt))
|
||||
print(f'Time window: {timewindow:g} s ({iterations} iterations)')
|
||||
print(f'Time step: {dt:g} s')
|
||||
|
||||
if fft:
|
||||
# FFT
|
||||
@@ -111,7 +111,9 @@ def mpl_plot(w, timewindow, dt, iterations, fft=False):
|
||||
pltrange = np.where(freqs > 4 * w.freq)[0][0]
|
||||
pltrange = np.s_[0:pltrange]
|
||||
|
||||
fig, (ax1, ax2) = plt.subplots(nrows=1, ncols=2, num=w.type, figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, (ax1, ax2) = plt.subplots(nrows=1, ncols=2, num=w.type,
|
||||
figsize=(20, 10), facecolor='w',
|
||||
edgecolor='w')
|
||||
|
||||
# Plot waveform
|
||||
ax1.plot(time, waveform, 'r', lw=2)
|
||||
@@ -119,7 +121,8 @@ def mpl_plot(w, timewindow, dt, iterations, fft=False):
|
||||
ax1.set_ylabel('Amplitude')
|
||||
|
||||
# Plot frequency spectra
|
||||
markerline, stemlines, baseline = ax2.stem(freqs[pltrange], power[pltrange], '-.')
|
||||
markerline, stemlines, baseline = ax2.stem(freqs[pltrange], power[pltrange],
|
||||
'-.', use_line_collection=True)
|
||||
plt.setp(baseline, 'linewidth', 0)
|
||||
plt.setp(stemlines, 'color', 'r')
|
||||
plt.setp(markerline, 'markerfacecolor', 'r', 'markeredgecolor', 'r')
|
||||
@@ -128,18 +131,23 @@ def mpl_plot(w, timewindow, dt, iterations, fft=False):
|
||||
ax2.set_ylabel('Power [dB]')
|
||||
|
||||
else:
|
||||
fig, ax1 = plt.subplots(num=w.type, figsize=(20, 10), facecolor='w', edgecolor='w')
|
||||
fig, ax1 = plt.subplots(num=w.type, figsize=(20, 10), facecolor='w',
|
||||
edgecolor='w')
|
||||
|
||||
# Plot waveform
|
||||
ax1.plot(time, waveform, 'r', lw=2)
|
||||
ax1.set_xlabel('Time [s]')
|
||||
ax1.set_ylabel('Amplitude')
|
||||
|
||||
[ax.grid(which='both', axis='both', linestyle='-.') for ax in fig.axes] # Turn on grid
|
||||
# Turn on grid
|
||||
[ax.grid(which='both', axis='both', linestyle='-.') for ax in fig.axes]
|
||||
|
||||
# Save a PDF/PNG of the figure
|
||||
# fig.savefig(os.path.dirname(os.path.abspath(__file__)) + os.sep + w.type + '.pdf', dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(os.path.dirname(os.path.abspath(__file__)) + os.sep + w.type + '.png', dpi=150, format='png', bbox_inches='tight', pad_inches=0.1)
|
||||
savefile = Path(__file__).parent / w.type
|
||||
# fig.savefig(savefile.with_suffix('.pdf'), dpi=None, format='pdf',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
# fig.savefig(savefile.with_suffix('.png'), dpi=150, format='png',
|
||||
# bbox_inches='tight', pad_inches=0.1)
|
||||
|
||||
return plt
|
||||
|
||||
@@ -158,7 +166,7 @@ if __name__ == "__main__":
|
||||
|
||||
# Check waveform parameters
|
||||
if args.type.lower() not in Waveform.types:
|
||||
raise CmdInputError('The waveform must have one of the following types {}'.format(', '.join(Waveform.types)))
|
||||
raise CmdInputError(f"The waveform must have one of the following types {', '.join(Waveform.types)}")
|
||||
if args.freq <= 0:
|
||||
raise CmdInputError('The waveform requires an excitation frequency value of greater than zero')
|
||||
|
||||
|
||||
在新工单中引用
屏蔽一个用户