rsqsim_api.io.read_utils ======================== .. py:module:: rsqsim_api.io.read_utils .. autoapi-nested-parse:: Utilities for reading RSQSim output files in text, binary, CSV, DXF, STL, and VTK formats. Attributes ---------- .. autoapisummary:: rsqsim_api.io.read_utils.catalogue_columns Functions --------- .. autoapisummary:: rsqsim_api.io.read_utils.read_text rsqsim_api.io.read_utils.read_binary rsqsim_api.io.read_utils.read_csv_and_array rsqsim_api.io.read_utils.read_earthquakes rsqsim_api.io.read_utils.read_earthquake_catalogue rsqsim_api.io.read_utils.read_ts_coords rsqsim_api.io.read_utils.read_dxf rsqsim_api.io.read_utils.read_stl rsqsim_api.io.read_utils.read_vtk Module Contents --------------- .. py:data:: catalogue_columns :value: ['t0', 'm0', 'mw', 'x', 'y', 'z', 'area', 'dt'] .. py:function:: read_text(file, format) Read scalar values from a text file output by RSQSim. Produced when RSQSim is compiled with serial (text) output mode. :param file: Path to the text file to read. :param format: Data type specifier: ``"d"`` for double (float64) or ``"i"`` for integer (int32). :returns: Flattened 1-D array of values read from the file. :rtype: numpy.ndarray :raises AssertionError: If ``format`` is not ``"d"`` or ``"i"``, or if ``file`` does not exist. .. py:function:: read_binary(file, format, endian = 'little') Read scalar values from a binary file output by RSQSim. :param file: Path to the binary file to read. :param format: Data type specifier: ``"d"`` for double (float64) or ``"i"`` for integer (int32). :param endian: Byte order of the file. ``"little"`` (default) is standard for most modern x86 systems; use ``"big"`` for non-standard platforms. :returns: Flattened 1-D array of values read from the file. :rtype: numpy.ndarray :raises AssertionError: If ``endian`` is not ``"little"`` or ``"big"``, if ``format`` is not ``"d"`` or ``"i"``, or if ``file`` does not exist. .. py:function:: read_csv_and_array(prefix, read_index = True) Read a catalogue CSV and its associated NumPy array files from a common prefix. Expects files named ``_catalogue.csv``, ``_events.npy``, ``_patches.npy``, ``_slip.npy``, and ``_slip_time.npy`` to exist on disk. :param prefix: File path prefix (with or without trailing underscore). :param read_index: If ``True`` (default), read the first column of the CSV as the DataFrame index. :type read_index: bool, optional :returns: ``[df, events, patches, slip, slip_time]`` where ``df`` is a ``pandas.DataFrame`` and the remaining elements are NumPy arrays loaded from the corresponding ``.npy`` files. :rtype: list :raises AssertionError: If ``prefix`` is an empty string. :raises FileNotFoundError: If any of the five expected files is missing. .. py:function:: read_earthquakes(earthquake_file, get_patch = False, eq_start_index = None, eq_end_index = None, endian = 'little') Read earthquake event data, inferring companion file names from the event file prefix. Based on R scripts by Keith Richards-Dinger. The earthquake file is expected to follow the naming convention ``eqs..out``. :param earthquake_file: Path to the RSQSim earthquake output file, usually with a ``.out`` suffix and named ``eqs..out``. :param get_patch: If ``True``, also read per-patch rupture data. Currently unused placeholder. :param eq_start_index: Zero-based index of the first earthquake to read. Both ``eq_start_index`` and ``eq_end_index`` must be provided together. :param eq_end_index: Zero-based index of the last earthquake to read (exclusive). :param endian: Byte order for binary companion files. Defaults to ``"little"``. :raises AssertionError: If ``endian`` is not ``"little"`` or ``"big"``, or if ``earthquake_file`` does not exist. :raises ValueError: If ``eq_start_index >= eq_end_index`` when both are provided. .. py:function:: read_earthquake_catalogue(catalogue_file) Read an RSQSim earthquake catalogue text file into a DataFrame. Parses the catalogue file produced by RSQSim, skipping the header block that ends with the line ``%%% end input files``, and loads the subsequent rows into a DataFrame with the standard catalogue columns ``t0, m0, mw, x, y, z, area, dt``. :param catalogue_file: Path to the RSQSim catalogue text file. :returns: DataFrame with columns ``["t0", "m0", "mw", "x", "y", "z", "area", "dt"]``, one row per earthquake event. :rtype: pandas.DataFrame :raises AssertionError: If ``catalogue_file`` does not exist. .. py:function:: read_ts_coords(filename) Read vertex and triangle data from a tsurf (``.ts``) file. Parses the SCEC Community Fault Model tsurf format, extracting vertex coordinates and triangle connectivity. Based on the MATLAB script ``ReadAndSaveCfm.m`` by Brendan Meade. :param filename: Path to the tsurf ``.ts`` file. :returns: * **vrtx** (*numpy.ndarray of shape (n_vertices, 4)*) -- Vertex data array. Each row is ``[vertex_id, x, y, z]``. * **trgl** (*numpy.ndarray of shape (n_triangles, 3), dtype int*) -- Triangle connectivity array. Each row gives the three vertex IDs forming a triangle. * **tri** (*numpy.ndarray of shape (n_triangles, 9)*) -- Triangle data array. Each row contains the (x, y, z) coordinates of all three corners concatenated: ``[x1,y1,z1, x2,y2,z2, x3,y3,z3]``. .. admonition:: Notes Copyright Paul Kaeufl, July 2014. Original MATLAB script: http://structure.rc.fas.harvard.edu/cfm/download/meade/ReadAndSaveCfm.m .. py:function:: read_dxf(dxf_file) Read a triangulated mesh and boundary polyline from a DXF file. Expects a DXF file exported from Move containing exactly one ``POLYLINE`` boundary and one or more ``3DFACE`` triangles. :param dxf_file: Path to the DXF file. :returns: * **triangle_array** (*numpy.ndarray of shape (n_triangles, 9)*) -- Each row contains the (x, y, z) coordinates of the three triangle corners concatenated: ``[x1,y1,z1, x2,y2,z2, x3,y3,z3]``. * **boundary_array** (*numpy.ndarray of shape (n_points, 3)*) -- (x, y, z) coordinates of the boundary polyline vertices. :raises AssertionError: If ``dxf_file`` does not exist, or if the file does not contain both ``3DFACE`` and ``POLYLINE`` entities. :raises ValueError: If the file contains more than one ``POLYLINE`` boundary. .. py:function:: read_stl(stl_file, min_point_sep=100.0) Read a triangulated surface mesh from an STL file. Merges near-duplicate vertices (within ``min_point_sep``) by averaging their positions and removing degenerate triangles that result from the merge. :param stl_file: Path to the STL mesh file. :param min_point_sep: Minimum separation in metres below which two vertices are considered duplicates and merged. Defaults to 100.0 m. :returns: Each row contains the (x, y, z) coordinates of the three triangle corners: ``[x1,y1,z1, x2,y2,z2, x3,y3,z3]``. :rtype: numpy.ndarray of shape (n_triangles, 9) :raises AssertionError: If ``stl_file`` does not exist or if the mesh contains no triangular cells. .. py:function:: read_vtk(vtk_file, min_point_sep=1.0) Read a triangulated surface mesh with slip and rake data from a VTK file. Merges near-duplicate vertices (within ``min_point_sep``) and removes degenerate triangles, preserving the associated slip and rake cell data. :param vtk_file: Path to the VTK mesh file. :param min_point_sep: Minimum separation in metres below which two vertices are considered duplicates and merged. Defaults to 1.0 m. :returns: * **mesh_as_array** (*numpy.ndarray of shape (n_triangles, 9)*) -- Each row contains the (x, y, z) coordinates of the three triangle corners: ``[x1,y1,z1, x2,y2,z2, x3,y3,z3]``. * **slip** (*numpy.ndarray of shape (n_triangles,)*) -- Slip magnitude (metres) for each triangle after duplicate removal. * **rake** (*numpy.ndarray of shape (n_triangles,)*) -- Rake angle (degrees) for each triangle after duplicate removal. :raises AssertionError: If ``vtk_file`` does not exist, if the mesh contains no triangular cells, or if the ``slip`` and ``rake`` cell data arrays are absent.