diff --git a/gprMax/__init__.py b/gprMax/__init__.py index 9010516f..1580cb9b 100644 --- a/gprMax/__init__.py +++ b/gprMax/__init__.py @@ -7,54 +7,32 @@ Electromagnetic wave propagation simulation software. """ -from ._version import __version__ -from .cmds_single_use import Discretisation -from .cmds_single_use import Domain -from .cmds_single_use import TimeWindow -from .cmds_single_use import Title -from .cmds_single_use import NumThreads -from .cmds_single_use import TimeStepStabilityFactor -from .cmds_single_use import PMLCells -from .cmds_single_use import SrcSteps -from .cmds_single_use import RxSteps -from .cmds_single_use import ExcitationFile -from .cmds_multiple import Waveform -from .cmds_multiple import VoltageSource -from .cmds_multiple import HertzianDipole -from .cmds_multiple import MagneticDipole -from .cmds_multiple import TransmissionLine -from .cmds_multiple import Rx -from .cmds_multiple import RxArray -from .cmds_multiple import Snapshot -from .cmds_multiple import Material -from .cmds_multiple import AddDebyeDispersion -from .cmds_multiple import AddLorentzDispersion -from .cmds_multiple import AddDrudeDispersion -from .cmds_multiple import SoilPeplinski -from .cmds_multiple import GeometryView -from .cmds_multiple import GeometryObjectsWrite -from .cmds_multiple import PMLCFS -from .subgrids.user_objects import SubGridHSG -from .subgrids.user_objects import ReferenceRx +import gprMax.config as config -from .cmds_geometry.edge import Edge -from .cmds_geometry.plate import Plate -from .cmds_geometry.triangle import Triangle +from ._version import __version__ +from .cmds_geometry.add_grass import AddGrass +from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness +from .cmds_geometry.add_surface_water import AddSurfaceWater from .cmds_geometry.box import Box from .cmds_geometry.cylinder import Cylinder from .cmds_geometry.cylindrical_sector import CylindricalSector -from .cmds_geometry.sphere import Sphere +from .cmds_geometry.edge import Edge from .cmds_geometry.fractal_box import FractalBox -from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness -from .cmds_geometry.add_surface_water import AddSurfaceWater -from .cmds_geometry.add_grass import AddGrass from .cmds_geometry.geometry_objects_read import GeometryObjectsRead - -from .hash_cmds_file import user_libs_fn_to_scene_obj - -from .scene import Scene +from .cmds_geometry.plate import Plate +from .cmds_geometry.sphere import Sphere +from .cmds_geometry.triangle import Triangle +from .cmds_multiple import (PMLCFS, AddDebyeDispersion, AddDrudeDispersion, + AddLorentzDispersion, GeometryObjectsWrite, + GeometryView, HertzianDipole, MagneticDipole, + Material, Rx, RxArray, Snapshot, SoilPeplinski, + TransmissionLine, VoltageSource, Waveform) +from .cmds_single_use import (Discretisation, Domain, ExcitationFile, + NumThreads, PMLCells, RxSteps, SrcSteps, + TimeStepStabilityFactor, TimeWindow, Title) from .gprMax import run as run - -import gprMax.config as config +from .hash_cmds_file import user_libs_fn_to_scene_obj +from .scene import Scene +from .subgrids.user_objects import ReferenceRx, SubGridHSG __name__ = 'gprMax' diff --git a/gprMax/cmds_geometry/add_grass.py b/gprMax/cmds_geometry/add_grass.py index 41220bde..741d0e4d 100644 --- a/gprMax/cmds_geometry/add_grass.py +++ b/gprMax/cmds_geometry/add_grass.py @@ -18,15 +18,14 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..exceptions import CmdInputError -from ..fractals import FractalSurface -from ..fractals import Grass -from ..materials import Material +from ..fractals import FractalSurface, Grass +from ..materials import DispersiveMaterial from ..utilities import round_value +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) @@ -175,15 +174,15 @@ class AddGrass(UserObjectGeometry): # Check to see if grass has been already defined as a material if not any(x.ID == 'grass' for x in grid.materials): - m = Material(len(grid.materials), 'grass') + m = DispersiveMaterial(len(grid.materials), 'grass') m.averagable = False m.type = 'builtin, debye' - m.er = Material.grasseri - m.deltaer.append(Material.grassdeltaer) - m.tau.append(Material.grasstau) + m.er = DispersiveMaterial.grasseri + m.deltaer.append(DispersiveMaterial.grassdeltaer) + m.tau.append(DispersiveMaterial.grasstau) grid.materials.append(m) - if Material.maxpoles == 0: - Material.maxpoles = 1 + if config.get_model_config().materials['maxpoles'] == 0: + config.get_model_config().materials['maxpoles'] = 1 # Check if time step for model is suitable for using grass grass = next((x for x in grid.materials if x.ID == 'grass')) diff --git a/gprMax/cmds_geometry/add_surface_roughness.py b/gprMax/cmds_geometry/add_surface_roughness.py index 83aec2e2..6b20c22c 100644 --- a/gprMax/cmds_geometry/add_surface_roughness.py +++ b/gprMax/cmds_geometry/add_surface_roughness.py @@ -18,14 +18,13 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry -from ..cython.geometry_primitives import build_box from ..exceptions import CmdInputError from ..fractals import FractalSurface from ..utilities import round_value +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/add_surface_water.py b/gprMax/cmds_geometry/add_surface_water.py index 2dc9d6f7..024135d2 100644 --- a/gprMax/cmds_geometry/add_surface_water.py +++ b/gprMax/cmds_geometry/add_surface_water.py @@ -19,10 +19,11 @@ import logging import gprMax.config as config -from .cmds_geometry import UserObjectGeometry + from ..exceptions import CmdInputError -from ..materials import Material +from ..materials import DispersiveMaterial from ..utilities import round_value +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) @@ -126,15 +127,15 @@ class AddSurfaceWater(UserObjectGeometry): # Check to see if water has been already defined as a material if not any(x.ID == 'water' for x in grid.materials): - m = Material(len(grid.materials), 'water') + m = DispersiveMaterial(len(grid.materials), 'water') m.averagable = False m.type = 'builtin, debye' - m.er = Material.watereri - m.deltaer.append(Material.waterdeltaer) - m.tau.append(Material.watertau) + m.er = DispersiveMaterial.watereri + m.deltaer.append(DispersiveMaterial.waterdeltaer) + m.tau.append(DispersiveMaterial.watertau) grid.materials.append(m) - if Material.maxpoles == 0: - Material.maxpoles = 1 + if config.get_model_config().materials['maxpoles'] == 0: + config.get_model_config().materials['maxpoles'] = 1 # Check if time step for model is suitable for using water water = next((x for x in grid.materials if x.ID == 'water')) diff --git a/gprMax/cmds_geometry/box.py b/gprMax/cmds_geometry/box.py index 06f19b9b..51b5cf28 100644 --- a/gprMax/cmds_geometry/box.py +++ b/gprMax/cmds_geometry/box.py @@ -18,13 +18,13 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..cython.geometry_primitives import build_box from ..exceptions import CmdInputError from ..materials import Material +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/build_templates.py b/gprMax/cmds_geometry/build_templates.py index c1ccad15..0c722459 100644 --- a/gprMax/cmds_geometry/build_templates.py +++ b/gprMax/cmds_geometry/build_templates.py @@ -16,9 +16,7 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from jinja2 import Environment -from jinja2 import PackageLoader -from jinja2 import select_autoescape +from jinja2 import Environment, PackageLoader, select_autoescape env = Environment( loader=PackageLoader(__name__, 'templates'), diff --git a/gprMax/cmds_geometry/cylinder.py b/gprMax/cmds_geometry/cylinder.py index f0bafebc..0ead2438 100644 --- a/gprMax/cmds_geometry/cylinder.py +++ b/gprMax/cmds_geometry/cylinder.py @@ -17,13 +17,14 @@ # along with gprMax. If not, see . import logging -import numpy as np import gprMax.config as config -from .cmds_geometry import UserObjectGeometry +import numpy as np + from ..cython.geometry_primitives import build_cylinder from ..exceptions import CmdInputError from ..materials import Material +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/cylindrical_sector.py b/gprMax/cmds_geometry/cylindrical_sector.py index 4ef0fcbd..4ffa965e 100644 --- a/gprMax/cmds_geometry/cylindrical_sector.py +++ b/gprMax/cmds_geometry/cylindrical_sector.py @@ -18,13 +18,13 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..cython.geometry_primitives import build_cylindrical_sector from ..exceptions import CmdInputError from ..materials import Material +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/edge.py b/gprMax/cmds_geometry/edge.py index 19a81c24..9359ef7b 100644 --- a/gprMax/cmds_geometry/edge.py +++ b/gprMax/cmds_geometry/edge.py @@ -19,11 +19,11 @@ import logging import gprMax.config as config -from .cmds_geometry import UserObjectGeometry -from ..cython.geometry_primitives import build_edge_x -from ..cython.geometry_primitives import build_edge_y -from ..cython.geometry_primitives import build_edge_z + +from ..cython.geometry_primitives import (build_edge_x, build_edge_y, + build_edge_z) from ..exceptions import CmdInputError +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/fractal_box.py b/gprMax/cmds_geometry/fractal_box.py index 965bff56..1c742abe 100644 --- a/gprMax/cmds_geometry/fractal_box.py +++ b/gprMax/cmds_geometry/fractal_box.py @@ -18,12 +18,12 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..exceptions import CmdInputError from ..fractals import FractalVolume +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/fractal_box_builder.py b/gprMax/cmds_geometry/fractal_box_builder.py index 9a329bd1..3f6ff9f7 100644 --- a/gprMax/cmds_geometry/fractal_box_builder.py +++ b/gprMax/cmds_geometry/fractal_box_builder.py @@ -16,13 +16,13 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry -from ..cython.geometry_primitives import build_voxels_from_array -from ..cython.geometry_primitives import build_voxels_from_array_mask +from ..cython.geometry_primitives import (build_voxels_from_array, + build_voxels_from_array_mask) from ..exceptions import CmdInputError +from .cmds_geometry import UserObjectGeometry class FractalBoxBuilder(UserObjectGeometry): diff --git a/gprMax/cmds_geometry/geometry_objects_read.py b/gprMax/cmds_geometry/geometry_objects_read.py index 51c25d84..ed8145e8 100644 --- a/gprMax/cmds_geometry/geometry_objects_read.py +++ b/gprMax/cmds_geometry/geometry_objects_read.py @@ -19,14 +19,14 @@ import logging from pathlib import Path +import gprMax.config as config import h5py -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..cython.geometry_primitives import build_voxels_from_array from ..exceptions import CmdInputError from ..hash_cmds_file import get_user_objects from ..utilities import round_value +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/plate.py b/gprMax/cmds_geometry/plate.py index 4c43a295..f12ed285 100644 --- a/gprMax/cmds_geometry/plate.py +++ b/gprMax/cmds_geometry/plate.py @@ -19,11 +19,11 @@ import logging import gprMax.config as config -from .cmds_geometry import UserObjectGeometry -from ..cython.geometry_primitives import build_face_yz -from ..cython.geometry_primitives import build_face_xz -from ..cython.geometry_primitives import build_face_xy + +from ..cython.geometry_primitives import (build_face_xy, build_face_xz, + build_face_yz) from ..exceptions import CmdInputError +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/sphere.py b/gprMax/cmds_geometry/sphere.py index 49a0f7c2..ab359fa5 100644 --- a/gprMax/cmds_geometry/sphere.py +++ b/gprMax/cmds_geometry/sphere.py @@ -18,13 +18,13 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..cython.geometry_primitives import build_sphere from ..exceptions import CmdInputError from ..materials import Material +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_geometry/triangle.py b/gprMax/cmds_geometry/triangle.py index ed1cbbb2..68046594 100644 --- a/gprMax/cmds_geometry/triangle.py +++ b/gprMax/cmds_geometry/triangle.py @@ -18,13 +18,13 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config -from .cmds_geometry import UserObjectGeometry from ..cython.geometry_primitives import build_triangle from ..exceptions import CmdInputError from ..materials import Material +from .cmds_geometry import UserObjectGeometry logger = logging.getLogger(__name__) diff --git a/gprMax/cmds_multiple.py b/gprMax/cmds_multiple.py index f0ededa4..86e5c767 100644 --- a/gprMax/cmds_multiple.py +++ b/gprMax/cmds_multiple.py @@ -18,23 +18,22 @@ import logging +import gprMax.config as config import numpy as np -import gprMax.config as config from .cmds_geometry.cmds_geometry import UserObjectGeometry -from .exceptions import CmdInputError +from .exceptions import CmdInputError, GeneralError from .geometry_outputs import GeometryObjects as GeometryObjectsUser -from .materials import Material as MaterialUser from .materials import DispersiveMaterial as DispersiveMaterialUser +from .materials import Material as MaterialUser from .materials import PeplinskiSoil as PeplinskiSoilUser -from .pml import CFSParameter -from .pml import CFS +from .pml import CFS, CFSParameter from .receivers import Rx as RxUser from .snapshots import Snapshot as SnapshotUser -from .sources import VoltageSource as VoltageSourceUser from .sources import HertzianDipole as HertzianDipoleUser from .sources import MagneticDipole as MagneticDipoleUser from .sources import TransmissionLine as TransmissionLineUser +from .sources import VoltageSource as VoltageSourceUser from .subgrids.base import SubGridBase from .utilities import round_value from .waveforms import Waveform as WaveformUser @@ -162,7 +161,7 @@ class VoltageSource(UserObjectMulti): # Check if there is a waveformID in the waveforms list if not any(x.ID == waveform_id for x in grid.waveforms): - raise CmdInputError(f"'{self.params_str()}' there is no waveform with the identifier {tmp[5]}") + raise CmdInputError(f"'{self.params_str()}' there is no waveform with the identifier {waveform_id}") v = VoltageSourceUser() v.polarisation = polarisation @@ -594,7 +593,7 @@ class RxArray(UserObjectMulti): r.zcoordorigin = z r.ID = r.__class__.__name__ + '(' + str(x) + ',' + str(y) + ',' + str(z) + ')' for key in RxUser.defaultoutputs: - r.outputs[key] = np.zeros(grid.iterations, dtype=config.dtypes['float_or_double']) + r.outputs[key] = np.zeros(grid.iterations, dtype=config.sim_config.dtypes['float_or_double']) logger.info(f" Receiver at {r.xcoord * grid.dx:g}m, {r.ycoord * grid.dy:g}m, {r.zcoord * grid.dz:g}m with output component(s) {', '.join(r.outputs)} created.") grid.rxs.append(r) @@ -1080,7 +1079,7 @@ class GeometryObjectsWrite(UserObjectMulti): try: p1 = self.kwargs['p1'] p2 = self.kwargs['p2'] - filename = self.kwargs['filename'] + basefilename = self.kwargs['filename'] except KeyError: raise CmdInputError(f"'{self.params_str()}' requires exactly seven parameters") @@ -1088,9 +1087,9 @@ class GeometryObjectsWrite(UserObjectMulti): x0, y0, z0 = p1 x1, y1, z1 = p2 - g = GeometryObjectsUser(x0, y0, z0, x1, y1, z1, filename) + g = GeometryObjectsUser(x0, y0, z0, x1, y1, z1, basefilename) - logger.info(f'Geometry objects in the volume from {p1[0] * grid.dx:g}m, {p1[1] * grid.dy:g}m, {p1[2] * grid.dz:g}m, to {p2[0] * grid.dx:g}m, {p2[1] * grid.dy:g}m, {p2[2] * grid.dz:g}m, will be written to {g.filename}, with materials written to {g.materialsfilename}') + logger.info(f'Geometry objects in the volume from {p1[0] * grid.dx:g}m, {p1[1] * grid.dy:g}m, {p1[2] * grid.dz:g}m, to {p2[0] * grid.dx:g}m, {p2[1] * grid.dy:g}m, {p2[2] * grid.dz:g}m, will be written to {g.filename_hdf5}, with materials written to {g.filename_materials}') # Append the new GeometryView object to the geometry objects to write list grid.geometryobjectswrite.append(g) diff --git a/gprMax/cmds_single_use.py b/gprMax/cmds_single_use.py index 0c26f773..10142ad2 100644 --- a/gprMax/cmds_single_use.py +++ b/gprMax/cmds_single_use.py @@ -20,19 +20,16 @@ import inspect import logging from pathlib import Path -from colorama import init -from colorama import Fore -from colorama import Style -init() -import numpy as np -from scipy.constants import c -from scipy import interpolate - import gprMax.config as config +import numpy as np +from colorama import Fore, Style, init +init() +from scipy import interpolate +from scipy.constants import c + from .exceptions import CmdInputError +from .utilities import round_value, set_omp_threads from .waveforms import Waveform -from .utilities import round_value -from .utilities import set_omp_threads logger = logging.getLogger(__name__) @@ -95,31 +92,31 @@ class Domain(UserObjectSingle): try: G.nx, G.ny, G.nz = uip.discretise_point(self.kwargs['p1']) except KeyError: - raise CmdInputError(f"'{self.params_str()}' please specify a point") + raise CmdInputError(self.__str__() + ' please specify a point') if G.nx == 0 or G.ny == 0 or G.nz == 0: - raise CmdInputError(f"'{self.params_str()}' requires at least one cell in every dimension") + raise CmdInputError(self.__str__() + ' requires at least one cell in every dimension') logger.info(f"Domain size: {self.kwargs['p1'][0]:g} x {self.kwargs['p1'][1]:g} x {self.kwargs['p1'][2]:g}m ({G.nx:d} x {G.ny:d} x {G.nz:d} = {(G.nx * G.ny * G.nz):g} cells)") # Calculate time step at CFL limit; switch off appropriate PMLs for 2D if G.nx == 1: - G.mode = '2D TMx' + config.get_model_config().mode = '2D TMx' G.pmlthickness['x0'] = 0 G.pmlthickness['xmax'] = 0 elif G.ny == 1: - G.mode = '2D TMy' + config.get_model_config().mode = '2D TMy' G.pmlthickness['y0'] = 0 G.pmlthickness['ymax'] = 0 elif G.nz == 1: - G.mode = '2D TMz' + config.get_model_config().mode = '2D TMz' G.pmlthickness['z0'] = 0 G.pmlthickness['zmax'] = 0 else: - G.mode = '3D' + config.get_model_config().mode = '3D' G.calculate_dt() - logger.info(f'Mode: {G.mode}') + logger.info(f'Mode: {config.get_model_config().mode}') logger.info(f'Time step (at CFL limit): {G.dt:g} secs') @@ -140,14 +137,14 @@ class Discretisation(UserObjectSingle): G.dl = np.array(self.kwargs['p1']) G.dx, G.dy, G.dz = self.kwargs['p1'] except KeyError: - raise CmdInputError(f"'{self.params_str()}' discretisation requires a point") + raise CmdInputError(self.__str__() + ' discretisation requires a point') if G.dl[0] <= 0: - raise CmdInputError(f"'{self.params_str()}' discretisation requires the x-direction spatial step to be greater than zero") + raise CmdInputError(self.__str__() + ' discretisation requires the x - direction spatial step to be greater than zero') if G.dl[1] <= 0: - raise CmdInputError(f"'{self.params_str()}' discretisation requires the y-direction spatial step to be greater than zero") + raise CmdInputError(self.__str__() + ' discretisation requires the y - direction spatial step to be greater than zero') if G.dl[2] <= 0: - raise CmdInputError(f"'{self.params_str()}' discretisation requires the z-direction spatial step to be greater than zero") + raise CmdInputError(self.__str__() + ' discretisation requires the z - direction spatial step to be greater than zero') logger.info(f'Spatial discretisation: {G.dl[0]:g} x {G.dl[1]:g} x {G.dl[2]:g}m') @@ -182,12 +179,12 @@ class TimeWindow(UserObjectSingle): G.timewindow = tmp G.iterations = int(np.ceil(tmp / G.dt)) + 1 else: - raise CmdInputError(f"'{self.params_str()}' must have a value greater than zero") + raise CmdInputError(self.__str__() + ' must have a value greater than zero') except KeyError: pass if not G.timewindow: - raise CmdInputError(f"'{self.params_str()}' specify a time or number of iterations") + raise CmdInputError(self.__str__() + ' specify a time or number of iterations') logger.info(f'Time window: {G.timewindow:g} secs ({G.iterations} iterations)') diff --git a/gprMax/config.py b/gprMax/config.py index 31bc69f1..a878bacb 100644 --- a/gprMax/config.py +++ b/gprMax/config.py @@ -17,23 +17,19 @@ # along with gprMax. If not, see . import logging -from pathlib import Path import sys +from pathlib import Path -from colorama import init -from colorama import Fore -from colorama import Style -init() import cython import numpy as np +from colorama import Fore, Style, init +init() from scipy.constants import c from scipy.constants import epsilon_0 as e0 from scipy.constants import mu_0 as m0 from .exceptions import GeneralError -from .utilities import detect_gpus -from .utilities import get_host_info -from .utilities import get_terminal_width +from .utilities import detect_gpus, get_host_info, get_terminal_width logger = logging.getLogger(__name__) diff --git a/gprMax/contexts.py b/gprMax/contexts.py index 0661fa2d..d51bf13f 100644 --- a/gprMax/contexts.py +++ b/gprMax/contexts.py @@ -21,15 +21,12 @@ import logging import sys import gprMax.config as config + from ._version import __version__, codename from .exceptions import GeneralError from .model_build_run import ModelBuildRun -from .solvers import create_solver -from .solvers import create_G -from .utilities import get_terminal_width -from .utilities import human_size -from .utilities import logo -from .utilities import timer +from .solvers import create_G, create_solver +from .utilities import get_terminal_width, human_size, logo, timer logger = logging.getLogger(__name__) diff --git a/gprMax/cython/fields_updates_hsg.pyx b/gprMax/cython/fields_updates_hsg.pyx index f212d628..4550b8a8 100644 --- a/gprMax/cython/fields_updates_hsg.pyx +++ b/gprMax/cython/fields_updates_hsg.pyx @@ -18,7 +18,6 @@ import numpy as np cimport numpy as np - from cython.parallel import prange diff --git a/gprMax/exceptions.py b/gprMax/exceptions.py index 2172a65c..28817a6c 100644 --- a/gprMax/exceptions.py +++ b/gprMax/exceptions.py @@ -18,8 +18,7 @@ import logging -from colorama import init -from colorama import Fore +from colorama import Fore, init init() logger = logging.getLogger(__name__) diff --git a/gprMax/fields_outputs.py b/gprMax/fields_outputs.py index 87b72289..3ea924ba 100644 --- a/gprMax/fields_outputs.py +++ b/gprMax/fields_outputs.py @@ -17,8 +17,8 @@ # along with gprMax. If not, see . import logging -from string import Template from pathlib import Path +from string import Template import h5py diff --git a/gprMax/fractals.py b/gprMax/fractals.py index 0fb8966b..14fee38b 100644 --- a/gprMax/fractals.py +++ b/gprMax/fractals.py @@ -16,12 +16,11 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . +import gprMax.config as config import numpy as np from scipy import fftpack -import gprMax.config as config -from .cython.fractals_generate import generate_fractal2D -from .cython.fractals_generate import generate_fractal3D +from .cython.fractals_generate import generate_fractal2D, generate_fractal3D from .utilities import round_value np.seterr(divide='raise') @@ -220,7 +219,7 @@ class Grass: self.numblades = numblades self.geometryparams = np.zeros((self.numblades, 6), - dtype=config.dtypes['float_or_double']) + dtype=config.sim_config.dtypes['float_or_double']) self.seed = None # Randomly defined parameters that will be used to calculate geometry diff --git a/gprMax/geometry_outputs.py b/gprMax/geometry_outputs.py index 57ff6fdd..a49b53b6 100644 --- a/gprMax/geometry_outputs.py +++ b/gprMax/geometry_outputs.py @@ -19,16 +19,15 @@ import logging import os from pathlib import Path - -import h5py -import numpy as np from struct import pack import gprMax.config as config +import h5py +import numpy as np from ._version import __version__ -from .cython.geometry_outputs import define_normal_geometry -from .cython.geometry_outputs import define_fine_geometry +from .cython.geometry_outputs import (define_fine_geometry, + define_normal_geometry) from .utilities import round_value logger = logging.getLogger(__name__) diff --git a/gprMax/gprMax.py b/gprMax/gprMax.py index b9dbaa15..76bbcc7e 100644 --- a/gprMax/gprMax.py +++ b/gprMax/gprMax.py @@ -20,8 +20,8 @@ import argparse import logging import gprMax.config as config -from .contexts import Context -from .contexts import MPIContext + +from .contexts import Context, MPIContext from .utilities import setup_logging logger = logging.getLogger(__name__) diff --git a/gprMax/grid.py b/gprMax/grid.py index 715efac6..7d41bc8c 100644 --- a/gprMax/grid.py +++ b/gprMax/grid.py @@ -16,22 +16,16 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from collections import OrderedDict import decimal as d - -from colorama import init -from colorama import Fore -from colorama import Style -init() -import numpy as np +from collections import OrderedDict import gprMax.config as config +import numpy as np +from colorama import Fore, Style, init +init() from .exceptions import GeneralError -from .pml import PML -from .pml import CFS -from .utilities import fft_power -from .utilities import human_size -from .utilities import round_value +from .pml import CFS, PML +from .utilities import fft_power, human_size, round_value np.seterr(invalid='raise') @@ -274,13 +268,13 @@ class FDTDGrid: def calculate_dt(self): """Calculate time step at the CFL limit.""" - if self.mode == '2D TMx': + if config.get_model_config().mode == '2D TMx': self.dt = 1 / (config.sim_config.em_consts['c'] * np.sqrt((1 / self.dy**2) + (1 / self.dz**2))) - elif self.mode == '2D TMy': + elif config.get_model_config().mode == '2D TMy': self.dt = 1 / (config.sim_config.em_consts['c'] * np.sqrt((1 / self.dx**2) + (1 / self.dz**2))) - elif self.mode == '2D TMz': + elif config.get_model_config().mode == '2D TMz': self.dt = 1 / (config.sim_config.em_consts['c'] * np.sqrt((1 / self.dx**2) + (1 / self.dy**2))) else: @@ -352,9 +346,12 @@ def dispersion_analysis(G): results (dict): Results from dispersion analysis """ - # Physical phase velocity error (percentage); grid sampling density; - # material with maximum permittivity; maximum significant frequency; error message - results = {'deltavp': False, 'N': False, 'material': False, 'maxfreq': [], 'error': ''} + # deltavp: physical phase velocity error (percentage) + # N: grid sampling density + # material: material with maximum permittivity + # maxfreq: maximum significant frequency + # error: error message + results = {'deltavp': None, 'N': None, 'material': None, 'maxfreq': [], 'error': ''} # Find maximum significant frequency if G.waveforms: diff --git a/gprMax/hash_cmds_file.py b/gprMax/hash_cmds_file.py index bc5900f1..b00843f2 100644 --- a/gprMax/hash_cmds_file.py +++ b/gprMax/hash_cmds_file.py @@ -16,13 +16,14 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from io import StringIO import logging import os -from pathlib import Path import sys +from io import StringIO +from pathlib import Path import gprMax.config as config + from .exceptions import CmdInputError from .hash_cmds_geometry import process_geometrycmds from .hash_cmds_multiuse import process_multicmds diff --git a/gprMax/hash_cmds_geometry.py b/gprMax/hash_cmds_geometry.py index 3b54b98f..275a123f 100644 --- a/gprMax/hash_cmds_geometry.py +++ b/gprMax/hash_cmds_geometry.py @@ -18,17 +18,18 @@ import numpy as np -from .cmds_geometry.edge import Edge -from .cmds_geometry.plate import Plate -from .cmds_geometry.triangle import Triangle +from .cmds_geometry.add_grass import AddGrass +from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness +from .cmds_geometry.add_surface_water import AddSurfaceWater from .cmds_geometry.box import Box from .cmds_geometry.cylinder import Cylinder from .cmds_geometry.cylindrical_sector import CylindricalSector +from .cmds_geometry.edge import Edge from .cmds_geometry.fractal_box import FractalBox +from .cmds_geometry.plate import Plate from .cmds_geometry.sphere import Sphere -from .cmds_geometry.add_surface_roughness import AddSurfaceRoughness -from .cmds_geometry.add_surface_water import AddSurfaceWater -from .cmds_geometry.add_grass import AddGrass +from .cmds_geometry.triangle import Triangle +from .exceptions import CmdInputError from .utilities import round_value @@ -181,15 +182,18 @@ def process_geometrycmds(geometry): # Isotropic case with no user specified averaging if len(tmp) == 10: - CylindricalSector(normal=normal, ctl1=ctl1, ctl2=ctl2, extent1=extent1, extent2=extent2, r=r, start=start, end=end, msterial_id=tmp[9]) + cylindrical_sector = CylindricalSector(normal=normal, ctr1=ctr1, ctr2=ctr2, extent1=extent1, + extent2=extent2, r=r, start=start, end=end, msterial_id=tmp[9]) # Isotropic case with user specified averaging elif len(tmp) == 11: - CylindricalSector(normal=normal, ctl1=ctl1, ctl2=ctl2, extent1=extent1, extent2=extent2, r=r, start=start, end=end, averaging=tmp[10], material_id=tmp[9]) + cylindrical_sector = CylindricalSector(normal=normal, ctr1=ctr1, ctr2=ctr2, extent1=extent1, extent2=extent2, + r=r, start=start, end=end, averaging=tmp[10], material_id=tmp[9]) # Uniaxial anisotropic case elif len(tmp) == 12: - CylindricalSector(normal=normal, ctl1=ctl1, ctl2=ctl2, extent1=extent1, extent2=extent2, r=r, start=start, end=end, material_ids=tmp[9:]) + cylindrical_sector = CylindricalSector(normal=normal, ctr1=ctr1, ctr2=ctr2, extent1=extent1, + extent2=extent2, r=r, start=start, end=end, material_ids=tmp[9:]) else: raise CmdInputError("'" + ' '.join(tmp) + "'" + ' too many parameters have been given') diff --git a/gprMax/hash_cmds_multiuse.py b/gprMax/hash_cmds_multiuse.py index 1944da3d..d4ccd501 100644 --- a/gprMax/hash_cmds_multiuse.py +++ b/gprMax/hash_cmds_multiuse.py @@ -16,21 +16,11 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from .cmds_multiple import Waveform -from .cmds_multiple import VoltageSource -from .cmds_multiple import HertzianDipole -from .cmds_multiple import MagneticDipole -from .cmds_multiple import TransmissionLine -from .cmds_multiple import Material -from .cmds_multiple import Snapshot -from .cmds_multiple import AddDebyeDispersion -from .cmds_multiple import AddLorentzDispersion -from .cmds_multiple import AddDrudeDispersion -from .cmds_multiple import SoilPeplinski -from .cmds_multiple import GeometryView -from .cmds_multiple import GeometryObjectsWrite -from .cmds_multiple import PMLCFS -from .cmds_multiple import Rx +from .cmds_multiple import (PMLCFS, AddDebyeDispersion, AddDrudeDispersion, + AddLorentzDispersion, GeometryObjectsWrite, + GeometryView, HertzianDipole, MagneticDipole, + Material, Rx, RxArray, Snapshot, SoilPeplinski, + TransmissionLine, VoltageSource, Waveform) from .exceptions import CmdInputError @@ -275,9 +265,9 @@ def process_multicmds(multicmds): if len(tmp) != 7: raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' requires exactly seven parameters') - p1 = float(tmp[0]), float(tmp[1]), float(tmp[2]) - p2 = float(tmp[3]), float(tmp[4]), float(tmp[5]) - gow = GeometryObjectsWrite(p1=p1, p2=p2, filename=tmp[6]) - scene_objects.append(gow) + p1 = float(tmp[0]), float(tmp[1]), float(tmp[2]) + p2 = float(tmp[3]), float(tmp[4]), float(tmp[5]) + gow = GeometryObjectsWrite(p1=p1, p2=p2, filename=tmp[6]) + scene_objects.append(gow) return scene_objects diff --git a/gprMax/hash_cmds_singleuse.py b/gprMax/hash_cmds_singleuse.py index 30b22db2..5d1727bd 100644 --- a/gprMax/hash_cmds_singleuse.py +++ b/gprMax/hash_cmds_singleuse.py @@ -16,17 +16,10 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from .cmds_single_use import Title -from .cmds_single_use import NumThreads -from .cmds_single_use import Discretisation -from .cmds_single_use import Domain -from .cmds_single_use import TimeStepStabilityFactor -from .cmds_single_use import TimeWindow -from .cmds_single_use import PMLCells -from .cmds_single_use import SrcSteps -from .cmds_single_use import RxSteps -from .cmds_single_use import ExcitationFile -from .cmds_single_use import OutputDir +from .cmds_single_use import (Discretisation, Domain, ExcitationFile, + NumThreads, OutputDir, PMLCells, RxSteps, + SrcSteps, TimeStepStabilityFactor, TimeWindow, + Title) from .exceptions import CmdInputError diff --git a/gprMax/model_build_run.py b/gprMax/model_build_run.py index 907c6203..8ab05dbc 100644 --- a/gprMax/model_build_run.py +++ b/gprMax/model_build_run.py @@ -19,37 +19,31 @@ import datetime import itertools import logging -from pathlib import Path import platform -import psutil import sys +from pathlib import Path -from colorama import init -from colorama import Fore -from colorama import Style -init() +import gprMax.config as config import numpy as np +import psutil +from colorama import Fore, Style, init +init() from terminaltables import SingleTable from tqdm import tqdm -import gprMax.config as config -from .cython.yee_cell_build import build_electric_components -from .cython.yee_cell_build import build_magnetic_components +from .cython.yee_cell_build import (build_electric_components, + build_magnetic_components) from .exceptions import GeneralError from .fields_outputs import write_hdf5_outputfile from .grid import dispersion_analysis from .hash_cmds_file import parse_hash_commands -from .materials import Material -from .materials import process_materials -from .pml import build_pml -from .pml import print_pml_info +from .materials import Material, process_materials +from .pml import build_pml, print_pml_info from .scene import Scene from .snapshots import Snapshot from .solvers import create_solver -from .utilities import get_terminal_width -from .utilities import human_size -from .utilities import mem_check_all -from .utilities import set_omp_threads +from .utilities import (get_terminal_width, human_size, mem_check_all, + set_omp_threads) logger = logging.getLogger(__name__) diff --git a/gprMax/mpi.py b/gprMax/mpi.py index 9771cd76..63b0afff 100644 --- a/gprMax/mpi.py +++ b/gprMax/mpi.py @@ -16,9 +16,9 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from enum import IntEnum import logging import time +from enum import IntEnum from mpi4py import MPI diff --git a/gprMax/pml.py b/gprMax/pml.py index 366ff3e5..0fa14e9c 100644 --- a/gprMax/pml.py +++ b/gprMax/pml.py @@ -18,9 +18,8 @@ from importlib import import_module -import numpy as np - import gprMax.config as config +import numpy as np class CFSParameter: diff --git a/gprMax/receivers.py b/gprMax/receivers.py index c4d0df9e..c34de534 100644 --- a/gprMax/receivers.py +++ b/gprMax/receivers.py @@ -16,9 +16,8 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -import numpy as np - import gprMax.config as config +import numpy as np class Rx: diff --git a/gprMax/scene.py b/gprMax/scene.py index fc215a97..9fa1aa03 100644 --- a/gprMax/scene.py +++ b/gprMax/scene.py @@ -21,11 +21,9 @@ import logging from .cmds_geometry.cmds_geometry import UserObjectGeometry from .cmds_geometry.fractal_box_builder import FractalBoxBuilder from .cmds_multiple import UserObjectMulti -from .cmds_single_use import UserObjectSingle -from .cmds_single_use import Domain -from .cmds_single_use import Discretisation -from .cmds_single_use import TimeWindow -from .exceptions import CmdInputError +from .cmds_single_use import (Discretisation, Domain, TimeWindow, + UserObjectSingle) +from .exceptions import CmdInputError, GeneralError from .materials import create_built_in_materials from .subgrids.user_objects import SubGridBase as SubGridUserBase from .user_inputs import create_user_input_points diff --git a/gprMax/snapshots.py b/gprMax/snapshots.py index 972e2231..1af7a57b 100644 --- a/gprMax/snapshots.py +++ b/gprMax/snapshots.py @@ -16,12 +16,12 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from struct import pack import sys - -import numpy as np +from struct import pack import gprMax.config as config +import numpy as np + from .cython.snapshots import calculate_snapshot_fields from .utilities import round_value diff --git a/gprMax/solvers.py b/gprMax/solvers.py index 16fa1d4f..31e71c57 100644 --- a/gprMax/solvers.py +++ b/gprMax/solvers.py @@ -17,11 +17,10 @@ # along with gprMax. If not, see . import gprMax.config as config -from .grid import FDTDGrid -from .grid import CUDAGrid + +from .grid import CUDAGrid, FDTDGrid from .subgrids.updates import create_updates as create_subgrid_updates -from .updates import CPUUpdates -from .updates import CUDAUpdates +from .updates import CPUUpdates, CUDAUpdates def create_G(): diff --git a/gprMax/sources.py b/gprMax/sources.py index eb2376e6..952ce9d3 100644 --- a/gprMax/sources.py +++ b/gprMax/sources.py @@ -18,12 +18,10 @@ from copy import deepcopy +import gprMax.config as config import numpy as np -import gprMax.config as config -from .grid import Ix -from .grid import Iy -from .grid import Iz +from .grid import Ix, Iy, Iz from .utilities import round_value diff --git a/gprMax/subgrids/subgrid_hsg.py b/gprMax/subgrids/subgrid_hsg.py index c85784fc..5758cb7b 100644 --- a/gprMax/subgrids/subgrid_hsg.py +++ b/gprMax/subgrids/subgrid_hsg.py @@ -16,16 +16,16 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from colorama import init -from colorama import Fore -from colorama import Style +import gprMax.config as config +from colorama import Fore, Style, init + +from ..cython.fields_updates_hsg import (cython_update_electric_os, + cython_update_is, + cython_update_magnetic_os) +from .base import SubGridBase + init() -import gprMax.config as config -from .base import SubGridBase -from ..cython.fields_updates_hsg import cython_update_is -from ..cython.fields_updates_hsg import cython_update_magnetic_os -from ..cython.fields_updates_hsg import cython_update_electric_os class SubGridHSG(SubGridBase): diff --git a/gprMax/subgrids/updates.py b/gprMax/subgrids/updates.py index 594445c1..01267032 100644 --- a/gprMax/subgrids/updates.py +++ b/gprMax/subgrids/updates.py @@ -19,10 +19,9 @@ import logging from ..exceptions import GeneralError -from .precursor_nodes import PrecursorNodes -from .precursor_nodes import PrecursorNodesFiltered -from .subgrid_hsg import SubGridHSG from ..updates import CPUUpdates +from .precursor_nodes import PrecursorNodes, PrecursorNodesFiltered +from .subgrid_hsg import SubGridHSG log = logging.getLogger(__name__) diff --git a/gprMax/subgrids/user_objects.py b/gprMax/subgrids/user_objects.py index 8e6182c5..aa75b8ae 100644 --- a/gprMax/subgrids/user_objects.py +++ b/gprMax/subgrids/user_objects.py @@ -19,11 +19,10 @@ from copy import copy import numpy as np - from gprMax import config + from ..cmds_geometry.cmds_geometry import UserObjectGeometry -from ..cmds_multiple import UserObjectMulti -from ..cmds_multiple import Rx +from ..cmds_multiple import Rx, UserObjectMulti from ..exceptions import CmdInputError from .multi import ReferenceRx as ReferenceRxUser from .subgrid_hsg import SubGridHSG as SubGridHSGUser diff --git a/gprMax/updates.py b/gprMax/updates.py index 2c48cfd5..fcb42451 100644 --- a/gprMax/updates.py +++ b/gprMax/updates.py @@ -16,27 +16,26 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from importlib import import_module import logging - -import numpy as np +from importlib import import_module import gprMax.config as config +import numpy as np + from .cuda.fields_updates import kernel_template_fields from .cuda.snapshots import kernel_template_store_snapshot from .cuda.source_updates import kernel_template_sources -from .cython.fields_updates_normal import update_electric as update_electric_cpu -from .cython.fields_updates_normal import update_magnetic as update_magnetic_cpu -from .fields_outputs import store_outputs as store_outputs_cpu +from .cython.fields_updates_normal import \ + update_electric as update_electric_cpu +from .cython.fields_updates_normal import \ + update_magnetic as update_magnetic_cpu +from .exceptions import GeneralError from .fields_outputs import kernel_template_store_outputs -from .receivers import htod_rx_arrays -from .receivers import dtoh_rx_array -from .snapshots import Snapshot -from .snapshots import htod_snapshot_array -from .snapshots import dtoh_snapshot_array +from .fields_outputs import store_outputs as store_outputs_cpu +from .receivers import dtoh_rx_array, htod_rx_arrays +from .snapshots import Snapshot, dtoh_snapshot_array, htod_snapshot_array from .sources import htod_src_arrays -from .utilities import round32 -from .utilities import timer +from .utilities import human_size, round32, timer class CPUUpdates: @@ -692,7 +691,7 @@ class CUDAUpdates: # Copy data from any snapshots back to correct snapshot objects if self.grid.snapshots and not config.get_model_config().cuda['snapsgpu2cpu']: for i, snap in enumerate(self.grid.snapshots): - dtoh_snapshot_arra(self.snapEx_gpu.get(), + dtoh_snapshot_array(self.snapEx_gpu.get(), self.snapEy_gpu.get(), self.snapEz_gpu.get(), self.snapHx_gpu.get(), diff --git a/gprMax/user_inputs.py b/gprMax/user_inputs.py index 25778548..facf6df6 100644 --- a/gprMax/user_inputs.py +++ b/gprMax/user_inputs.py @@ -18,19 +18,17 @@ import logging -from colorama import init -from colorama import Fore -from colorama import Style -init() -import numpy as np - import gprMax.config as config +import numpy as np +from colorama import Fore, Style, init +init() from .exceptions import CmdInputError from .subgrids.base import SubGridBase from .utilities import round_value logger = logging.getLogger(__name__) + """Module contains classes to handle points supplied by a user. The classes implement a common interface such that geometry building objects such as box or triangle do not need to have any knowledge which grid to diff --git a/gprMax/utilities.py b/gprMax/utilities.py index da6457a2..180f248e 100644 --- a/gprMax/utilities.py +++ b/gprMax/utilities.py @@ -16,34 +16,30 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from contextlib import contextmanager import codecs import decimal as d import logging import os import platform -import psutil import re import subprocess -from shutil import get_terminal_size import sys import textwrap +from contextlib import contextmanager +from shutil import get_terminal_size + +import gprMax.config as config +import numpy as np +import psutil +from colorama import Fore, Style, init +init() +from .exceptions import GeneralError try: from time import thread_time as timer_fn except ImportError: from time import perf_counter as timer_fn - -from colorama import init -from colorama import Fore -from colorama import Style -init() -import numpy as np - -import gprMax.config as config -from .exceptions import GeneralError - logger = logging.getLogger(__name__) @@ -525,7 +521,7 @@ def detect_gpus(): else: deviceIDsavail = range(drv.Device.count()) - # Gather information about selected/detected GPUs + # Gather information about detected GPUs gpus = [] for ID in deviceIDsavail: gpu = GPU(deviceID=ID) @@ -535,24 +531,24 @@ def detect_gpus(): return gpus -def check_gpus(gpus): - """Check if requested Nvidia GPU(s) deviceID(s) exist. +# def check_gpus(gpus): +# """Check if requested Nvidia GPU(s) deviceID(s) exist. - Args: - gpus (list): List of GPU object(s). - """ +# Args: +# gpus (list): List of GPU object(s). +# """ - # Check if requested device ID(s) exist - for ID in deviceIDs: - if ID not in deviceIDsavail: - raise GeneralError(f'GPU with device ID {ID} does not exist') +# # Check if requested device ID(s) exist +# for ID in deviceIDs: +# if ID not in deviceIDsavail: +# raise GeneralError(f'GPU with device ID {ID} does not exist') - # Gather information about selected/detected GPUs - gpus = [] - for ID in deviceIDsavail: - gpu = GPU(deviceID=ID) - gpu.get_gpu_info(drv) - gpus.append(gpu) +# # Gather information about selected/detected GPUs +# gpus = [] +# for ID in deviceIDsavail: +# gpu = GPU(deviceID=ID) +# gpu.get_gpu_info(drv) +# gpus.append(gpu) def timer(): diff --git a/gprMax/waveforms.py b/gprMax/waveforms.py index 0bbcbf88..36eb643f 100644 --- a/gprMax/waveforms.py +++ b/gprMax/waveforms.py @@ -16,8 +16,11 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . +import logging import numpy as np +logger = logging.getLogger(__name__) + class Waveform: """Definitions of waveform shapes that can be used with sources.""" @@ -129,5 +132,5 @@ class Waveform: return ampvalue def __str__(self): - log.debug('Do we need this?') + logger.debug('Do we need this?') return f'Waveform: ID={self.ID}, type={self.type}, amp{self.amp}, freq={self.freq}' diff --git a/tests/analytical_solutions.py b/tests/analytical_solutions.py index 0ef302a1..4ce8610a 100644 --- a/tests/analytical_solutions.py +++ b/tests/analytical_solutions.py @@ -1,7 +1,6 @@ import numpy as np -from gprMax.config import c -from gprMax.config import e0 +import gprMax.config as config from gprMax.waveforms import Waveform @@ -52,52 +51,41 @@ def hertzian_dipole_fs(iterations, dt, dxdydz, rx): x = rx[0] y = rx[1] z = rx[2] - if z == 0: - sign_z = 1 - else: - sign_z = np.sign(z) # Coordinates of Rx for Ex FDTD component Ex_x = x + 0.5 * dx Ex_y = y Ex_z = z - 0.5 * dz Er_x = np.sqrt((Ex_x**2 + Ex_y**2 + Ex_z**2)) - tau_Ex = Er_x / c + tau_Ex = Er_x / config.sim_config.em_consts['c'] # Coordinates of Rx for Ey FDTD component Ey_x = x Ey_y = y + 0.5 * dy Ey_z = z - 0.5 * dz Er_y = np.sqrt((Ey_x**2 + Ey_y**2 + Ey_z**2)) - tau_Ey = Er_y / c + tau_Ey = Er_y / config.sim_config.em_consts['c'] # Coordinates of Rx for Ez FDTD component Ez_x = x Ez_y = y Ez_z = z Er_z = np.sqrt((Ez_x**2 + Ez_y**2 + Ez_z**2)) - tau_Ez = Er_z / c + tau_Ez = Er_z / config.sim_config.em_consts['c'] # Coordinates of Rx for Hx FDTD component Hx_x = x Hx_y = y + 0.5 * dy Hx_z = z Hr_x = np.sqrt((Hx_x**2 + Hx_y**2 + Hx_z**2)) - tau_Hx = Hr_x / c + tau_Hx = Hr_x / config.sim_config.em_consts['c'] # Coordinates of Rx for Hy FDTD component Hy_x = x + 0.5 * dx Hy_y = y Hy_z = z Hr_y = np.sqrt((Hy_x**2 + Hy_y**2 + Hy_z**2)) - tau_Hy = Hr_y / c - - # Coordinates of Rx for Hz FDTD component - Hz_x = x + 0.5 * dx - Hz_y = y + 0.5 * dy - Hz_z = z - 0.5 * dz - Hr_z = np.sqrt((Hz_x**2 + Hz_y**2 + Hz_z**2)) - tau_Hz = Hr_z / c + tau_Hy = Hr_y / config.sim_config.em_consts['c'] # Initialise fields fields = np.zeros((iterations, 6)) @@ -118,30 +106,24 @@ def hertzian_dipole_fs(iterations, dt, dxdydz, rx): f_Ez = w.calculate_value((timestep * dt) - tau_Ez, dt) * dl fdot_Ez = wdot.calculate_value((timestep * dt) - tau_Ez, dt) * dl - fint_Hx = wint.calculate_value((timestep * dt) - tau_Hx, dt) * dl f_Hx = w.calculate_value((timestep * dt) - tau_Hx, dt) * dl fdot_Hx = wdot.calculate_value((timestep * dt) - tau_Hx, dt) * dl - fint_Hy = wint.calculate_value((timestep * dt) - tau_Hy, dt) * dl f_Hy = w.calculate_value((timestep * dt) - tau_Hy, dt) * dl fdot_Hy = wdot.calculate_value((timestep * dt) - tau_Hy, dt) * dl - fint_Hz = wint.calculate_value((timestep * dt) - tau_Hz, dt) * dl - f_Hz = w.calculate_value((timestep * dt) - tau_Hz, dt) * dl - fdot_Hz = wdot.calculate_value((timestep * dt) - tau_Hz, dt) * dl - # Ex - fields[timestep, 0] = ((Ex_x * Ex_z) / (4 * np.pi * e0 * Er_x**5)) * (3 * (fint_Ex + (tau_Ex * f_Ex)) + (tau_Ex**2 * fdot_Ex)) + fields[timestep, 0] = ((Ex_x * Ex_z) / (4 * np.pi * config.sim_config.em_consts['e0'] * Er_x**5)) * (3 * (fint_Ex + (tau_Ex * f_Ex)) + (tau_Ex**2 * fdot_Ex)) # Ey try: tmp = Ey_y / Ey_x except ZeroDivisionError: tmp = 0 - fields[timestep, 1] = tmp * ((Ey_x * Ey_z) / (4 * np.pi * e0 * Er_y**5)) * (3 * (fint_Ey + (tau_Ey * f_Ey)) + (tau_Ey**2 * fdot_Ey)) + fields[timestep, 1] = tmp * ((Ey_x * Ey_z) / (4 * np.pi * config.sim_config.em_consts['e0'] * Er_y**5)) * (3 * (fint_Ey + (tau_Ey * f_Ey)) + (tau_Ey**2 * fdot_Ey)) # Ez - fields[timestep, 2] = (1 / (4 * np.pi * e0 * Er_z**5)) * ((2 * Ez_z**2 - (Ez_x**2 + Ez_y**2)) * (fint_Ez + (tau_Ez * f_Ez)) - (Ez_x**2 + Ez_y**2) * tau_Ez**2 * fdot_Ez) + fields[timestep, 2] = (1 / (4 * np.pi * config.sim_config.em_consts['e0'] * Er_z**5)) * ((2 * Ez_z**2 - (Ez_x**2 + Ez_y**2)) * (fint_Ez + (tau_Ez * f_Ez)) - (Ez_x**2 + Ez_y**2) * tau_Ez**2 * fdot_Ez) # Hx fields[timestep, 3] = - (Hx_y / (4 * np.pi * Hr_x**3)) * (f_Hx + (tau_Hx * fdot_Hx)) diff --git a/tests/test_experimental.py b/tests/test_experimental.py index 99235748..acc96189 100644 --- a/tests/test_experimental.py +++ b/tests/test_experimental.py @@ -17,13 +17,12 @@ # along with gprMax. If not, see . import argparse -from pathlib import Path import sys +from pathlib import Path import h5py -import numpy as np import matplotlib.pyplot as plt - +import numpy as np from gprMax.exceptions import CmdInputError """Plots a comparison of fields between given simulation output and experimental data files.""" diff --git a/tests/test_models.py b/tests/test_models.py index cc541b7b..eb580372 100644 --- a/tests/test_models.py +++ b/tests/test_models.py @@ -16,21 +16,21 @@ # You should have received a copy of the GNU General Public License # along with gprMax. If not, see . -from pathlib import Path import sys +from pathlib import Path -from colorama import init, Fore, Style -init() +import gprMax import h5py -import numpy as np import matplotlib.pyplot as plt +import numpy as np +from colorama import Fore, Style, init +init() +from gprMax.exceptions import GeneralError +from tests.analytical_solutions import hertzian_dipole_fs if sys.platform == 'linux': plt.switch_backend('agg') -import gprMax -from gprMax.exceptions import GeneralError -from tests.analytical_solutions import hertzian_dipole_fs """Compare field outputs diff --git a/tools/HPC_scripts/gprmax_omp_mpi.sh b/tools/HPC_scripts/gprmax_omp_mpi.sh index dc1d8599..3c511ed1 100644 --- a/tools/HPC_scripts/gprmax_omp_mpi.sh +++ b/tools/HPC_scripts/gprmax_omp_mpi.sh @@ -16,7 +16,7 @@ #$ -pe mpi 176 ### Job script name: -#$ -N gprmax_omp_mpi.sh +#$ -N gprmax_omp_mpi_no_spawn.sh ##################################################################################### ### Initialise environment module @@ -34,4 +34,4 @@ export OMP_NUM_THREADS=16 ### Run gprMax with input file cd $HOME/gprMax -python -m gprMax mymodel.in -n 10 -mpi 11 +mpirun -n 11 python -m gprMax mymodel.in -n 10 -mpi diff --git a/tools/HPC_scripts/gprmax_omp_mpi_no_spawn.sh b/tools/HPC_scripts/gprmax_omp_mpi_no_spawn.sh deleted file mode 100644 index 5cd91c3e..00000000 --- a/tools/HPC_scripts/gprmax_omp_mpi_no_spawn.sh +++ /dev/null @@ -1,37 +0,0 @@ -#!/bin/sh -##################################################################################### -### Change to current working directory: -#$ -cwd - -### Specify runtime (hh:mm:ss): -#$ -l h_rt=01:00:00 - -### Email options: -#$ -m ea -M joe.bloggs@email.com - -### Resource reservation: -#$ -R y - -### Parallel environment ($NSLOTS): -#$ -pe mpi 176 - -### Job script name: -#$ -N gprmax_omp_mpi_no_spawn.sh -##################################################################################### - -### Initialise environment module -. /etc/profile.d/modules.sh - -### Load and activate Anaconda environment for gprMax, i.e. Python 3 and required packages -module load anaconda -source activate gprMax - -### Load OpenMPI -module load openmpi - -### Set number of OpenMP threads per MPI task (each gprMax model) -export OMP_NUM_THREADS=16 - -### Run gprMax with input file -cd $HOME/gprMax -mpirun -n 11 python -m gprMax mymodel.in -n 10 --mpi-no-spawn diff --git a/tools/inputfile_old2new.py b/tools/inputfile_old2new.py index 0fca1f32..fd02a9f3 100644 --- a/tools/inputfile_old2new.py +++ b/tools/inputfile_old2new.py @@ -223,7 +223,7 @@ while(lindex < len(inputlines)): lindex += 1 elif cmdname == '#bowtie': - print("Command '{}', is no longer supported. You can create the bowtie shape using two triangle commands.".format(inputlines[lindex], replacement)) + print("Command '{}', is no longer supported. You can create the bowtie shape using two triangle commands.".format(inputlines[lindex])) inputlines.pop(lindex) elif cmdname == '#cylinder': diff --git a/tools/outputfiles_merge.py b/tools/outputfiles_merge.py index 475cade3..bfdaacde 100644 --- a/tools/outputfiles_merge.py +++ b/tools/outputfiles_merge.py @@ -23,8 +23,8 @@ from pathlib import Path import h5py import numpy as np - from gprMax._version import __version__ +from gprMax.exceptions import CmdInputError def get_output_data(filename, rxnumber, rxcomponent): diff --git a/tools/plot_Ascan.py b/tools/plot_Ascan.py index 6dbe200b..9a6a4520 100644 --- a/tools/plot_Ascan.py +++ b/tools/plot_Ascan.py @@ -20,10 +20,9 @@ import argparse from pathlib import Path import h5py -import numpy as np -import matplotlib.pyplot as plt import matplotlib.gridspec as gridspec - +import matplotlib.pyplot as plt +import numpy as np from gprMax.exceptions import CmdInputError from gprMax.receivers import Rx from gprMax.utilities import fft_power diff --git a/tools/plot_Bscan.py b/tools/plot_Bscan.py index 11ea0578..336cc0d0 100644 --- a/tools/plot_Bscan.py +++ b/tools/plot_Bscan.py @@ -20,10 +20,10 @@ import argparse from pathlib import Path import h5py -import numpy as np import matplotlib.pyplot as plt - +import numpy as np from gprMax.exceptions import CmdInputError + from .outputfiles_merge import get_output_data diff --git a/tools/plot_antenna_params.py b/tools/plot_antenna_params.py index 92558467..83f93a4d 100644 --- a/tools/plot_antenna_params.py +++ b/tools/plot_antenna_params.py @@ -20,10 +20,9 @@ import argparse from pathlib import Path import h5py -import numpy as np -import matplotlib.pyplot as plt import matplotlib.gridspec as gridspec - +import matplotlib.pyplot as plt +import numpy as np from gprMax.exceptions import CmdInputError @@ -412,10 +411,10 @@ 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 - savename1 = file.stem + '_tl_params' - savename1 = file.parent / savename1 - savename2 = file.stem + '_ant_params' - savename2 = file.parent / savename2 + savename1 = filename.stem + '_tl_params' + savename1 = filename.parent / savename1 + savename2 = filename.stem + '_ant_params' + savename2 = filename.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', diff --git a/tools/plot_antenna_params_UtEl.py b/tools/plot_antenna_params_UtEl.py new file mode 100644 index 00000000..5f24e7cb --- /dev/null +++ b/tools/plot_antenna_params_UtEl.py @@ -0,0 +1,429 @@ +# Copyright (C) 2015-2020: The University of Edinburgh +# Authors: Craig Warren and Antonis Giannopoulos +# +# This file is part of gprMax. +# +# gprMax is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# gprMax is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with gprMax. If not, see . + +import argparse +from pathlib import Path + +import h5py +import matplotlib.gridspec as gridspec +import matplotlib.pyplot as plt +import numpy as np +from gprMax.exceptions import CmdInputError +from scipy.io import loadmat + + +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. + + Args: + filename (string): Filename (including path) of output file. + tltxnumber (int): Transmitter antenna - transmission line number + tlrxnumber (int): Receiver antenna - transmission line number + rxnumber (int): Receiver antenna - output number + rxcomponent (str): Receiver antenna - output electric field component + + Returns: + antennaparams (dict): Antenna parameters. + """ + + # Open output file and read some attributes + file = Path(filename) + f = h5py.File(file, 'r') + dxdydz = f.attrs['dx_dy_dz'] + dt = f.attrs['dt'] + iterations = f.attrs['Iterations'] + + # Calculate time array and frequency bin spacing + time = np.linspace(0, (iterations - 1) * dt, num=iterations) + df = 1 / np.amax(time) + + 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) + '/' + + # Incident voltages/currents + Vinc = f[tltxpath + 'Vinc'][:] + Iinc = f[tltxpath + 'Iinc'][:] + + # Total (incident + reflected) voltages/currents + Vtotal = f[tltxpath + 'Vtotal'][:] + Itotal = f[tltxpath + 'Itotal'][:] + + # Reflected voltages/currents + Vref = Vtotal - Vinc + Iref = Itotal - Iinc + + # If a receiver antenna is used (with a transmission line or receiver), get received voltage for s21 + if tlrxnumber: + tlrxpath = '/tls/tl' + str(tlrxnumber) + '/' + Vrec = f[tlrxpath + 'Vtotal'][:] + + elif rxnumber: + rxpath = '/rxs/rx' + str(rxnumber) + '/' + availableoutputs = list(f[rxpath].keys()) + + if rxcomponent not in availableoutputs: + raise CmdInputError(f"{rxcomponent} output requested, but the available output for receiver {rxnumber} is {', '.join(availableoutputs)}") + + rxpath += rxcomponent + + # Received voltage + if rxcomponent == 'Ex': + Vrec = f[rxpath][:] * -1 * dxdydz[0] + elif rxcomponent == 'Ey': + Vrec = f[rxpath][:] * -1 * dxdydz[1] + elif rxcomponent == 'Ez': + Vrec = f[rxpath][:] * -1 * dxdydz[2] + f.close() + + # Frequency bins + freqs = np.fft.fftfreq(Vinc.size, d=dt) + + # Delay correction - current lags voltage, so delay voltage to match current timestep + delaycorrection = np.exp(1j * 2 * np.pi * freqs * (dt / 2)) + + # Calculate s11 and (optionally) s21 + s11 = np.abs(np.fft.fft(Vref) / np.fft.fft(Vinc)) + if tlrxnumber or rxnumber: + s21 = np.abs(np.fft.fft(Vrec) / np.fft.fft(Vinc)) + + # Calculate input impedance + zin = (np.fft.fft(Vtotal) * delaycorrection) / np.fft.fft(Itotal) + + # Calculate input admittance + yin = np.fft.fft(Itotal) / (np.fft.fft(Vtotal) * delaycorrection) + + # Convert to decibels (ignore warning from taking a log of any zero values) + with np.errstate(divide='ignore'): + Vincp = 20 * np.log10(np.abs((np.fft.fft(Vinc) * delaycorrection))) + Iincp = 20 * np.log10(np.abs(np.fft.fft(Iinc))) + Vrefp = 20 * np.log10(np.abs((np.fft.fft(Vref) * delaycorrection))) + Irefp = 20 * np.log10(np.abs(np.fft.fft(Iref))) + Vtotalp = 20 * np.log10(np.abs((np.fft.fft(Vtotal) * delaycorrection))) + Itotalp = 20 * np.log10(np.abs(np.fft.fft(Itotal))) + s11 = 20 * np.log10(s11) + + # Replace any NaNs or Infs from zero division + Vincp[np.invert(np.isfinite(Vincp))] = 0 + Iincp[np.invert(np.isfinite(Iincp))] = 0 + Vrefp[np.invert(np.isfinite(Vrefp))] = 0 + Irefp[np.invert(np.isfinite(Irefp))] = 0 + Vtotalp[np.invert(np.isfinite(Vtotalp))] = 0 + Itotalp[np.invert(np.isfinite(Itotalp))] = 0 + s11[np.invert(np.isfinite(s11))] = 0 + + # Create dictionary of antenna parameters + antennaparams = {'time': time, 'freqs': freqs, 'Vinc': Vinc, 'Vincp': Vincp, 'Iinc': Iinc, 'Iincp': Iincp, + 'Vref': Vref, 'Vrefp': Vrefp, 'Iref': Iref, 'Irefp': Irefp, + 'Vtotal': Vtotal, 'Vtotalp': Vtotalp, 'Itotal': Itotal, 'Itotalp': Itotalp, + 's11': s11, 'zin': zin, 'yin': yin} + if tlrxnumber or rxnumber: + with np.errstate(divide='ignore'): # Ignore warning from taking a log of any zero values + s21 = 20 * np.log10(s21) + s21[np.invert(np.isfinite(s21))] = 0 + antennaparams['s21'] = s21 + + return antennaparams + + +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. + + Args: + filename (string): Filename (including path) of output file. + time (array): Simulation time. + freq (array): Frequencies for FFTs. + Vinc, Vincp, Iinc, Iincp (array): Time and frequency domain representations of incident voltage and current. + Vref, Vrefp, Iref, Irefp (array): Time and frequency domain representations of reflected voltage and current. + Vtotal, Vtotalp, Itotal, Itotalp (array): Time and frequency domain representations of total voltage and current. + s11, s21 (array): s11 and, optionally, s21 parameters. + zin, yin (array): Input impedance and input admittance parameters. + + Returns: + plt (object): matplotlib plot object. + """ + + # Set plotting range + pltrangemin = 1 + # 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 > 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(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') + # gs1 = gridspec.GridSpec(4, 2, hspace=0.7) + # ax = plt.subplot(gs1[0, 0]) + # ax.plot(time, Vinc, 'r', lw=2, label='Vinc') + # ax.set_title('Incident voltage') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Voltage [V]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot frequency spectra of incident voltage + # ax = plt.subplot(gs1[0, 1]) + # 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') + # ax.plot(freqs[pltrange], Vincp[pltrange], 'r', lw=2) + # ax.set_title('Incident voltage') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot incident current + # ax = plt.subplot(gs1[1, 0]) + # ax.plot(time, Iinc, 'b', lw=2, label='Vinc') + # ax.set_title('Incident current') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Current [A]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot frequency spectra of incident current + # ax = plt.subplot(gs1[1, 1]) + # 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') + # ax.plot(freqs[pltrange], Iincp[pltrange], 'b', lw=2) + # ax.set_title('Incident current') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot total voltage + # ax = plt.subplot(gs1[2, 0]) + # ax.plot(time, Vtotal, 'r', lw=2, label='Vinc') + # ax.set_title('Total (incident + reflected) voltage') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Voltage [V]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot frequency spectra of total voltage + # ax = plt.subplot(gs1[2, 1]) + # 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') + # ax.plot(freqs[pltrange], Vtotalp[pltrange], 'r', lw=2) + # ax.set_title('Total (incident + reflected) voltage') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot total current + # ax = plt.subplot(gs1[3, 0]) + # ax.plot(time, Itotal, 'b', lw=2, label='Vinc') + # ax.set_title('Total (incident + reflected) current') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Current [A]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot frequency spectra of total current + # ax = plt.subplot(gs1[3, 1]) + # 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') + # ax.plot(freqs[pltrange], Itotalp[pltrange], 'b', lw=2) + # ax.set_title('Total (incident + reflected) current') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot reflected (reflected) voltage + # ax = plt.subplot(gs1[4, 0]) + # ax.plot(time, Vref, 'r', lw=2, label='Vref') + # ax.set_title('Reflected voltage') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Voltage [V]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot frequency spectra of reflected voltage + # ax = plt.subplot(gs1[4, 1]) + # 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') + # ax.plot(freqs[pltrange], Vrefp[pltrange], 'r', lw=2) + # ax.set_title('Reflected voltage') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot reflected (reflected) current + # ax = plt.subplot(gs1[5, 0]) + # ax.plot(time, Iref, 'b', lw=2, label='Iref') + # ax.set_title('Reflected current') + # ax.set_xlabel('Time [s]') + # ax.set_ylabel('Current [A]') + # ax.set_xlim([0, np.amax(time)]) + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot frequency spectra of reflected current + # ax = plt.subplot(gs1[5, 1]) + # 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') + # ax.plot(freqs[pltrange], Irefp[pltrange], 'b', lw=2) + # ax.set_title('Reflected current') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Power [dB]') + # ax.grid(which='both', axis='both', linestyle='-.') + + # Figure 2 + # Load measured s11 + s11_meas = loadmat('/Volumes/Users/cwarren/Dropbox (Northumbria University)/Research/GPR-antennas/Utsi Electronics/500MHz_Vivaldi_antenna/measured/s11_measured.mat') + + # Plot frequency spectra of s11 + fig2, ax = plt.subplots(num='Antenna parameters', figsize=(20, 12), facecolor='w', edgecolor='w') + gs2 = gridspec.GridSpec(1, 1, hspace=0.3) + ax = plt.subplot(gs2[0, 0]) + markerline, stemlines, baseline = ax.stem(freqs[pltrange], s11[pltrange], '-.', use_line_collection=True) + markerline2, stemlines2, baseline2 = ax.stem(s11_meas['f'], s11_meas['s11_t'], '-.', use_line_collection=True) + plt.setp(baseline, 'linewidth', 0) + plt.setp(stemlines, 'color', 'w') + plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g') + plt.setp(baseline2, 'linewidth', 0) + plt.setp(stemlines2, 'color', 'w') + plt.setp(markerline2, 'markerfacecolor', 'b', 'markeredgecolor', 'b') + ax.plot(freqs[pltrange], s11[pltrange], 'g', lw=2, label='Simulation (Debye, 50$\Omega$)') + ax.plot(s11_meas['f'], s11_meas['s11_t'], 'b', lw=2, label='Measured') + ax.set_title('s11') + ax.set_xlabel('Frequency [Hz]') + ax.set_ylabel('Power [dB]') + ax.set_xlim([0, 3e9]) + ax.set_ylim([-40, 0]) + ax.grid(which='both', axis='both', linestyle='-.') + ax.legend(fontsize=12, frameon=False) + + # 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], '-.', use_line_collection=True) + plt.setp(baseline, 'linewidth', 0) + plt.setp(stemlines, 'color', 'g') + plt.setp(markerline, 'markerfacecolor', 'g', 'markeredgecolor', 'g') + ax.plot(freqs[pltrange], s21[pltrange], 'g', lw=2) + ax.set_title('s21') + ax.set_xlabel('Frequency [Hz]') + ax.set_ylabel('Power [dB]') + # ax.set_xlim([0.88e9, 1.02e9]) + # ax.set_ylim([-25, 50]) + ax.grid(which='both', axis='both', linestyle='-.') + + # Plot input resistance (real part of impedance) + # ax = plt.subplot(gs2[1, 0]) + # 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') + # ax.plot(freqs[pltrange], zin[pltrange].real, 'g', lw=2) + # ax.set_title('Input impedance (resistive)') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Resistance [Ohms]') + # # ax.set_xlim([0.88e9, 1.02e9]) + # ax.set_ylim(bottom=0) + # # ax.set_ylim([0, 300]) + # ax.grid(which='both', axis='both', linestyle='-.') + # + # # Plot input reactance (imaginery part of impedance) + # ax = plt.subplot(gs2[1, 1]) + # 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') + # ax.plot(freqs[pltrange], zin[pltrange].imag, 'g', lw=2) + # ax.set_title('Input impedance (reactive)') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Reactance [Ohms]') + # # ax.set_xlim([0.88e9, 1.02e9]) + # # ax.set_ylim([-300, 300]) + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot input admittance (magnitude) + # ax = plt.subplot(gs2[2, 0]) + # 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') + # ax.plot(freqs[pltrange], np.abs(yin[pltrange]), 'g', lw=2) + # ax.set_title('Input admittance (magnitude)') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Admittance [Siemens]') + # ax.set_xlim([0.88e9, 1.02e9]) + # ax.set_ylim([0, 0.035]) + # ax.grid(which='both', axis='both', linestyle='-.') + + # Plot input admittance (phase) + # ax = plt.subplot(gs2[2, 1]) + # 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') + # ax.plot(freqs[pltrange], np.angle(yin[pltrange], deg=True), 'g', lw=2) + # ax.set_title('Input admittance (phase)') + # ax.set_xlabel('Frequency [Hz]') + # ax.set_ylabel('Phase [degrees]') + # ax.set_xlim([0.88e9, 1.02e9]) + # ax.set_ylim([-40, 100]) + # 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(filename.with_suffix('.pdf'), dpi=None, format='pdf', bbox_inches='tight', pad_inches=0.1) + + return plt + + +if __name__ == "__main__": + + # Parse command line arguments + parser = argparse.ArgumentParser(description='Plots antenna parameters (s11, s21 parameters and input impedance) from an output file containing a transmission line source.', usage='cd gprMax; python -m tools.plot_antenna_params outputfile') + parser.add_argument('outputfile', help='name of output file including path') + parser.add_argument('--tltx-num', default=1, type=int, help='transmitter antenna - transmission line number') + parser.add_argument('--tlrx-num', type=int, help='receiver antenna - transmission line number') + parser.add_argument('--rx-num', type=int, help='receiver antenna - output number') + parser.add_argument('--rx-component', type=str, help='receiver antenna - output electric field component', choices=['Ex', 'Ey', 'Ez']) + args = parser.parse_args() + + antennaparams = calculate_antenna_params(args.outputfile, args.tltx_num, args.tlrx_num, args.rx_num, args.rx_component) + plthandle = mpl_plot(args.outputfile, **antennaparams) + plthandle.show() diff --git a/tools/plot_source_wave.py b/tools/plot_source_wave.py index 5e51b7bd..7d7299cd 100644 --- a/tools/plot_source_wave.py +++ b/tools/plot_source_wave.py @@ -19,12 +19,10 @@ import argparse from pathlib import Path -import numpy as np import matplotlib.pyplot as plt - +import numpy as np from gprMax.exceptions import CmdInputError -from gprMax.utilities import fft_power -from gprMax.utilities import round_value +from gprMax.utilities import fft_power, round_value from gprMax.waveforms import Waveform diff --git a/user_libs/antenna_patterns/initial_save.py b/user_libs/antenna_patterns/initial_save.py index 38538d9c..b0fc312d 100644 --- a/user_libs/antenna_patterns/initial_save.py +++ b/user_libs/antenna_patterns/initial_save.py @@ -9,13 +9,10 @@ import argparse import os import sys +import gprMax.config as config import h5py -import numpy as np import matplotlib.pyplot as plt - -from gprMax.config import c -from gprMax.config import z0 - +import numpy as np # Parse command line arguments parser = argparse.ArgumentParser(description='Calculate and store (in a Numpy file) field patterns from a simulation with receivers positioned in circles around an antenna.', usage='cd gprMax; python -m user_libs.antenna_patterns.initial_save outputfile') @@ -55,9 +52,9 @@ traceno = np.s_[:] # All traces # Critical angle and velocity if epsr: mr = 1 - z1 = np.sqrt(mr / epsr) * z0 - v1 = c / np.sqrt(epsr) - thetac = np.round(np.arcsin(v1 / c) * (180 / np.pi)) + z1 = np.sqrt(mr / epsr) * config.sim_config.em_consts['z0'] + v1 = config.sim_config.em_consts['c'] / np.sqrt(epsr) + thetac = np.round(np.arcsin(v1 / config.sim_config.em_consts['c']) * (180 / np.pi)) wavelength = v1 / f # Print some useful information @@ -164,8 +161,8 @@ for radius in range(0, len(radii)): Ethetasum[index] = np.sum(Etheta[:, index]**2) / z1 Hthetasum[index] = np.sum(Htheta[:, index]**2) / z1 else: - Ethetasum[index] = np.sum(Etheta[:, index]**2) / z0 - Hthetasum[index] = np.sum(Htheta[:, index]**2) / z0 + Ethetasum[index] = np.sum(Etheta[:, index]**2) / config.sim_config.em_consts['z0'] + Hthetasum[index] = np.sum(Htheta[:, index]**2) / config.sim_config.em_consts['z0'] index += 1 diff --git a/user_libs/antenna_patterns/plot_fields.py b/user_libs/antenna_patterns/plot_fields.py index f78247ad..e3ff3351 100644 --- a/user_libs/antenna_patterns/plot_fields.py +++ b/user_libs/antenna_patterns/plot_fields.py @@ -9,12 +9,9 @@ import argparse import os import sys -import numpy as np +import gprMax.config as config import matplotlib.pyplot as plt - -from gprMax.config import c -from gprMax.config import z0 - +import numpy as np # Parse command line arguments parser = argparse.ArgumentParser(description='Plot field patterns from a simulation with receivers positioned in circles around an antenna. This module should be used after the field pattern data has been processed and stored using the initial_save.py module.', usage='cd gprMax; python -m user_libs.antenna_patterns.plot_fields numpyfile') @@ -52,9 +49,9 @@ step = 12 # Critical angle and velocity if epsr: mr = 1 - z1 = np.sqrt(mr / epsr) * z0 - v1 = c / np.sqrt(epsr) - thetac = np.round(np.rad2deg(np.arcsin(v1 / c))) + z1 = np.sqrt(mr / epsr) * config.sim_config.em_consts['z0'] + v1 = config.sim_config.em_consts['c'] / np.sqrt(epsr) + thetac = np.round(np.rad2deg(np.arcsin(v1 / config.sim_config.em_consts['c']))) wavelength = v1 / f # Print some useful information