From af5ba53d058904c5ef9d72571f8447dfcc0de71c Mon Sep 17 00:00:00 2001 From: Dongdong Tian Date: Wed, 2 Oct 2024 14:40:04 +0800 Subject: [PATCH] Add private _FocalMechanismConvention class and refactor Figure.meca --- pygmt/src/_common.py | 212 ++++++++++++++++++++++++++++++++++++++++- pygmt/src/meca.py | 221 +++++-------------------------------------- 2 files changed, 233 insertions(+), 200 deletions(-) diff --git a/pygmt/src/_common.py b/pygmt/src/_common.py index 96ae9791b71..cd35e08fc70 100644 --- a/pygmt/src/_common.py +++ b/pygmt/src/_common.py @@ -2,9 +2,11 @@ Common functions used in multiple PyGMT functions/methods. """ +from collections.abc import Sequence from pathlib import Path -from typing import Any +from typing import Any, ClassVar, Literal +from pygmt.exceptions import GMTInvalidInput from pygmt.src.which import which @@ -39,3 +41,211 @@ def _data_geometry_is_point(data: Any, kind: str) -> bool: except FileNotFoundError: pass return False + + +class _FocalMechanismConvention: + """ + Class to handle focal mechanism convention, code, and associated parameters. + + Parameters + ---------- + convention + The focal mechanism convention. Valid values are: + + - ``"aki"``: Aki and Richards convention. + - ``"gcmt"``: Global CMT (Centroid Moment Tensor) convention. + - ``"partial"``: Partial focal mechanism convention. + - ``"mt"``: Moment Tensor convention. + - ``"principal_axis"``: Principal axis convention. + component + The component of the seismic moment tensor to plot. Valid values are: + + - ``"full"``: the full tensor seismic moment tensor + - ``"dc"``: the closest double coupe defined from the moment tensor (zero trace + and zero determinant) + - ``"deviatoric"``: deviatoric part of the moment tensor (zero trace) + + Only valid for conventions ``"mt"`` and ``"principal_axis"``. + + Attributes + ---------- + convention + The focal mechanism convention. + params + List of parameters associated with the focal mechanism convention. + code + The single-letter code that can be used in meca/coupe's -S option. + + Examples + -------- + >>> from pygmt.src._common import _FocalMechanismConvention + + >>> conv = _FocalMechanismConvention("aki") + >>> conv.convention, conv.code + ('aki', 'a') + >>> conv.params + ['strike', 'dip', 'rake', 'magnitude'] + + >>> conv = _FocalMechanismConvention("gcmt") + >>> conv.convention, conv.code + ('gcmt', 'c') + >>> conv.params + ['strike1', 'dip1', 'rake1', 'strike2', 'dip2', 'rake2', 'mantissa', 'exponent'] + + >>> conv = _FocalMechanismConvention("partial") + >>> conv.convention, conv.code + ('partial', 'p') + >>> conv.params + ['strike1', 'dip1', 'strike2', 'fault_type', 'magnitude'] + + >>> conv = _FocalMechanismConvention("mt", component="dc") + >>> conv.convention, conv.code + ('mt', 'd') + ['mrr', 'mtt', 'mff', 'mrt', 'mrf', 'mtf', 'exponent'] + + >>> conv = _FocalMechanismConvention("principal_axis", component="deviatoric") + >>> conv.convention, conv.code + ('principal_axis', 't') + + >>> conv = _FocalMechanismConvention("a") + >>> conv.convention, conv.code + ('aki', 'a') + >>> conv.params + ['strike', 'dip', 'rake', 'magnitude'] + + >>> conv = _FocalMechanismConvention.from_params( + ... ["strike", "dip", "rake", "magnitude"] + ... ) + >>> conv.convention, conv.code + ('aki', 'a') + + >>> conv = _FocalMechanismConvention(convention="invalid") + Traceback (most recent call last): + ... + GMTInvalidInput: Invalid focal mechanism convention 'invalid'. + + >>> conv = _FocalMechanismConvention("mt", component="invalid") + Traceback (most recent call last): + ... + GMTInvalidInput: Invalid component 'invalid' for focal mechanism convention 'mt'. + + >>> _FocalMechanismConvention.from_params(["strike", "dip", "rake"]) + Traceback (most recent call last): + ... + GMTInvalidInput: Fail to determine focal mechanism convention from ... + """ + + # Mapping of focal mechanism conventions to their single-letter codes. + _conventions: ClassVar = { + "aki": "a", + "gcmt": "c", + "partial": "p", + "mt": {"full": "m", "deviatoric": "z", "dc": "d"}, + "principal_axis": {"full": "x", "deviatoric": "t", "dc": "y"}, + } + + # Mapping of single-letter codes to focal mechanism convention names + _codes: ClassVar = { + "a": "aki", + "c": "gcmt", + "p": "partial", + "m": "mt", + "z": "mt", + "d": "mt", + "x": "principal_axis", + "t": "principal_axis", + "y": "principal_axis", + } + + # Mapping of focal mechanism conventions to their parameters. + _params: ClassVar = { + "aki": ["strike", "dip", "rake", "magnitude"], + "gcmt": [ + "strike1", + "dip1", + "rake1", + "strike2", + "dip2", + "rake2", + "mantissa", + "exponent", + ], + "partial": ["strike1", "dip1", "strike2", "fault_type", "magnitude"], + "mt": ["mrr", "mtt", "mff", "mrt", "mrf", "mtf", "exponent"], + "principal_axis": [ + "t_value", + "t_azimuth", + "t_plunge", + "n_value", + "n_azimuth", + "n_plunge", + "p_value", + "p_azimuth", + "p_plunge", + "exponent", + ], + } + + def __init__( + self, + convention: Literal["aki", "gcmt", "partial", "mt", "principal_axis"], + component: Literal["full", "deviatoric", "dc"] = "full", + ): + """ + Initialize the FocalMechanismConvention object. + """ + if convention in self._conventions: + # Convention is given via 'convention' and 'component' parameters. + if component not in {"full", "deviatoric", "dc"}: + msg = ( + f"Invalid component '{component}' for focal mechanism convention " + f"'{convention}'." + ) + raise GMTInvalidInput(msg) + + self.convention = convention + self.code = self._conventions[convention] + if isinstance(self.code, dict): + self.code = self.code[component] + elif convention in self._codes: + # Convention is given as a single-letter code. + self.code = convention + self.convention = self._codes[convention] + else: + msg = f"Invalid focal mechanism convention '{convention}'." + raise GMTInvalidInput(msg) + self.params = self._params[self.convention] + + @staticmethod + def from_params( + params: Sequence[str], component: Literal["full", "deviatoric", "dc"] = "full" + ) -> "_FocalMechanismConvention": + """ + Create a FocalMechanismConvention object from a sequence of parameters. + + The method checks if the given parameters are a superset of a known focal + mechanism convention to determine the convention. If the parameters are not + sufficient to determine the convention, an exception is raised. + + Parameters + ---------- + params + Sequence of parameters to determine the focal mechanism convention. The + order of the parameters does not matter. + + Returns + ------- + _FocalMechanismConvention + The FocalMechanismConvention object. + + Raises + ------ + GMTInvalidInput + If the focal mechanism convention cannot be determined from the given + parameters + """ + for convention, param_list in _FocalMechanismConvention._params.items(): + if set(param_list).issubset(set(params)): + return _FocalMechanismConvention(convention, component=component) + msg = "Fail to determine focal mechanism convention from the data column names." + raise GMTInvalidInput(msg) diff --git a/pygmt/src/meca.py b/pygmt/src/meca.py index 9822dae1f23..4d5b704c74d 100644 --- a/pygmt/src/meca.py +++ b/pygmt/src/meca.py @@ -5,183 +5,9 @@ import numpy as np import pandas as pd from pygmt.clib import Session -from pygmt.exceptions import GMTError, GMTInvalidInput +from pygmt.exceptions import GMTInvalidInput from pygmt.helpers import build_arg_list, fmt_docstring, kwargs_to_strings, use_alias - - -def convention_code(convention, component="full"): - """ - Determine the convention code for focal mechanisms. - - The convention code can be used in meca's -S option. - - Parameters - ---------- - convention : str - The focal mechanism convention. Can be one of the following: - - - ``"aki"``: Aki and Richards - - ``"gcmt"``: Global Centroid Moment Tensor - - ``"partial"``: Partial focal mechanism - - ``"mt"``: Moment tensor - - ``"principal_axis"``: Principal axis - - Single letter convention codes like ``"a"`` and ``"c"`` are also - supported but undocumented. - - component : str - The component of the focal mechanism. Only used when ``convention`` is - ``"mt"`` or ``"principal_axis"``. Can be one of the following: - - - ``"full"``: Full moment tensor - - ``"deviatoric"``: Deviatoric moment tensor - - ``"dc"``: Double couple - - Returns - ------- - str - The single-letter convention code used in meca's -S option. - - Examples - -------- - >>> convention_code("aki") - 'a' - >>> convention_code("gcmt") - 'c' - >>> convention_code("partial") - 'p' - - >>> convention_code("mt", component="full") - 'm' - >>> convention_code("mt", component="deviatoric") - 'z' - >>> convention_code("mt", component="dc") - 'd' - >>> convention_code("principal_axis", component="full") - 'x' - >>> convention_code("principal_axis", component="deviatoric") - 't' - >>> convention_code("principal_axis", component="dc") - 'y' - - >>> for code in ["a", "c", "m", "d", "z", "p", "x", "y", "t"]: - ... assert convention_code(code) == code - - >>> convention_code("invalid") - Traceback (most recent call last): - ... - pygmt.exceptions.GMTInvalidInput: Invalid convention 'invalid'. - - >>> convention_code("mt", "invalid") # doctest: +NORMALIZE_WHITESPACE - Traceback (most recent call last): - ... - pygmt.exceptions.GMTInvalidInput: - Invalid component 'invalid' for convention 'mt'. - """ - # Codes for focal mechanism formats determined by "convention" - codes1 = {"aki": "a", "gcmt": "c", "partial": "p"} - # Codes for focal mechanism formats determined by both "convention" and - # "component" - codes2 = { - "mt": {"deviatoric": "z", "dc": "d", "full": "m"}, - "principal_axis": {"deviatoric": "t", "dc": "y", "full": "x"}, - } - - if convention in codes1: - return codes1[convention] - if convention in codes2: - if component not in codes2[convention]: - raise GMTInvalidInput( - f"Invalid component '{component}' for convention '{convention}'." - ) - return codes2[convention][component] - if convention in {"a", "c", "m", "d", "z", "p", "x", "y", "t"}: - return convention - raise GMTInvalidInput(f"Invalid convention '{convention}'.") - - -def convention_name(code): - """ - Determine the name of a focal mechanism convention from its code. - - Parameters - ---------- - code : str - The single-letter convention code. - - Returns - ------- - str - The name of the focal mechanism convention. - - Examples - -------- - >>> convention_name("a") - 'aki' - >>> convention_name("aki") - 'aki' - """ - name = { - "a": "aki", - "c": "gcmt", - "p": "partial", - "z": "mt", - "d": "mt", - "m": "mt", - "x": "principal_axis", - "y": "principal_axis", - "t": "principal_axis", - }.get(code) - return name if name is not None else code - - -def convention_params(convention): - """ - Return the list of focal mechanism parameters for a given convention. - - Parameters - ---------- - convention : str - The focal mechanism convention. Can be one of the following: - - - ``"aki"``: Aki and Richards - - ``"gcmt"``: Global Centroid Moment Tensor - - ``"partial"``: Partial focal mechanism - - ``"mt"``: Moment tensor - - ``"principal_axis"``: Principal axis - - Returns - ------- - list - The list of focal mechanism parameters. - """ - return { - "aki": ["strike", "dip", "rake", "magnitude"], - "gcmt": [ - "strike1", - "dip1", - "rake1", - "strike2", - "dip2", - "rake2", - "mantissa", - "exponent", - ], - "mt": ["mrr", "mtt", "mff", "mrt", "mrf", "mtf", "exponent"], - "partial": ["strike1", "dip1", "strike2", "fault_type", "magnitude"], - "principal_axis": [ - "t_value", - "t_azimuth", - "t_plunge", - "n_value", - "n_azimuth", - "n_plunge", - "p_value", - "p_azimuth", - "p_plunge", - "exponent", - ], - }[convention] +from pygmt.src._common import _FocalMechanismConvention @fmt_docstring @@ -204,7 +30,7 @@ def convention_params(convention): t="transparency", ) @kwargs_to_strings(R="sequence", c="sequence_comma", p="sequence") -def meca( # noqa: PLR0912, PLR0913, PLR0915 +def meca( # noqa: PLR0912, PLR0913 self, spec, scale, @@ -401,17 +227,10 @@ def meca( # noqa: PLR0912, PLR0913, PLR0915 # Convert spec to pandas.DataFrame unless it's a file if isinstance(spec, dict | pd.DataFrame): # spec is a dict or pd.DataFrame - # determine convention from dict keys or pd.DataFrame column names - for conv in ["aki", "gcmt", "mt", "partial", "principal_axis"]: - if set(convention_params(conv)).issubset(set(spec.keys())): - convention = conv - break - else: - if isinstance(spec, dict): - msg = "Keys in dict 'spec' do not match known conventions." - else: - msg = "Column names in pd.DataFrame 'spec' do not match known conventions." - raise GMTError(msg) + # Determine convention from dict keys or pd.DataFrame column names + _convention = _FocalMechanismConvention.from_params( + spec.keys(), component=component + ) # convert dict to pd.DataFrame so columns can be reordered if isinstance(spec, dict): @@ -423,13 +242,16 @@ def meca( # noqa: PLR0912, PLR0913, PLR0915 ) elif isinstance(spec, np.ndarray): # spec is a numpy array if convention is None: - raise GMTInvalidInput("'convention' must be specified for an array input.") - # make sure convention is a name, not a code - convention = convention_name(convention) + msg = "'convention' must be specified for an array input." + raise GMTInvalidInput(msg) + + _convention = _FocalMechanismConvention( + convention=convention, component=component + ) # Convert array to pd.DataFrame and assign column names spec = pd.DataFrame(np.atleast_2d(spec)) - colnames = ["longitude", "latitude", "depth", *convention_params(convention)] + colnames = ["longitude", "latitude", "depth", *_convention.params] # check if spec has the expected number of columns ncolsdiff = len(spec.columns) - len(colnames) if ncolsdiff == 0: @@ -441,10 +263,15 @@ def meca( # noqa: PLR0912, PLR0913, PLR0915 elif ncolsdiff == 3: colnames += ["plot_longitude", "plot_latitude", "event_name"] else: - raise GMTInvalidInput( + msg = ( f"Input array must have {len(colnames)} to {len(colnames) + 3} columns." ) + raise GMTInvalidInput(msg) spec.columns = colnames + else: + _convention = _FocalMechanismConvention( + convention=convention, component=component + ) # Now spec is a pd.DataFrame or a file if isinstance(spec, pd.DataFrame): @@ -472,7 +299,7 @@ def meca( # noqa: PLR0912, PLR0913, PLR0915 # expected columns are: # longitude, latitude, depth, focal_parameters, # [plot_longitude, plot_latitude] [event_name] - newcols = ["longitude", "latitude", "depth", *convention_params(convention)] + newcols = ["longitude", "latitude", "depth", *_convention.params] if "plot_longitude" in spec.columns and "plot_latitude" in spec.columns: newcols += ["plot_longitude", "plot_latitude"] if kwargs.get("A") is None: @@ -483,11 +310,7 @@ def meca( # noqa: PLR0912, PLR0913, PLR0915 if spec.columns.tolist() != newcols: spec = spec.reindex(newcols, axis=1) - # determine data_format from convention and component - data_format = convention_code(convention=convention, component=component) - - # Assemble -S flag - kwargs["S"] = f"{data_format}{scale}" + kwargs["S"] = f"{_convention.code}{scale}" with Session() as lib: with lib.virtualfile_in(check_kind="vector", data=spec) as vintbl: lib.call_module(module="meca", args=build_arg_list(kwargs, infile=vintbl))