ogstools.meshlib.mesh_series module#
- class ogstools.meshlib.mesh_series.MeshSeries[source]#
Bases:
object
A wrapper around pyvista and meshio for reading of pvd and xdmf timeseries.
Initialize a MeshSeries object
- param filepath:
Path to the PVD or XDMF file.
- param time_unit:
Data unit of the timevalues.
- param data_length_unit:
Length unit of the mesh data.
- param output_length_unit:
Length unit in plots.
- returns:
A MeshSeries object
- __init__(filepath, time_unit='s', spatial_unit='m', spatial_output_unit='m')[source]#
Initialize a MeshSeries object
- param filepath:
Path to the PVD or XDMF file.
- param time_unit:
Data unit of the timevalues.
- param data_length_unit:
Length unit of the mesh data.
- param output_length_unit:
Length unit in plots.
- returns:
A MeshSeries object
- data(variable_name)[source]#
Returns an DataItems object, that allows array indexing. To get “geometry”/”points” or “topology”/”cells” read the first time step and use pyvista functionality Selection example: ms = MeshSeries() temp = ms.data(“temperature”) time_step1_temps = temp[1,:] temps_at_some_points = temp[:,1:3] :param variable_name: Name the variable (e.g.”temperature”) :returns: Returns an objects that allows array indexing.
- Return type:
- aggregate(variable, np_func, axis, mask=None)[source]#
Aggregate data of all timesteps using a specified function.
- Parameters:
variable (Variable | str) – The mesh variable to be aggregated.
func – The aggregation function to apply.
- Returns:
A numpy.ndarray of the same length as the timesteps if axis=0 or of the same length as the data if axis=1.
- Return type:
ndarray
- aggregate_over_time(variable, func)[source]#
Aggregate data over all timesteps using a specified function.
- closest_timestep(timevalue)[source]#
Return the corresponding timestep from a timevalue.
- Return type:
int
- closest_timevalue(timevalue)[source]#
Return the closest timevalue to a timevalue.
- Return type:
float
- mesh(timestep, lazy_eval=True)[source]#
Selects mesh at given timestep all data function.
- Return type:
- rawdata_file()[source]#
Checks, if working with the raw data is possible. For example, OGS Simulation results with XDMF support efficient raw data access via h5py
- Returns:
The location of the file containing the raw data. If it does not support efficient read (e.g., no efficient slicing), it returns None.
- Return type:
Path | None
- read_interp(timevalue, lazy_eval=True)[source]#
Return the temporal interpolated mesh for a given timevalue.
- Return type:
- property timesteps: range#
Return the timesteps of the timeseries data.
- timevalues(time_unit=None)[source]#
Return the timevalues, optionally converted to another time unit.
- Return type:
ndarray
- values(data_name, selection=None)[source]#
Get the data in the MeshSeries for all timesteps.
- Parameters:
data_name (str) – Name of the data in the MeshSeries.
selection (slice | ndarray | None) –
Can limit the data to be read. - Time is always the first dimension. - If None, it takes the selection that is defined in the xdmf file. - If a tuple or np.ndarray: see how h5py uses Numpy array indexing. - If a slice: see Python slice reference. - If a string: see example:
Example:
"|0 0 0:1 1 1:1 190 3:97 190 3"
This represents the selection
[(offset(0,0,0): step(1,1,1) : end(1,190,3) : of_data_with_size(97,190,30))]
.
- Returns:
A numpy array of the requested data for all timesteps.
- Return type:
ndarray
- time_of_min(variable)[source]#
Returns a Mesh with the time of the variable minimum as data.
- Return type:
- time_of_max(variable)[source]#
Returns a Mesh with the time of the variable maximum as data.
- Return type:
- aggregate_over_domain(variable, func, mask=None)[source]#
Aggregate data over domain per timestep using a specified function.
- Parameters:
variable (Variable | str) – The mesh variable to be aggregated.
func (Literal['min', 'max', 'mean', 'median', 'sum', 'std', 'var']) – The aggregation function to apply.
mask (ndarray | None) – A numpy array as a mask for the domain.
- Returns:
A numpy array with aggregated data.
- Return type:
ndarray
- plot_domain_aggregate(variable, func, timesteps=None, time_unit='s', mask=None, ax=None, **kwargs)[source]#
Plot the transient aggregated data over the domain per timestep.
- Parameters:
variable (Variable | str) – The mesh variable to be aggregated.
func (Literal['min', 'max', 'mean', 'median', 'sum', 'std', 'var']) – The aggregation function to apply.
timesteps (slice | None) – A slice to select the timesteps. Default: all.
time_unit (str | None) – Output unit of the timevalues.
mask (ndarray | None) – A numpy array as a mask for the domain.
ax (Axes | None) – matplotlib axis to use for plotting
kwargs (Any) – Keyword args passed to matplotlib’s plot function.
- Returns:
A matplotlib Figure or None if plotting on existing axis.
- Return type:
Figure | None
- probe(points, data_name, interp_method='linear')[source]#
Probe the MeshSeries at observation points.
- Parameters:
points (ndarray) – The observation points to sample at.
data_name (str) – Name of the data to sample.
interp_method (Literal['nearest', 'linear']) – Interpolation method, defaults to linear
- Returns:
numpy array of interpolated data at observation points.
- Return type:
ndarray
- plot_probe(points, variable, variable_abscissa=None, labels=None, time_unit='s', interp_method='linear', colors=None, linestyles=None, ax=None, fill_between=False, **kwargs)[source]#
Plot the transient variable on the observation points in the MeshSeries.
- param points:
The points to sample at.
- param variable:
The variable to be sampled.
- param labels:
The labels for each observation point.
- param time_unit:
Output unit of the timevalues.
- param interp_method:
Choose the interpolation method, defaults to linear for xdmf MeshSeries and probefilter for pvd MeshSeries.
- param interp_backend:
Interpolation backend for PVD MeshSeries.
- param kwargs:
Keyword arguments passed to matplotlib’s plot function.
- returns:
A matplotlib Figure
- Return type:
Figure | None
- animate(variable, timesteps=None, mesh_func=lambda mesh: ..., plot_func=lambda *_: ..., **kwargs)[source]#
Create an animation for a variable with given timesteps.
- Parameters:
variable (Variable) – the field to be visualized on all timesteps
timesteps (Sequence | None) – if sequence of int: the timesteps to animate if sequence of float: the timevalues to animate
mesh_func (Callable[[Mesh], Mesh]) – A function which expects to read a mesh and return a mesh. Useful, for slicing / clipping / thresholding the meshseries for animation.
plot_func (Callable[[Axes, float], None]) – A function which expects to read a matplotlib Axes and the time value of the current frame. Useful to customize the plot in the animation.
- Keyword Arguments:
- Return type:
FuncAnimation
- plot_time_slice(variable, points, y_axis='auto', interpolate=True, time_unit='s', time_logscale=False, fig=None, ax=None, cbar=True, **kwargs)[source]#
- Parameters:
variable (Variable | str) – The variable to be visualized.
points (ndarray) – The points along which the data is sampled over time.
y_axis (Literal['x', 'y', 'z', 'dist', 'auto']) – The component of the sampling points which labels the y-axis. By default, if only one coordinate of the points is changing, this axis is taken, otherwise the distance along the line is taken.
interpolate (bool) – Smoothen the result be interpolation.
time_unit (str) – Time unit displayed on the x-axis.
time_logscale (bool) – Should log-scaling be applied to the time-axis?
fig (Figure | None) – matplotlib figure to use for plotting.
ax (Axes | None) – matplotlib axis to use for plotting.
cb_loc – Colorbar location. If None, omit colorbar.
- Keyword Arguments:
cb_labelsize: colorbar labelsize
cb_loc: colorbar location (‘left’ or ‘right’)
cb_pad: colorbar padding
cmap: colormap
vmin: minimum value for colorbar
vmax: maximum value for colorbar
num_levels: number of levels for colorbar
figsize: figure size
dpi: resolution
- Return type:
Figure