# -*- coding: utf-8 -*-
"""
@author: Rochlitz.R
"""
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as clrs
from custEM.misc import pyhed_calculations as phc
from custEM.post import plot_utils as pu
from custEM.post.plot_utils import PlotBase
import custEM.misc
import sys
import os
from matplotlib import cm
from matplotlib.colors import LogNorm
from matplotlib import rcParams
[docs]
class PlotFD(PlotBase):
"""
Plot class for visualization of custEM results.
class internal functions:
-------------------------
- import_line_data()
load data interpolated on lines.
- import_slice_data()
load data interpolated on slices.
- load_default_line_data()
import all data interpolated on default coordinate-axes.
- load_default_slice_data()
import all data interpolated on default slices orthogonal to the
coordinate axes.
- plot_line_data()
plot either total or secondary **E** or **H** fields component-wise
with real and imaginary parts on an observation line.
- plot_line_errors()
plot mismatches between two datasets interpolated on the same line
or between data and reference solutions component-wise or magnitude
with real and imaginary parts on an observation line.
- plot_slice_data()
plot either total or secondary **E** or **H** fields component-wise
with real and imaginary parts on an observation slice.
- plot_slice_errors()
plot mismatches between two datasets interpolated on the same slice
or between data and reference solutions component-wise or magnitude
with real and imaginary parts on an observation slice.
- add_to_slice_plot()
plot 2D data in an existing (sub)figure for custom illustrations.
- add_errors_to_slice_plot()
plot mismatches of components or field magnitudes of 2D data
in an existing (sub)figure for custom illustrations.
- add_3D_slice_plot()
adds visulaization of 2D data (e.g., with topography) to an
existing (sub)figure as 3D plot.
class internal utility functions:
---------------------------------
- init_main_parameters()
utility function for the import-functions.
- init_component_integers()
utility function for evaulating the correct coordinate directions.
- calc_pyhed_reference()
deprecated, might not work anymore, might be re-implemented soon.
- eval_properties()
evaluate plot-properties.
- reduce_slice_data()
if slice data very interpolated on a grid that is too fine, one
could apply a stride (default *1*) to reduce the amount of data,
otherwise, 2D visualization can become very slow or even crash.
- print_import_error()
print warning, if specified data could not be imported
- get_coords()
evaluate correct coordinates for 1D or 2D visualization
- get_coord_names()
utility function for **get_coords()**
- init_cmap()
initalize colormap (adjust colormap type and range)
- add_cbar()
add a proper colorbar to a figure
- rel_errors(), abs_errors(), ratio()
calculate relative, absolute errors or ratio between two datasets
using the absolute ('abs') values of the latter
- mean_errors(), median_errors()
calculate **mean* or *median* of an error distribution
using the absolute ('abs') values of the latter
- save_plot()
save a plot in the save-directory
- further class-external utility functions are defined at the bottom!
see description at the bottom
"""
def __init__(self, mod, mesh, approach, **plot_kwargs):
"""
Class for importing, comparing and illustrating interpolated FEM
results obtained with custEM in 1D, 2D and 3D. So far, lines and slices
parallel or orthogonal to the coordinate-system axes are supported,
even with topography. However, support for arbitrary directed
interpolation paths or planes still needs to be implemented.
Note:
-----
- Field quantities are referred to as **E_t**, **H_t**, **E_s**
or **H_s**. Potential illustration is not supported yet.
Required arguments:
-------------------
- mod_name, type str
name of the model
- mesh_name, type str
name of input mesh without suffix
- approach, type str
either **E_t**, **H_t**, **Am_t**, **An_t** or
**E_s**, **H_s**, **Am_s**, **An_s**. Note that **H_t** and
**An_t** are theoretically considered but not working yet.
Keyword arguments:
------------------
- r_dir = '../results/', type str
specify custom modeling results directory if required.
- s_dir = '.', type str
specify custom save-directory for figures if required.
- fig_size = None, type list with len = 2
global adjustment of figure sizes [inch].
- dpi = 300, type int
dpi value for saving figures in "png" format
- label_color = '#000000', type str
string of hexahedral color code to define axis and label colors
for plots
- fs = 12, type int
default font size for tick-labels, legends etc.
- dg_space = None, type bool
set *True* if data are loaded from discontinuous interpolation
to automatically remove multiple values at the same positions and
sort the data. Please note, this approach is experimental and
should be used with caution
"""
super(self.__class__, self).__init__()
self.mod = mod
self.mesh = mesh
self.approach = approach
self.load_parameters = True
for key in plot_kwargs:
if key not in self.__dict__:
print('!!! Warning: Unknown Plot class argument set:', key)
self.__dict__.update(**plot_kwargs)
rcParams['axes.edgecolor'] = self.label_color
rcParams['xtick.color'] = self.label_color
rcParams['ytick.color'] = self.label_color
rcParams['axes.labelcolor'] = self.label_color
rcParams['text.color'] = self.label_color
rcParams['axes.unicode_minus'] = False
if not os.path.isdir(self.s_dir):
os.makedirs(self.s_dir)
rcParams['axes.unicode_minus'] = False
if self.fig_size is not None:
rcParams['figure.figsize'] = self.fig_size[0], self.fig_size[1]
if self.load_parameters:
self.load_model_parameters()
###########################################################################
# # # import section # # #
###########################################################################
[docs]
def import_line_data(self, linE, mod=None, mesh=None, approach=None,
EH='EH', sf=None, path=None, key=None, stride=1,
tx=None, freq=None):
"""
import line data from regular line (line-name must end with either
"_x", "_y" or "_z" for x-, y-, z-directed lines, respectively).
Required arguments:
-------------------
- linE, type str
name of the line-mesh on which the required data were interpolated
- mod = None, type str
specify if a different model than the overall defined **self.mod**
(input when initiallizing Plot instance) should be imported
- mesh = None, type str
specify if a different mesh than the overall defined **self.mesh**
(Plot instance) should be imported
- approach = None, type str
specify if a different approach than the overall defined
**self.approach** (Plot instance) should be imported
- EH='EH', type str
Alternatively **E** or **H**, of solely electric or magentic fields
should be imported
- sf=False, type bool
Set **True**, if secondary fields should be imported instead of
total fields
- path = None, type str
if non-default interpolation directorys were chosen, specify the
path in which the interpolation results are located
- key = None, type str
specify a **key** which can be used all over the **Plot**
class functions to identify the imported data with this name;
Otherwise, the name will generated by default
*(mod + '_E_t_on_' + linE)*
- stride = 1, type int
downsample observation points; return only every *stride*-th point
- tx = None, type int
if multiple Tx have been computed in the same model, specify which
one to import
- freq = None, type int
if multiple frequencies have been computed in the same model,
specify which one to import
"""
mod, mesh, approach, path = self.init_main_parameters(
mod, mesh, approach, path)
comp = self.init_component_integers(linE[-1])
if freq is not None:
path += 'f_' + str(freq) + '_'
if tx is not None:
path += 'tx_' + str(tx) + '_'
if 'E' in EH:
self.general_import(path, 'E_t', linE, key, comp, stride)
if 'H' in EH:
self.general_import(path, 'H_t', linE, key, comp, stride)
if '_s' in approach and sf is None:
sf = True
if sf:
if 'E' in EH:
self.general_import(path, 'E_s', linE, key, comp, stride)
if 'H' in EH:
self.general_import(path, 'H_s', linE, key, comp, stride)
[docs]
def import_point_data(self, points, mod=None, mesh=None, approach=None,
EH='EH', sf=None, path=None, key=None,
tx=None, freq=None):
"""
import line data from regular line (line-name must end with either
"_x", "_y" or "_z" for x-, y-, z-directed lines, respectively).
Required arguments:
-------------------
- points, type str
name of the path-mesh on which the required data were interpolated
- mod = None, type str
specify if a different model than the overall defined **self.mod**
(input when initiallizing Plot instance) should be imported
- mesh = None, type str
specify if a different mesh than the overall defined **self.mesh**
(Plot instance) should be imported
- approach = None, type str
specify if a different approach than the overall defined
**self.approach** (Plot instance) should be imported
- EH='EH', type str
Alternatively **E** or **H**, of solely electric or magentic fields
should be imported
- sf=False, type bool
Set **True**, if secondary fields should be imported instead of
total fields
- path = None, type str
if non-default interpolation directorys were chosen, specify the
path in which the interpolation results are located
- key = None, type str
specify a **key** which can be used all over the **Plot**
class functions to identify the imported data with this name;
Otherwise, the name will generated by default
*(mod + '_E_t_on_' + linE)*
- tx = None, type int
if multiple Tx have been computed in the same model, specify which
one to import
- freq = None, type int
if multiple frequencies have been computed in the same model,
specify which one to import
"""
mod, mesh, approach, path = self.init_main_parameters(
mod, mesh, approach, path)
if freq is not None:
path += 'f_' + str(freq) + '_'
if tx is not None:
path += 'tx_' + str(tx) + '_'
if 'E' in EH:
self.general_import(path, 'E_t', points, key, None, stride)
if 'H' in EH:
self.general_import(path, 'H_t', points, key, None, stride)
if '_s' in approach and sf is None:
sf = True
if sf:
if 'E' in EH:
self.general_import(path, 'E_s', points, key, None, stride)
if 'H' in EH:
self.general_import(path, 'H_s', points, key, None, stride)
[docs]
def import_slice_data(self, slicE, mod=None, mesh=None, approach=None,
EH='EH', sf=None, path=None, key=None,
coord_key=None, stride=1, tx=None, freq=None):
"""
import slice data from regular slices (slice-name must end with either
"_x", "_y" or "_z" for x-, y-, z-directed slice-normals, respectively).
Required arguments:
-------------------
- slicE, type str
name of the slice-mesh on which required data were interpolated
- mod = None, type str
specify if a different model than the overall defined **self.mod**
(input when initiallizing Plot instance) should be imported
- mesh = None, type str
specify if a different mesh than the overall defined **self.mesh**
(Plot instance) should be imported
- approach = None, type str
specify if a different approach than the overall defined
**self.approach** (Plot instance) should be imported
- EH='EH', type str
Alternatively **E** or **H**, of solely electric or magentic fields
should be imported
- sf=False, type bool
Set **True**, if secondary fields should be imported instead of
total fields
- path = None, type str
if non-default interpolation directorys were chosen, specify the
path in which the interpolation results are located
- key = None, type str
specify a **key** which can be used all over the **Plot**
class functions to identify the imported data with this name;
Otherwise, the name will generated by default:
*(mod + '_E_t_on_' + linE)*
- coord_key = None, type str
specify a **coord_key**, if interpolated data are imported several
times with a different **stride**, e.g. for 'contour' and 'vector-
arorow '(plt.contourf and plt.quiver) visualization at the same
time. Later on, the underlying data coordinates can be referenced
for plotting by these **coord_key**s
- stride = 1, type int
specify, if only every *stride*th (e.g., 10th) data point from a
dense slice-interpolation mesh should be imported (stride equal
in x- and y- direction). Too large 2D
datasets make visualization VERY slow and even lead to crashes
- tx = None, type int
if multiple Tx have been computed in the same model, specify which
one to import
- freq = None, type int
if multiple frequencies have been computed in the same model,
specify which one to import
"""
mod, mesh, approach, path = self.init_main_parameters(
mod, mesh, approach, path)
comp1, comp2 = self.init_component_integers(slicE[-1], line=False)
if freq is not None:
path += 'f_' + str(freq) + '_'
if tx is not None:
path += 'tx_' + str(tx) + '_'
if coord_key is None:
coord_key = slicE
if 'E' in EH:
self.general_import(path, 'E_t', slicE, key, [comp1, comp2],
stride, coord_key)
if 'H' in EH:
self.general_import(path, 'H_t', slicE, key, [comp1, comp2],
stride, coord_key)
if '_s' in approach and sf is None:
sf = True
if sf:
if 'E' in EH:
self.general_import(path, 'E_s', slicE, key, [comp1, comp2],
stride, coord_key)
if 'H' in EH:
self.general_import(path, 'H_s', slicE, key, [comp1, comp2],
stride, coord_key)
[docs]
def load_default_line_data(self, interp_meshes='small', **kwargs):
"""
laod all or specific line data from a default set of
line-interpolation meshes.
Keyword arguments:
------------------
- interp_meshes = 'small', type str
name of the default set of line-interpolation meshes
- further keyword arguments listed in the **import_line_data()** docs
"""
if interp_meshes == 'small':
self.import_line_data('default_small_line_x', **kwargs)
self.import_line_data('default_small_line_y', **kwargs)
self.import_line_data('default_small_line_z', **kwargs)
elif interp_meshes == 'large':
self.import_line_data('default_large_line_x', **kwargs)
self.import_line_data('default_large_line_y', **kwargs)
self.import_line_data('default_large_line_z', **kwargs)
else:
print('Error!, Invalid interpolation mesh type!')
raise SystemExit
if '_s' in self.approach:
if interp_meshes == 'small':
self.import_line_data('default_small_line_x', **kwargs)
self.import_line_data('default_small_line_y', **kwargs)
self.import_line_data('default_small_line_z', **kwargs)
elif interp_meshes == 'large':
self.import_line_data('default_large_line_x', **kwargs)
self.import_line_data('default_large_line_y', **kwargs)
self.import_line_data('default_large_line_z', **kwargs)
[docs]
def load_default_slice_data(self, interp_meshes='small', **kwargs):
"""
laod all or specific line data from a default set of
slice-interpolation meshes.
Keyword arguments:
------------------
- interp_meshes = 'small', type str
name of the default set of slice-interpolation meshes
- further keyword arguments listed in the **import_slice_data()** docs
"""
if interp_meshes == 'small':
self.import_slice_data('default_small_slice_x', **kwargs)
self.import_slice_data('default_small_slice_y', **kwargs)
self.import_slice_data('default_small_slice_z', **kwargs)
elif interp_meshes == 'large':
self.import_slice_data('default_large_slice_x', **kwargs)
self.import_slice_data('default_large_slice_y', **kwargs)
self.import_slice_data('default_large_slice_z', **kwargs)
else:
print('Error!, Invalid interpolation mesh type!')
raise SystemExit
if '_s' in self.approach:
if interp_meshes == 'small':
self.import_slice_data('default_small_slice_x', **kwargs)
self.import_slice_data('default_small_slice_y', **kwargs)
self.import_slice_data('default_small_slice_z', **kwargs)
elif interp_meshes == 'large':
self.import_slice_data('default_large_slice_x', **kwargs)
self.import_slice_data('default_large_slice_y', **kwargs)
self.import_slice_data('default_large_slice_z', **kwargs)
###########################################################################
# # # line plot section # # #
###########################################################################
[docs]
def plot_point_data(self, mesh=None, EH='EH', mod=None, s_name=None,
ylim=None, grid=True, ap=False, sf=False, new=True,
label=None, legend=True, title=None, sharey=False,
key=None, stations=None, fs=12, save=True, **kwargs):
"""
work in priogress, see also *plot_line_data* docs
"""
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, points=True)
if stations is None:
stations = coords
xlim = [0 - 0.2, len(stations) - 0.8]
for k, key in enumerate(keys):
if new:
setattr(self, EH[k] + '_Lax', pu.make_subfigure_box(
var=EH[k], fs=fs, ylim=ylim, xlim=xlim, sf=sf,
grid=grid, cc=cc, ap=ap, sharey=sharey,
add_km=False, clr=self.label_color)[1])
if ap:
d1, d2 = self.amp_phase_comp(self.point_data[key])
print(d1.shape, d2.shape)
else:
d1 = np.abs(self.point_data[key].real)
d2 = np.abs(self.point_data[key].imag)
for j in range(3):
getattr(self, EH[k] + '_Lax')[j, 0].plot(
coords[:len(stations)], d1[:, j][stations], **kwargs)
tmp, = getattr(self, EH[k] + '_Lax')[j, 1].plot(
coords[:len(stations)], d2[:, j][stations], **kwargs)
if j == 2:
tmp.set_label(label)
if legend:
getattr(self, EH[k] + '_Lax')[2, 1].legend(loc='best')
self.save_plot(EH[k], title, save, s_name, fs, mesh)
[docs]
def plot_point_errors(self, mod=None, mod2=None, mesh=None, EH='EH',
err_type='rel', key=None, key2=None,
ylim=[-1e2, 1e2], new=True, sf=False, sf2=False,
log_scale=False, pyhed_ref=False, stations=None,
label=None, xlim=[-5., 5.], legend=True, fs=12,
title=None, save=True, s_name=None, **kwargs):
"""
Plot mismatch between two EM datasets, separated by
real and imaginary parts as well as the three vector components.
Keyword arguments:
------------------
- Most keyword arguments are identical to the ones documented in the
**plot_line_data()** function. Therefore, only differing arguments
are listed here.
- mod2 = None, type str
specify the 2nd model for mismatch calculation.
- key2 = None, type str
access the second model via the **key** if required. This
overwrites the **mod2** argument.
- err_type = 'rel', type str
specify which type of errors between both datasets should be used:
Either **'rel'** or **'abs'** or **'ratio'**.
- log_scale = False, type bool
set **True** for logarithmic error representation.
"""
if label is None:
label = err_type + ' [%]'
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, points=True)
keys2, label, EH, coords, cc = self.eval_properties(
mod2, key2, label, EH, sf2, mesh, points=True)
if stations is None:
stations = coords
xlim = [0 - 0.2, len(stations) - 0.8]
for k, key in enumerate(keys):
key2 = [kk for kk in keys2][k]
if pyhed_ref:
self.line_data[key2] = self.calc_pyhed_reference(
key2, EH=EH[k])
if err_type == 'rel':
diff, mag = self.rel_errors(self.point_data[key],
self.point_data[key2], store=False)
elif err_type == 'abs':
diff, mag = self.abs_errors(self.point_data[key],
self.point_data[key2], store=False)
elif err_type == 'ratio':
diff, mag = self.ratio(self.point_data[key],
self.point_data[key2], store=False)
else:
print('Error! "err_type" mus be "rel", "abs", or "ratio"!')
raise SystemExit
if new:
setattr(self, EH[k] + '_err_Lax', pu.make_subfigure_box(
sf=sf, height=4, var=EH[k], ylim=ylim, fs=fs,
xlim=xlim, log_scale=log_scale, err_plot='diff',
cc=cc, add_km=False, clr=self.label_color)[1])
if log_scale:
for j in range(3):
getattr(self, EH[k] + '_err_Lax')[j, 0].plot(
coords[:len(stations)], np.log(np.abs(
diff[:, j].real))[:len(stations)], **kwargs)
getattr(self, EH[k] + '_err_Lax')[j, 1].plot(
coords[:len(stations)], np.log(np.abs(
diff[:, j].imag))[:len(stations)], **kwargs)
getattr(self, EH[k] + '_err_Lax')[3, 0].plot(
coords[:len(stations)], np.log(
np.abs(mag.real))[:len(stations)], **kwargs)
tmp, = getattr(self, EH[k] + '_err_Lax')[3, 1].plot(
coords[:len(stations)], np.log(np.abs(
mag.imag))[:len(stations)], label=label, **kwargs)
else:
for j in range(3):
getattr(self, EH[k] + '_err_Lax')[j, 0].plot(
coords[:len(stations)], np.abs(
diff[:, j].real)[:len(stations)], **kwargs)
getattr(self, EH[k] + '_err_Lax')[j, 1].plot(
coords[:len(stations)], np.abs(
diff[:, j].imag)[:len(stations)], **kwargs)
getattr(self, EH[k] + '_err_Lax')[3, 0].plot(
coords[:len(stations)], np.abs(
mag.real)[:len(stations)], **kwargs)
tmp, = getattr(self, EH[k] + '_err_Lax')[3, 1].plot(
coords[:len(stations)], np.abs(
mag.imag)[:len(stations)], label=label,
**kwargs)
tmp.set_label(label)
if legend is not False:
getattr(self, EH[k] + '_err_Lax')[3, 1].legend(loc='best')
self.save_plot(EH[k], title, save, s_name, fs, mesh)
[docs]
def plot_line_data(self, mesh=None, EH='EH', mod=None, s_name=None,
km=1e3, ylim=None, xlim=[-5., 5.], grid=True,
sf=False, label=None, legend=True, title=None,
sharey=False, ap=False,
key=None, new=True, fs=12, save=True, **kwargs):
"""
Plot EM-data (if topo: projected) on 1D straight lines, separated by
real and imaginary parts as well as the three vector components.
The y-axis is always scaled logarithmically.
Keyword arguments:
------------------
- mesh = None, type str
specify if data are not based on the same coordinates as the ones
of the first dataset imported.
- EH = 'EH', type str
specify if either both (default) or only **E** or **H** fields
should be plotted.
- mod = None, type str
model name to access data from another model that **self.mod**
which was imported in the same **Plot** instance for comparison.
- s_name = None, type str
specify a custom name for saving plot(s) in **self.s_dir**.
default is **title + '_' + mesh + '.pdf'**
- km = 1e3, type float
Automatically convert coordinates from 'm' to km', set to **1.**
if 'm' should be kept or something else.
- ylim = None, type list with two arguments 'y_min' and 'y_max'
specify custom amplitude limits if required.
- xlim = [-5., 5.], type list with two arguments 'x_min' and 'x_max'
specify custom coordinate limits if required.
- grid = True, type bool
switch on grid lines or not.
- sf = False, type bool
set **True** If secondary fields and no total fields of the given
model should be visualized.
- label = None, type str
specify custom line label if required.
- legend = True, type bool,
switch on legend or not.
- title = None, type str
specify figure title if required. If **s_name** is **None**, title-
and mesh-name will be used to autmatically generate a save-name.
- key = None, type str
access a certain model via the **key** if required. This
overwrites the **mod** argument.
- new = True, type bool
set to False, if the current data should be added to the last
line plot, hence, no new figure will be created.
- fs = 12, type int
font size for legend and ticklabels.
- save = True, type bool
flag which controls if the figure is saved.
- **kwargs, type keyword arguments
additional keyword arguments for plots, e.g., 'lw' or 'color'.
"""
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh)
for k, key in enumerate(keys):
if new:
setattr(self, EH[k] + '_Lax', pu.make_subfigure_box(
var=EH[k], fs=fs, ylim=ylim, xlim=xlim, sf=sf,
grid=grid, cc=cc, ap=ap, sharey=sharey,
clr=self.label_color)[1])
if ap:
d1, d2 = self.amp_phase_comp(self.line_data[key])
else:
d1 = np.abs(self.line_data[key].real)
d2 = np.abs(self.line_data[key].imag)
for j in range(3):
getattr(self, EH[k] + '_Lax')[j, 0].plot(
coords/km, d1[:, j], **kwargs)
tmp, = getattr(self, EH[k] + '_Lax')[j, 1].plot(
coords/km, d2[:, j], **kwargs)
if j == 2:
tmp.set_label(label)
if legend:
getattr(self, EH[k] + '_Lax')[2, 1].legend(loc='best')
self.save_plot(EH[k], title, save, s_name, fs, mesh)
[docs]
def plot_line_errors(self, mod=None, mod2=None, mesh=None, EH='EH',
err_type='rel', key=None, key2=None, km=1e3,
ylim=[1e-2, 1e2], new=True, sf=False, sf2=None,
log_scale=True, pyhed_ref=False, magnitude=False,
label=None, xlim=[-5., 5.], legend=True, fs=12,
title=None, save=True, s_name=None, **kwargs):
"""
Plot mismatch between two EM datasets, separated by
real and imaginary parts as well as the three vector components.
Keyword arguments:
------------------
- Most keyword arguments are identical to the ones documented in the
**plot_line_data()** function. Therefore, only differing arguments
are listed here.
- mod2 = None, type str
specify the 2nd model for mismatch calculation.
- key2 = None, type str
access the second model via the **key** if required. This
overwrites the **mod2** argument.
- err_type = 'rel', type str
specify which type of errors between both datasets should be used:
Either **'rel'** or **'abs'** or **'ratio'**.
- log_scale = False, type bool
set **True** for logarithmic error representation.
"""
if sf2 is None:
sf2 = sf
if label is None:
label = err_type + ' [%]'
keys, label, EH_dummy, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh)
keys2, label, EH, coords, cc = self.eval_properties(
mod2, key2, label, EH, sf2, mesh)
for k, key in enumerate(keys):
key2 = [kk for kk in keys2][k]
if pyhed_ref:
self.line_data[key2] = self.calc_pyhed_reference(
key2, EH=EH[k])
if err_type == 'rel':
diff, mag = self.rel_errors(self.line_data[key],
self.line_data[key2], store=False)
elif err_type == 'abs':
diff, mag = self.abs_errors(self.line_data[key],
self.line_data[key2], store=False)
elif err_type == 'ratio':
diff, mag = self.ratio(self.line_data[key],
self.line_data[key2], store=False)
else:
print('Error! "err_type" mus be "rel", "abs", or "ratio"!')
raise SystemExit
if new:
if magnitude:
setattr(self, EH[k] + '_err_Lax', pu.make_subfigure_box(
sf=sf, height=4, var=EH[k], ylim=ylim, fs=fs,
xlim=xlim, log_scale=log_scale, err_plot='diff',
cc=cc, clr=self.label_color)[1])
else:
setattr(self, EH[k] + '_err_Lax', pu.make_subfigure_box(
sf=sf, height=3, var=EH[k], ylim=ylim, fs=fs,
xlim=xlim, log_scale=log_scale, err_plot='diff',
cc=cc, clr=self.label_color)[1])
if log_scale:
for j in range(2):
getattr(self, EH[k] + '_err_Lax')[j, 0].plot(
coords/km, np.abs(diff[:, j].real), **kwargs)
getattr(self, EH[k] + '_err_Lax')[j, 1].plot(
coords/km, np.abs(diff[:, j].imag), **kwargs)
getattr(self, EH[k] + '_err_Lax')[2, 0].plot(
coords/km, np.abs(diff[:, 2].real), **kwargs)
tmp, = getattr(self, EH[k] + '_err_Lax')[2, 1].plot(
coords/km, np.abs(diff[:, 2].imag), label=label, **kwargs)
else:
for j in range(2):
getattr(self, EH[k] + '_err_Lax')[j, 0].plot(
coords/km, np.abs(diff[:, j].real), **kwargs)
getattr(self, EH[k] + '_err_Lax')[j, 1].plot(
coords/km, np.abs(diff[:, j].imag), **kwargs)
getattr(self, EH[k] + '_err_Lax')[2, 0].plot(
coords/km, np.abs(mag.real), **kwargs)
tmp, = getattr(self, EH[k] + '_err_Lax')[2, 1].plot(
coords/km, np.abs(mag.imag), label=label,
**kwargs)
tmp.set_label(label)
if legend is not False:
getattr(self, EH[k] + '_err_Lax')[2, 1].legend(loc='best')
self.save_plot(EH[k], title, save, s_name, fs, mesh)
[docs]
def add_to_line_plot(self, ax, pos, errors=False, mod=None, mod2=None,
comp='x', log_scale=None, tcolor='#002C72', ap=False,
key=None, key2=None, sf=False, sf2=False, grid=True,
ri='real', EH='H', mesh=None, label=None, km=1e3,
llabel=None, fs=12, ylim=[0, 1e2], xlim=[-5., 5.],
err_type='rel', **kwargs):
"""
Add either line data or errors to a specific subplot.
Keyword arguments:
------------------
- All keyword arguments are identical to the ones documented in the
previous methods.
"""
if mesh is None:
mesh = [key for key in self.line_coords][0]
l_hacked_flag = True
if label is None:
l_hacked_flag = False
keys, label, EH_dummy, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh)
keys2, label, EH, coords, cc = self.eval_properties(
mod2, key2, label, EH, sf2, mesh)
for k, key in enumerate(keys):
key2 = [kk for kk in keys2][k]
if not errors and ri == 'real':
if comp in ['x', 'y', 'z']:
d = np.abs(self.line_data[key][:, self.c_idx[comp]].real)
elif comp == 'mag':
d = np.linalg.norm(self.line_data[key].real, axis=1)
elif not errors and ri == 'imag':
if comp in ['x', 'y', 'z']:
d = np.abs(self.line_data[key][:, self.c_idx[comp]].imag)
elif comp == 'mag':
d = np.linalg.norm(self.line_data[key].imag, axis=1)
elif not errors and ri == 'amp':
d, dummy = self.amp_phase_comp(self.line_data[key])
if comp in ['x', 'y', 'z']:
d = d[:, self.c_idx[comp]]
elif comp == 'mag':
d = np.linalg.norm(d, axis=1)
elif not errors and ri == 'phase':
dummy, d = self.amp_phase_comp(self.line_data[key])
if comp in ['x', 'y', 'z']:
d = d[:, self.c_idx[comp]]
elif comp == 'mag':
d = np.linalg.norm(d, axis=1)
ddr = np.linalg.norm(self.line_data[key].real, axis=1)
ddi = np.linalg.norm(self.line_data[key].imag, axis=1)
dd = np.zeros(len(ddr), dtype=complex)
dd.real = ddr
dd.imag = ddi
ddd = np.angle(dd, deg=True)
d = ddd
elif errors:
if not ap:
if err_type == 'rel':
diff, mag = self.rel_errors(self.line_data[key],
self.line_data[key2],
store=False)
elif err_type == 'abs':
diff, mag = self.abs_errors(self.line_data[key],
self.line_data[key2],
store=False)
elif err_type == 'ratio':
diff, mag = self.ratio(self.line_data[key],
self.line_data[key2],
store=False)
if ri == 'real':
if comp in ['x', 'y', 'z']:
d = np.abs(diff[:, self.c_idx[comp]].real)
elif comp == 'mag':
d = mag.real
elif ri == 'imag':
if comp in ['x', 'y', 'z']:
d = np.abs(diff[:, self.c_idx[comp]].imag)
elif comp == 'mag':
d = mag.imag
else:
if ri == 'amp':
d1, dummy = self.amp_phase_comp(self.line_data[key])
d2, dummy = self.amp_phase_comp(self.line_data[key2])
elif ri == 'phase':
dummy, d1 = self.amp_phase_comp(self.line_data[key])
dummy, d2 = self.amp_phase_comp(self.line_data[key2])
if err_type == 'ratio':
if ri == 'amp':
if comp in ['x', 'y', 'z']:
d = (d2[:, self.c_idx[comp]] /
d1[:, self.c_idx[comp]]) * 100.
elif comp == 'mag':
d = np.linalg.norm(d2, axis=1) / \
np.linalg.norm(d1, axis=1) * 100
elif ri == 'phase':
d = d1
else:
print('Not implmeneted yet')
raise SystemExit
if ri != 'phase':
data = np.abs(d)
else:
data = d
if llabel is None:
self.tmp_plot, = ax[pos[0], pos[1]].plot(coords/km, data,
**kwargs)
else:
self.tmp_plot, = ax[pos[0], pos[1]].plot(
coords/km, data, label=llabel, **kwargs)
if l_hacked_flag:
ax[pos[0], pos[1]].set_title(label, fontsize=fs, color=tcolor)
ax[pos[0], pos[1]].set_ylim(ylim)
ax[pos[0], pos[1]].set_xlim(xlim)
if grid:
ax[pos[0], pos[1]].grid(which='major', color='0.8', ls=':')
if log_scale:
ax[pos[0], pos[1]].set_yscale('log')
elif log_scale is False:
ax[pos[0], pos[1]].set_yscale('linear')
###########################################################################
# # # slice plot section # # #
###########################################################################
[docs]
def plot_slice_data(self, mesh=None, mod=None, sf=False, km=1., ap=False,
n_colors=21, EH='EH', label=None, fs=12, c_lim=None,
cmap=None, equal_axis=True, xlim=None, ylim=None,
key=None, title=None, save=True, s_name=None,
shift=0., var=None, dummy=True, kk=0):
"""
Plot 2D EM-data (if topo: projected) on straight slices, separated by
real and imaginary parts as well as the three vector components.
The y-axis is always scaled logarithmically. Amplitudes are represented
by the color map.
Keyword arguments:
------------------
Most keyword arguments are identical to the ones documented in the
**plot_line_data()** and **plot_line_errors()** functions.
Therefore, only differing arguments are listed here.
- n_colors = 101, type int
number of different colors for the colormap
- cmap = None, type str
specify a custom matplotlib colormap via string-name. Currently
supported: **'magma'**, **'viridis'** and **'RdBu_r'**
by default: plt.cm.magma is used for *E* fields,
plt.cm.virids is used for *H* fields,
plt.cm.RdBu_r is used for error plots
- c_lim = None, type list with two arguments 'c_min' and 'c_max'
specify custom amplitude limits if required
- equal_axis = True, type bool
set to **False** for non-equal x- and y-axis aspect ratio
"""
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, line=False)
coords0 = coords[0] * km
coords1 = coords[1] * km
for k, key in enumerate(keys):
if var is None:
var = [EH[k], EH[k], EH[k]]
dummy = False
kk = k
if k == 1 and dummy:
pass
else:
f, ax = pu.make_subfigure_box(
var=var, log_scale=False, cc=cc,
fs=fs, sf=sf, sliceplot=True,
clr=self.label_color, ap=ap)
setattr(self, EH[k] + '_Sax', ax)
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
self.slice_data[key], n_colors, c_lim)
if ap:
phase_colors = np.load(os.path.dirname(custEM.misc.__file__) +
'\hpluv_rgb_colors.npy')
self.phase_map = clrs.ListedColormap(phase_colors)
d1, d2 = self.amp_phase_comp(self.slice_data[key], shift=shift)
else:
d1 = np.abs(self.slice_data[key].real)
d2 = np.abs(self.slice_data[key].imag)
if EH[k] == 'E':
ll = r'$\mathbf{E}$ $(V/m)$'
else:
ll = r'$\mathbf{H}$ $(A/m)$'
for j in range(3):
if var[j] != EH[k]:
continue
if dummy and k == 1:
if var == 'E':
continue
cm = self.init_cmap(cmap, var[j])
d = [d1[:, j].reshape(len(coords0), len(coords1))]
d.append(d2[:, j].reshape(len(coords0), len(coords1)))
if ap:
self.tmp_plot = getattr(
self, EH[kk] + '_Sax')[j, 0].contourf(
coords0, coords1, d[0].T, levels=colors,
norm=clrs.LogNorm(), cmap=cm)
self.tmp_plot.set_clim(c1, c2)
if j == 2:
self.add_cbar(f, [c1, c2], fs=fs, label=ll,
pos=[0.8, 0.12, 0.03, 0.75],
additional=1, cmap2=dummy)
self.tmp_plot = getattr(
self, EH[kk] + '_Sax')[j, 1].pcolormesh(
coords0, coords1, d[1].T, cmap=self.phase_map)
self.tmp_plot.set_clim(-180, 180)
if j == 2:
self.add_phase_bar(f, fs=fs, dummy=dummy)
else:
for jj in range(2):
self.tmp_plot = getattr(
self, EH[kk] + '_Sax')[j, jj].contourf(
coords0, coords1, d[jj].T, levels=colors,
norm=clrs.LogNorm(), cmap=cm)
self.tmp_plot.set_clim(c1, c2)
if j == 2:
self.add_cbar(f, [c1, c2], fs=fs, label=ll,
additional=1, cmap2=dummy)
if var == ['E', 'E', 'E']:
var = None
pu.adjust_axes(getattr(self, EH[kk] + '_Sax'), equal_axis,
xlim, ylim)
pu.adjust_ticks_and_labels(ax)
self.save_plot(EH[kk], title, save, s_name, fs, mesh, png=True)
[docs]
def plot_slice_errors(self, mod=None, mod2=None, mesh=None, EH='EH',
err_type='rel', key=None, key2=None, ap=False,
n_colors=101, c_lim=[1e-1, 1e2], p_lim=1.,
sf=False, fs=12, pyhed_ref=False, mask=None,
label='err [%]', cmap=None, equal_axis=True,
xlim=None, ylim=None, magnitude=False, shift=0.,
title=None, save=True, s_name=None):
"""
Plot mismatches between two 2D EM-datasets, separated by
real and imaginary parts as well as the three vector components.
The y-axis is always scaled logarithmically. Amplitudes are represented
by the color map.
Keyword arguments:
------------------
- All keyword arguments are identical to the ones documented in the
**plot_line_data()**, **plot_line_errors()** and
**plot_slice_data()** functions. It is referred to the latter.
"""
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, line=False)
keys2, label, EH, coords, cc = self.eval_properties(
mod2, key2, label, EH, sf, mesh, line=False)
for k, key in enumerate(keys):
key2 = [kk for kk in keys2][k]
if pyhed_ref:
self.slice_data[key2] = self.calc_pyhed_reference(
key2, EH=EH[k])
if ap:
d11, d12 = self.amp_phase_comp(self.slice_data[key],
shift=shift)
d21, d22 = self.amp_phase_comp(self.slice_data[key2],
shift=shift)
if err_type == 'rel':
diff, mag = self.rel_errors(d11, d21, store=False)
diff2, mag2 = self.abs_errors(d12, d22, store=False)
diff.imag = diff2.real
mag.imag = mag2.real
else:
print('Error! "err_type" must be "rel" for "ap=True" !')
raise SystemExit
else:
d1 = self.slice_data[key]
d2 = self.slice_data[key2]
if err_type == 'rel':
diff, mag = self.rel_errors(d1, d2, store=False)
elif err_type == 'abs':
diff, mag = self.abs_errors(d1, d2, store=False)
elif err_type == 'ratio':
diff, mag = self.ratio(d1, d2, store=False)
else:
print('Error! "err_type" mus be "rel", "abs", or "ratio"!')
raise SystemExit
cm = self.init_cmap(cmap=cmap, var='err')
if magnitude:
h = 4
else:
h = 3
f, ax = pu.make_subfigure_box(
height=h, var=EH[k], log_scale=False, cc=cc, ap=ap,
fs=fs, err_plot='diff', sharey=True, sf=sf, sliceplot=True,
clr=self.label_color)
setattr(self, EH[k] + '_err_Sax', ax)
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
np.hstack((diff, mag.reshape(-1, 1))),
n_colors, c_lim, quant=90, symlog=True)
for j in range(h):
if j != 3:
d = [diff[:, j].real.reshape(len(coords[0]),
len(coords[1]))]
d.append(diff[:, j].imag.reshape(len(coords[0]),
len(coords[1])))
if mask is not None:
d[0][mask[0][: , j].reshape(d[0].shape)==False] = 0.
d[1][mask[1][: , j].reshape(d[1].shape)==False] = 0.
else:
d = [mag.real.reshape(len(coords[0]), len(coords[1]))]
d.append(mag.imag.reshape(len(coords[0]), len(coords[1])))
for jj in range(2):
if ap and jj == 1:
d[jj][d[jj] < -p_lim] = -p_lim
d[jj][d[jj] > p_lim] = p_lim
self.tmp_plot = getattr(
self, EH[k] + '_err_Sax')[j, jj].pcolormesh(
coords[0], coords[1], d[jj].T,
vmin=-p_lim, vmax=p_lim, cmap='PuOr')
else:
self.tmp_plot = getattr(
self, EH[k] + '_err_Sax')[j, jj].pcolormesh(
coords[0], coords[1], d[jj].T,
norm=clrs.SymLogNorm(linthresh=c1,
vmin=-10**maxlog,
vmax=10**maxlog), cmap=cm)
if ap and j == 2:
if jj == 0:
self.add_cbar(f, [c1, c2], fs=fs, additional=1,
cmap='RdBu_r', symlog=True,
pos=[0.8, 0.12, 0.03, 0.75],
label=r'$\mathbf{\epsilon}$ $(\%)$')
elif jj == 1:
self.add_phase_bar(f, fs=fs, delta=True,
p_lim=p_lim)
else:
if j == 2:
self.add_cbar(f, [c1, c2], fs=fs, additional=1,
cmap='RdBu_r', symlog=True,
label=r'$\mathbf{\epsilon}$ $(\%)$')
pu.adjust_axes(getattr(self, EH[k] + '_err_Sax'), equal_axis,
xlim, ylim, size=h)
pu.adjust_ticks_and_labels(ax)
if title is not None:
plt.suptitle(title, fontsize=fs)
self.save_plot(EH[k], title, save, s_name, fs, mesh, png=True)
[docs]
def plot_magnitude_angle_misfits(self, mod=None, mod2=None,
mesh=None, EH='EH', cmap_width=0.9,
err_type='rel', key=None, key2=None,
n_colors=101, e_lim=[1e-1, 1e2],
p_lim=10., sf=False, fs=12,
a_lim=[1e-12, 1e-0],
label='err [%]', cmap=None, shift=0.,
equal_axis=True, xlim=None, ylim=None,
title=None, save=True, s_name=None):
"""
Description!
Keyword arguments:
------------------
"""
keys, label, EH, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, line=False)
keys2, label, EH, coords, cc = self.eval_properties(
mod2, key2, label, EH, sf, mesh, line=False)
for k, key in enumerate(keys):
key2 = [kk for kk in keys2][k]
if err_type == 'rel':
diff, mag = self.rel_errors(self.slice_data[key],
self.slice_data[key2], store=False)
da, mag_a = self.abs_errors(self.slice_data[key],
self.slice_data[key2], store=False)
angles = self.abs_angles(self.slice_data[key],
self.slice_data[key2], store=False)
else:
print('Error! "err_type" must be "rel" for "ap=True" !')
raise SystemExit
cm = self.init_cmap(cmap=cmap, var='err')
cm2 = self.init_cmap(cmap=cmap, var=EH[k])
if EH[k] == 'E':
ll = r'$\mathbf{E}$ $(V/m)$'
else:
ll = r'$\mathbf{H}$ $(A/m)$'
f, ax = pu.make_subfigure_box(
height=3, width=3, var=EH[k], log_scale=False, cc=cc,
ap=True, fs=fs, err_plot='diff', sharey=True, sf=sf,
sliceplot=True, clr=self.label_color)
setattr(self, EH[k] + '_mag_amp_ax', ax)
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
np.hstack((diff, mag.reshape(-1, 1))),
n_colors, e_lim, quant=90, symlog=True)
minlog2, maxlog2, colors2, c12, c22 = pu.eval_colors(
np.hstack((da, mag_a.reshape(-1, 1))),
n_colors, a_lim, quant=90, symlog=False)
mag_a2 = np.abs(mag_a)
v1, _ = self.amp_phase_mag(self.slice_data[key])
v2, _ = self.amp_phase_mag(self.slice_data[key2])
rel_e = ((v2 - v1) / v1) * 100.
cos_a_real = np.sum(np.abs(self.slice_data[key] *
self.slice_data[key2]), axis=1)
scale_real = (np.linalg.norm(np.abs(
self.slice_data[key]), axis=1) *
np.linalg.norm(np.abs(
self.slice_data[key2]), axis=1))
angles2 = np.rad2deg(np.arccos(cos_a_real / scale_real))
pu.adjust_axes(getattr(self, EH[k] + '_mag_amp_ax'), equal_axis,
xlim, ylim, size=3, width=3)
p0 = getattr(
self, EH[k] + '_mag_amp_ax')[0, 2].get_position(
).get_points().flatten()
p1 = getattr(
self, EH[k] + '_mag_amp_ax')[1, 2].get_position(
).get_points().flatten()
p2 = getattr(
self, EH[k] + '_mag_amp_ax')[2, 2].get_position(
).get_points().flatten()
# abs misfits
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[0, 0].contourf(
coords[0], coords[1], mag_a.reshape(
len(coords[0]), len(coords[1])).real.T,
levels=colors2, norm=clrs.LogNorm(), cmap=cm2,
vmin=a_lim[0], vmax=a_lim[1])
for c in self.tmp_plot.collections:
c.set_edgecolor("face")
getattr(self, EH[k] + '_mag_amp_ax')[0, 0].set_title(
r'$\mathbf{\Delta}_{||\Re(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[0, 1].contourf(
coords[0], coords[1], mag_a.reshape(
len(coords[0]), len(coords[1])).imag.T,
levels=colors2, norm=clrs.LogNorm(), cmap=cm2,
vmin=a_lim[0], vmax=a_lim[1])
for c in self.tmp_plot.collections:
c.set_edgecolor("face")
getattr(self, EH[k] + '_mag_amp_ax')[0, 1].set_title(
r'$\mathbf{\Delta}_{||\Im(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[0, 2].contourf(
coords[0], coords[1], mag_a2.reshape(
len(coords[0]), len(coords[1])).T,
levels=colors2, norm=clrs.LogNorm(), cmap=cm2,
vmin=a_lim[0], vmax=a_lim[1])
for c in self.tmp_plot.collections:
c.set_edgecolor("face")
getattr(self, EH[k] + '_mag_amp_ax')[0, 2].set_title(
r'$\mathbf{\Delta}_{||\mathbf{' + EH[k] + '}||}$')
self.add_cbar(f, [c12, c22], fs=fs, additional=1,
pos=[0.9, p0[1], 0.02, (p0[1]-p1[1]) * cmap_width],
label=ll,
diff=0.05, label_pos=[0.1, p0[1] + 0.4])
# relative misfits
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[1, 0].pcolormesh(
coords[0], coords[1], mag.reshape(
len(coords[0]), len(coords[1])).real.T,
norm=clrs.SymLogNorm(linthresh=c1,
vmin=-10**maxlog,
vmax=10**maxlog),
cmap=cm, edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[1, 0].set_title(
r'$\mathbf{\epsilon}_{||\Re(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[1, 1].pcolormesh(
coords[0], coords[1], mag.reshape(
len(coords[0]), len(coords[1])).imag.T,
norm=clrs.SymLogNorm(linthresh=c1,
vmin=-10**maxlog,
vmax=10**maxlog),
cmap=cm, edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[1, 1].set_title(
r'$\mathbf{\epsilon}_{||\Im(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[1, 2].pcolormesh(
coords[0], coords[1], rel_e.reshape(
len(coords[0]), len(coords[1])).T,
norm=clrs.SymLogNorm(linthresh=c1,
vmin=-10**maxlog,
vmax=10**maxlog),
cmap=cm, edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[1, 2].set_title(
r'$\mathbf{\epsilon}_{||\mathbf{' + EH[k] + '}||}$')
self.add_cbar(f, [c1, c2], fs=fs, additional=1, symlog=True,
pos=[0.9, p1[1], 0.02, (p0[1]-p1[1]) * cmap_width],
label=r'$\mathbf{\epsilon}$ $(\%)$', diff=0.05,
label_pos=[1.4, p1[1] + 0.5])
# angles between magnitude vectors
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[2, 0].pcolormesh(
coords[0], coords[1], angles.reshape(
len(coords[0]), len(coords[1])).real.T,
vmin=0., vmax=p_lim, cmap='Purples', edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[2, 0].set_title(
r'$\mathbf{\alpha}_{||\Re(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[2, 1].pcolormesh(
coords[0], coords[1], angles.reshape(
len(coords[0]), len(coords[1])).imag.T,
vmin=0., vmax=p_lim, cmap='Purples', edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[2, 1].set_title(
r'$\mathbf{\alpha}_{||\Im(\mathbf{' + EH[k] + '})||}$')
self.tmp_plot = getattr(
self, EH[k] + '_mag_amp_ax')[2, 2].pcolormesh(
coords[0], coords[1], angles2.reshape(
len(coords[0]), len(coords[1])).T,
vmin=0., vmax=p_lim, cmap='Purples', edgecolors='face')
getattr(self, EH[k] + '_mag_amp_ax')[2, 2].set_title(
r'$\mathbf{\alpha}_{||\mathbf{' + EH[k] + '}||}$')
self.add_phase_bar(
f, fs=fs, absdelta=True, p_lim=p_lim,
pos=[0.9, p2[1], 0.02, (p0[1]-p1[1]) * cmap_width],
label=r'$\mathbf{\alpha}$ (°)',
label_pos=[1.4, p2[1] + 0.7])
pu.adjust_ticks_and_labels(ax, remove_top_right=True)
if title is not None:
plt.suptitle(title, fontsize=fs)
self.save_plot(EH[k], title, save, s_name, fs, mesh, png=True)
[docs]
def add_to_slice_plot(self, ax, pos, comp='mag', ri='real', n_colors=101,
mod=None, mesh=None, key=None, EH='E',
label=None, fs=12, c_lim=None, sf=False,
tcolor='#002C72', horz_only=False,
cmap=None, equal_axis=True, quiver=True,
title=None, save=True, s_name=None):
"""
Add a selected dataset, i.e., a specific component, real or imag part,
to a subplot specified by a given figure axes object and position.
In development (quiver option needs to get variable and no hard coded
keyword arguments!)
Required arguments:
-------------------
ax, type matplotlib figure axes object
axes object of a matplotlib figure which might be separated into
numerous subplots.
pos, type list of two entried
specify the subplot position of the **ax** object to plot in.
Keyword arguments:
------------------
- Most keyword arguments are identical to the ones documented in the
**plot_line_data()**, **plot_line_errors()** and
**plot_slice_data()** functions. Therefore, only new keyword
arguments are listed below:
- comp = 'mag', type str
specify if either the **'x'**, **'y'**, **'z'** component or the
magnitude (**'mag'**) should be visualized.
- ri = 'real', type str
plot either **'real'** or **'imag'** parts of the data.
- quiver = True, type bool
plot vector, this option is not optimized yet and work in progress.
"""
cmap = self.init_cmap(cmap=cmap, var=EH)
keys, dummy, dummy2, coords, cc = self.eval_properties(
mod, key, label, EH, sf, mesh, line=False)
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
self.slice_data[keys[0]], n_colors, c_lim)
if ri == 'real':
label_ri = r'$\Re$'
if comp == 'mag':
d = (np.linalg.norm(self.slice_data[keys[0]].real, axis=1).
reshape(len(coords[0]), len(coords[1])))
if horz_only:
d = np.linalg.norm(self.slice_data[keys[0]][:, :2].real,
axis=1).reshape(len(coords[0]),
len(coords[1]))
else:
d = np.abs(self.slice_data[keys[0]][
:, self.c_idx[comp]].real).reshape(
len(coords[0]), len(coords[1]))
elif ri == 'imag':
label_ri = r'$\Im$'
if comp == 'mag':
d = (np.linalg.norm(self.slice_data[keys[0]].imag, axis=1).
reshape(len(coords[0]), len(coords[1])))
if horz_only:
d = np.linalg.norm(self.slice_data[keys[0]][:, :2].imag,
axis=1).reshape(len(coords[0]),
len(coords[1]))
else:
d = np.abs(self.slice_data[keys[0]][
:, self.c_idx[comp]].imag).reshape(
len(coords[0]), len(coords[1]))
else:
('Warning! Choose either real or imag for "ri" variable!')
raise SystemExit
if label is None:
label = label_ri + '($\mathbf{' + EH + '}_' + comp + '$) FEM'
self.tmp_plot = ax[pos[0], pos[1]].contourf(
coords[0], coords[1], d.T, levels=colors,
norm=clrs.LogNorm(), cmap=cmap)
for c in self.tmp_plot.collections:
c.set_edgecolor("face")
ax[pos[0], pos[1]].set_title(label, fontsize=fs, color=tcolor)
self.tmp_plot.set_clim(10.**minlog, 10.**maxlog)
if equal_axis:
ax[pos[0], pos[1]].axis('equal')
try:
ax[pos[0], pos[1]].set_adjustable('box')
except exception as e:
ax[pos[0], pos[1]].set_adjustable('box-forced')
[docs]
def add_errors_to_slice_plot(self, ax, pos, err_key=None, key=None,
key2=None, comp='mag', mask=None,
ri='real', EH='E', mesh=None, n_colors=101,
label=None, fs=12, e_lim=[1e-1, 1e2],
tcolor='#002C72',
cmap=None, equal_axis=True, err_type='rel'):
"""
Add a selected dataset, i.e., a specific component, real or imag part,
to a subplot specified by a given figure axes object and position.
In development (quiver option needs to get variable and no hard coded
keyword arguments!)
Required arguments:
-------------------
ax, type matplotlib figure axes object
axes object of a matplotlib figure which might be separated into
numerous subplots.
pos, type list of two entried
specify the subplot position of the **ax** object to plot in.
err_key, type str
specify error data to plot - valid key for one of the dictionaries:
**self.rel_mag_errors**, **self.rel_comp_errors**,
**self.abs_mag_errors**, **self.abs_comp_errors**,
**self.ratio_mag_errors**, or **self.ratio_comp_errors**.
Keyword arguments:
------------------
- All keyword arguments are identical to the ones documented in the
**plot_line_data()**, **plot_line_errors()**, **plot_slice_data()**
and **add_to_slice_plot()** functions. Therefore, only new keyword
arguments are listed below:
"""
if mesh is None:
mesh = [key for key in self.slice_coords][0]
coords0, coords1 = self.get_coords(mesh)
coords = self.get_coords(mesh)
# account for pcolormesh coordinate shift
coords0 -= np.diff(coords0[:2])[0] / 2.
coords0 = np.append(coords0[:],
coords0[-1] + np.diff(coords0[:2]))
coords1 -= np.diff(coords1[:2])[0] / 2.
coords1 = np.append(coords1,
coords1[-1] + np.diff(coords1[:2]))
a, b = len(coords[0]), len(coords[1])
if err_type != 'rel' and err_type != 'abs' and err_type != 'ratio':
print('Error! "err_type" mus be "rel", "abs", or "ratio"!')
raise SystemExit
if err_type != 'abs':
var = 'err'
else:
var = EH
if comp == 'x':
c = 0
elif comp == 'y':
c = 1
elif comp == 'z':
c = 2
cm = self.init_cmap(cmap=cmap, var=var)
if ri == 'real':
if comp == 'mag':
if err_type == 'rel':
d = self.rel_mag_errors[err_key].real.reshape(a, b)
elif err_type == 'abs':
d = np.abs(self.abs_mag_errors[
err_key].real.reshape(a, b))
elif err_type == 'ratio':
d = self.ratio_mag_errors[err_key].real.reshape(a, b)
else:
if err_type == 'rel':
d = self.rel_comp_errors[
err_key][:, self.c_idx[comp]].real.reshape(a, b)
elif err_type == 'abs':
d = np.abs(self.abs_comp_errors[
err_key][:, self.c_idx[comp]].real.reshape(a, b))
elif err_type == 'ratio':
d = self.ratio_comp_errors[
err_key][:, self.c_idx[comp]].real.reshape(a, b)
elif ri == 'imag':
if comp == 'mag':
if err_type == 'rel':
d = self.rel_mag_errors[err_key].imag.reshape(a, b)
elif err_type == 'abs':
d = np.abs(self.abs_mag_errors[
err_key].imag.reshape(a, b))
elif err_type == 'ratio':
d = self.ratio_mag_errors[err_key].imag.reshape(a, b)
else:
if err_type == 'rel':
d = self.rel_comp_errors[
err_key][:, self.c_idx[comp]].imag.reshape(a, b)
elif err_type == 'abs':
d = np.abs(self.abs_comp_errors[
err_key][:, self.c_idx[comp]].imag.reshape(a, b))
elif err_type == 'ratio':
d = self.ratio_comp_errors[
err_key][:, self.c_idx[comp]].imag.reshape(a, b)
elif ri == 'abs':
if comp == 'mag':
v1, _ = self.amp_phase_mag(self.slice_data[
key + '_' + EH + '_t'])
v2, _ = self.amp_phase_mag(self.slice_data[
key2 + '_' + EH + '_t'])
d = (((v2 - v1) / v1) * 100.).reshape(len(coords[0]),
len(coords[1]))
else:
print('Error! Component-wise relative misfits of absolute '
'parts not implemented yet. Abprting ...')
elif ri == 'angles':
cos_a_real = np.sum(np.abs(
self.slice_data[key + '_' + EH + '_t'] *
self.slice_data[key2 + '_' + EH + '_t']), axis=1)
scale_real = (np.linalg.norm(np.abs(self.slice_data[
key + '_' + EH + '_t']), axis=1) *
np.linalg.norm(np.abs(self.slice_data[
key2 + '_' + EH + '_t']), axis=1))
d = np.rad2deg(np.arccos(cos_a_real / scale_real)).reshape(
len(coords[0]), len(coords[1]))
else:
('Typing error! "ri" variable must be "real", "imag", "abs" or '
'"angles"!')
raise SystemExit
if mask is not None:
d[mask[:, c].reshape(d.shape)==False] = 0.
if err_type != 'abs':
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
d, n_colors, e_lim, quant=90, symlog=True)
norm = clrs.SymLogNorm(linthresh=c1, vmin=-10**maxlog,
vmax=10**maxlog)
else:
minlog, maxlog, colors, c1, c2 = pu.eval_colors(
d, n_colors, e_lim, quant=90, symlog=False)
norm = clrs.LogNorm()
if ri != 'angles':
# self.tmp_plot = ax[pos[0], pos[1]].contourf(
# coords[0], coords[1], d.T, levels=colors,
# norm=norm, cmap=cm)
# for c in self.tmp_plot.collections:
# c.set_edgecolor("face")
self.tmp_plot = ax[pos[0], pos[1]].pcolormesh(
coords0, coords1, d.T,
edgecolors='face',
norm=clrs.SymLogNorm(linthresh=c1,
vmin=-10**maxlog,
vmax=10**maxlog), cmap=cm)
else:
self.tmp_plot = ax[pos[0], pos[1]].pcolormesh(
coords0, coords1, d.T, vmin=e_lim[0],
vmax=e_lim[1], cmap='Purples', edgecolors='face')
if err_type == 'abs':
self.tmp_plot.set_clim(10.**minlog, 10.**maxlog)
ax[pos[0], pos[1]].set_title(label, fontsize=fs, color=tcolor)
if equal_axis:
ax[pos[0], pos[1]].axis('equal')
try:
ax[pos[0], pos[1]].set_adjustable('box')
except exception as e:
ax[pos[0], pos[1]].set_adjustable('box-forced')
[docs]
def add_3D_slice_plot(self, fig, full_name, slicE, sp_pos=111, fs=12,
q_slicE=None, q_name=None, ri='real', cbar=True,
EH='E', mesh=None, n_colors=101, q_length=0.2,
label=None, c_lim=None, cmap='magma',
equal_axis=True,
err_type='rel'):
"""
Plot 2D EM-data with real 3D geometry (topography) at specified
positions in a given figure. In development (there are still
some hard coded keyword arguments which need to be made variable!)
Required arguments:
-------------------
- fig, type matplotlib figure object
figure to plot in
- full_name, type str
full name: a valid key for the dictionary **self.slice_data**,
(characterization via **mod** or **key** args. in progress).
- sliceE, type str
slice name: a valid key for the dictionary **self.slice_coords**
Keyword arguments:
------------------
- Most keyword arguments are identical to the ones documented in the
**plot_line_data()**, **plot_line_errors()**, **plot_slice_data()**
and **add_to_slice_plot()** functions. Therefore, only new keyword
arguments are listed below:
- q_name = None, type str
valid key for the dictionary **self.slice_data**, used to select
a dataset for plotting vector arrows with a lower grid/data density
than the model specified via **full_name** for the contour plot.
- q_sliceE, type str
valid key for the dictionary **self.slice_coord**, used to select
a dataset for plotting vector arrows with a lower grid/data density
than the model specified via **slicE** for the contour plot.
- q_length = 0.2, type float
specify a custom length for all vector arrows.
"""
from mpl_toolkits.mplot3d import axes3d
from scipy.interpolate import griddata
axes3d.Axes3D
ax = fig.add_subplot(sp_pos, projection='3d')
X = self.slice_coords[slicE].real[:, 0]/1e3
Y = self.slice_coords[slicE].real[:, 1]/1e3
Z = self.slice_coords[slicE].real[:, 2]/1e3
size = int(np.sqrt(len(X)))
if ri == 'real':
Zvals = np.linalg.norm(self.slice_data[full_name].real, axis=1
).reshape(size, size)
elif ri == 'imag':
Zvals = np.linalg.norm(self.slice_data[full_name].imag, axis=1
).reshape(size, size)
xi = np.linspace(X.min(), X.max(), size)
yi = np.linspace(Y.min(), Y.max(), size)
zi = griddata((X, Y), Z, (xi[None, :], yi[:, None]), method='cubic')
if c_lim is None:
origmin = np.min(Zvals)
origmax = np.max(Zvals)
else:
origmin = c_lim[0]
origmax = c_lim[1]
minlog = np.log(origmin * (1. / origmin))
maxlog = np.log(origmax * (1. / origmin))
colval = np.zeros((len(Zvals), len(Zvals)))
plotcol = np.zeros((len(Zvals), len(Zvals), 4))
cmap = self.init_cmap(cmap)
for k in range(len(Zvals)):
for j in range(len(Zvals)):
colval[k, j] = (1 * np.log(Zvals[k, j] * (1 / origmin)) -
minlog) / (maxlog - minlog)
plotcol[k, j, :] = cmap(colval[k, j])
xig, yig = np.meshgrid(xi, yi)
ax.plot_surface(xig, yig, zi, rstride=1, cstride=1, facecolors=plotcol,
cmap=cmap, norm=LogNorm(), linewidth=0,
antialiased=False, alpha=1.)
max_range = np.array([X.max()-X.min(), Y.max()-Y.min(),
Z.max()-Z.min()]).max() / 2.0
mid_x = (X.max()+X.min()) * 0.5
mid_y = (Y.max()+Y.min()) * 0.5
mid_z = (Z.max()+Z.min()) * 0.5
ax.set_xlim(mid_x - max_range, mid_x + max_range)
ax.set_ylim(mid_y - max_range, mid_y + max_range)
ax.set_zlim(mid_z - max_range, mid_z + max_range)
ax.set_xlabel('x (km)', fontsize=fs)
ax.set_ylabel('y (km)', fontsize=fs)
ax.set_zlabel('elevation (km)', fontsize=fs)
if q_slicE is not None:
X = self.slice_coords[q_slicE].real[:, 0]/1e3
Y = self.slice_coords[q_slicE].real[:, 1]/1e3
Z = self.slice_coords[q_slicE].real[:, 2]/1e3
if ri == 'real':
U = self.slice_data[q_name][:, 0].real
V = self.slice_data[q_name][:, 1].real
W = self.slice_data[q_name][:, 2].real
elif ri == 'imag':
U = self.slice_data[q_name][:, 0].imag
V = self.slice_data[q_name][:, 1].imag
W = self.slice_data[q_name][:, 2].imag
ax.quiver(X, Y, Z, U, V, W, pivot='tail', color='k',
normalize=True, length=q_length)
ax.view_init(elev=40.0, azim=125.)
ax.set_xlim3d((mid_x - max_range)*0.95, (mid_x + max_range)*0.95)
ax.set_ylim3d((mid_x - max_range)*0.95, (mid_x + max_range)*0.95)
ax.set_zlim3d((mid_x - max_range)*0.95, (mid_x + max_range)*0.95)
if cbar:
m = cm.ScalarMappable(cmap=cmap,
norm=LogNorm(vmin=c_lim[0], vmax=c_lim[1]))
m.set_array(np.array(c_lim))
fig.subplots_adjust(right=0.88)
cbar_ax = fig.add_axes([0.9, 0.12, 0.03, 0.75])
cb = fig.colorbar(m, cax=cbar_ax, label=label)
cb.ax.tick_params(labelsize=fs)
pu.adjust_log_ticks(cb, np.log10(origmin), np.log10(origmax),
additional=1)
[docs]
def add_cbar(self, fig, c_lim, cmap='magma', cmap2=False,
pos=[0.9, 0.12, 0.03, 0.75], horizontal=False,
pos2=[0.05, 0.12, 0.03, 0.75],
diff=0.05, n_colors=101, fs=12, symlog=False, label=None,
label_pos=[0., 1.02], additional=0, where='left'):
"""
Add colorbar to a 2D / 3D data plot at a custom position with custom
style options.
"""
cmap = self.init_cmap(cmap)
if cmap == 'RdBu_r':
symlog = True
minlog, maxlog, colors, c1, c2 = pu.eval_colors(None, n_colors, c_lim)
if where == 'left' or where == 'both':
fig.subplots_adjust(right=pos[0] - diff)
if cmap2:
pos[1] -= 0.02
pos[2] *= 0.5
cbar_ax = fig.add_axes(pos)
pos[0] += pos[2]
cbar_ax2 = fig.add_axes(pos)
else:
cbar_ax = fig.add_axes(pos)
elif where == 'right' or where == 'both':
fig.subplots_adjust(left=pos2[0] + diff)
cbar_ax = fig.add_axes(pos2)
elif where == 'top':
fig.subplots_adjust(top=pos[1] - diff)
cbar_ax = fig.add_axes(pos)
elif where == 'bottom':
fig.subplots_adjust(bottom=pos[1] + diff)
cbar_ax = fig.add_axes(pos)
if cmap2:
tmp_f, bx = plt.subplots(1)
bla = bx.contourf(np.linspace(0., 1., 3), np.linspace(0., 1., 3),
np.linspace(c1, c2, 9).reshape(3, 3),
levels=np.logspace(minlog, maxlog, n_colors),
cmap='viridis', norm=clrs.LogNorm())
bla.set_clim(c1, c2)
cb = tmp_f.colorbar(bla, cax=cbar_ax2)
plt.close(tmp_f)
tmp_f, bx = plt.subplots(1)
bla = bx.contourf(np.linspace(0., 1., 3), np.linspace(0., 1., 3),
np.linspace(c1, c2, 9).reshape(3, 3),
levels=np.logspace(minlog, maxlog, n_colors),
cmap='magma', norm=clrs.LogNorm())
bla.set_clim(c1, c2)
cb2 = tmp_f.colorbar(bla, cax=cbar_ax)
cb2.ax.set_yticklabels([])
plt.close(tmp_f)
ll = r'$\mathbf{E}$ $(V/m)$ /' + '\n' + r'$\mathbf{H}$ $(A/m)$'
if label is not None:
cbar_ax.text(label_pos[0], label_pos[1], ll, fontsize=fs,
color=self.label_color)
else:
if label is not None:
cbar_ax.text(label_pos[0], label_pos[1], label, fontsize=fs,
color=self.label_color)
if horizontal:
cb = fig.colorbar(self.tmp_plot, cax=cbar_ax, cmap=cmap,
orientation='horizontal')
else:
cb = fig.colorbar(self.tmp_plot, cax=cbar_ax)
pu.adjust_log_ticks(cb, minlog, maxlog, symlog=symlog,
additional=additional)
cb.ax.tick_params(labelsize=fs)
[docs]
def add_phase_bar(self, fig, pos=[0.9, 0.12, 0.03, 0.75], dummy=False,
n_colors=3, fs=12, label='$\phi$ (°)', p_lim=1.,
horizontal=False, where='left', diff=0.05,
label_pos=[0., 1.02], delta=False, absdelta=False):
"""
Add colorbar to a 2D / 3D data plot at a custom position with custom
style options.
"""
if delta:
label = '$\Delta\phi$ (°)'
if dummy:
pos[1] -= 0.02
if where == 'left' or where == 'both':
fig.subplots_adjust(right=pos[0] - diff)
cbar_ax = fig.add_axes(pos)
if label is not None:
cbar_ax.text(label_pos[0], label_pos[1], label, fontsize=fs,
color=self.label_color)
if horizontal:
cb = fig.colorbar(self.tmp_plot, cax=cbar_ax,
orientation='horizontal')
else:
cb = fig.colorbar(self.tmp_plot, cax=cbar_ax)
if delta:
tickrange = np.linspace(-p_lim, p_lim, n_colors)
elif absdelta:
tickrange = np.linspace(0., p_lim, n_colors)
else:
tickrange = np.linspace(-180., 180., 9, dtype=int)
cb.set_ticks(tickrange)
tlabels = [str(tlabel) for tlabel in tickrange]
cb.set_ticklabels(tlabels)
cb.ax.tick_params(labelsize=fs)
[docs]
def amp_phase_comp(self, data, name=None, store=True, shift=0.):
"""
Calculate component-wise amplitudes and phases.
"""
amp, phase = np.zeros((len(data), 3)), np.zeros((len(data), 3))
amp[:, 0] = np.absolute(data[:, 0])
amp[:, 1] = np.absolute(data[:, 1])
amp[:, 2] = np.absolute(data[:, 2])
phase[:, 0] = np.angle(data[:, 0], deg=True) + shift
phase[:, 1] = np.angle(data[:, 1], deg=True) + shift
phase[:, 2] = np.angle(data[:, 2], deg=True) + shift
if shift != 0.:
phase[phase > 180.] -= 360.
return(amp, phase)
[docs]
def amp_phase_mag(self, data, name=None, store=True):
"""
Calculate amplitude and phase of vector magnitudes.
"""
vec = np.zeros((len(data)), dtype=complex)
vec.real = np.linalg.norm(data.real, axis=1)
vec.imag = np.linalg.norm(data.imag, axis=1)
return(np.abs(vec), np.angle(vec, deg=True))
[docs]
def rel_errors(self, data1, data2, name=None, store=True):
"""
Calculate relative errors / mismatch between two datasets.
"""
rel_comp_diff = np.zeros((len(data1), 3), dtype=complex)
rel_magnitude_diff = np.zeros(len(data1), dtype=complex)
rel_magnitude_diff.real = \
100. * ((np.linalg.norm(data2.real, axis=1) -
np.linalg.norm(data1.real, axis=1)) /
np.linalg.norm(data1.real, axis=1))
rel_magnitude_diff.imag = \
100. * ((np.linalg.norm(data2.imag, axis=1) -
np.linalg.norm(data1.imag, axis=1)) /
np.linalg.norm(data1.imag, axis=1))
rel_comp_diff.real = 100. * ((data2.real - data1.real) /
np.abs(data1.real))
rel_comp_diff.imag = 100. * ((data2.imag - data1.imag) /
np.abs(data1.imag))
if name is None:
name = str(len(self.rel_comp_errors))
if store:
self.rel_comp_errors[name] = rel_comp_diff
self.rel_mag_errors[name] = rel_magnitude_diff
else:
return(rel_comp_diff, rel_magnitude_diff)
[docs]
def abs_angles(self, data1, data2, name=None, store=True):
"""
Calculate angles between vectors
"""
angles = np.zeros((len(data1)), dtype=complex)
cos_a_real = np.sum(data1.real * data2.real, axis=1)
cos_a_imag = np.sum(data1.imag * data2.imag, axis=1)
scale_real = (np.linalg.norm(data1.real, axis=1) *
np.linalg.norm(data2.real, axis=1))
scale_imag = (np.linalg.norm(data1.imag, axis=1) *
np.linalg.norm(data2.imag, axis=1))
angles.real = np.rad2deg(np.arccos(cos_a_real / scale_real))
angles.imag = np.rad2deg(np.arccos(cos_a_imag / scale_imag))
if np.mean(angles.real) > 90.:
angles.real = 180. - angles.real
if np.mean(angles.imag) > 90.:
angles.imag = 180. - angles.imag
if name is None:
name = 'None'
if store:
print('Warning, storing angles not implemented so far!')
else:
return(angles)
[docs]
def abs_errors(self, data1, data2, name=None, store=True):
"""
Calculate absolute errors / mismatch between two datasets.
"""
abs_comp_diff = data2 - data1
abs_magnitude_diff = np.zeros((len(data1)), dtype=complex)
abs_magnitude_diff.real = np.linalg.norm(abs_comp_diff.real, axis=1)
abs_magnitude_diff.imag = np.linalg.norm(abs_comp_diff.imag, axis=1)
# abs_comp_diff.real /= abs_magnitude_diff.real.reshape(-1, 1)
# abs_comp_diff.imag /= abs_magnitude_diff.imag.reshape(-1, 1)
if name is None:
name = str(len(self.abs_comp_errors))
if store:
self.abs_comp_errors[name] = abs_comp_diff
self.abs_mag_errors[name] = abs_magnitude_diff
else:
return(abs_comp_diff, abs_magnitude_diff)
[docs]
def mean_errors(self, key, comp='mag', err_type='rel'):
"""
Calculate Mean of abs(relative errors between two datasets).
"""
if err_type == 'rel':
data = self.rel_mag_errors
real = np.mean(np.abs(data[key].real))
imag = np.mean(np.abs(data[key].imag))
return(real, imag)
[docs]
def ratio(self, data1, data2, name=None, store=True):
"""
Calculate ratio between two datasets, e.g.,
primary field magnitudes / secondary field magnitudes.
"""
ratio_comp_diff = np.zeros((len(data1), 3), dtype=complex)
ratio_magnitude_diff = np.zeros(len(data1), dtype=complex)
ratio_comp_diff.real = data2.real / data1.real
ratio_comp_diff.imag = data2.imag / data1.imag
ratio_magnitude_diff.real = (np.linalg.norm(data2.real, axis=1) /
np.linalg.norm(data1.real, axis=1))
ratio_magnitude_diff.imag = (np.linalg.norm(data2.imag, axis=1) /
np.linalg.norm(data1.imag, axis=1))
if name is None:
name = str(len(self.ratio_comp_errors))
if store:
self.ratio_comp_errors[name] = 100 * ratio_comp_diff
self.ratio_mag_errors[name] = 100 * ratio_magnitude_diff
else:
return(100 * ratio_comp_diff, 100. * ratio_magnitude_diff)
[docs]
def save_plot(self, EH, title, save, s_name, fs=10, mesh=None, png=False):
if mesh is None:
mesh = ''
if title is not None and save:
plt.suptitle(title, fontsize=fs)
if png:
plt.savefig(self.s_dir + '/' + EH + '_' + title + '_' + mesh +
'.png', bbox_inches='tight',
dpi=self.dpi, transparent=True)
else:
plt.savefig(self.s_dir + '/' + EH + '_' + title + '_' + mesh +
'.pdf', bbox_inches='tight')
if s_name is not None and save:
if png:
plt.savefig(self.s_dir + '/' + EH + '_' + s_name + '.png',
bbox_inches='tight', dpi=self.dpi,
transparent=True)
else:
plt.savefig(self.s_dir + '/' + EH + '_' + s_name + '.pdf',
bbox_inches='tight')
if title is None and s_name is None and save:
print('Warning!, specify either a *title* or *s_name* to save '
'figures!')