Skip to content

Datasets module

Package Module scubas.datasets Class Layer, Site, PROFILES

Why this module matters

Most SCUBAS workflows pass Site objects into conductivity, model, and cable utilities. Defining these correctly improves downstream stability.

scubas.datasets now exposes lightweight data classes powered by dataclass and offers convenience helpers for synthesising layered conductivity profiles. The modernised layout makes it easier to construct and inspect a Site, while remaining compatible with existing utilities such as the conductivity and cable modules.

Highlights:

  • Layer is an immutable dataclass that computes resistivity lazily via a property.
  • Site accepts any iterable of Layer objects and provides helpers for thickness/conductivity/resistivity extraction as lists or individual values.
  • The Site.init static method builds a Site from parallel sequences of conductivities, thicknesses, and names—a pattern used widely in PROFILES.

Quick start

from scubas.datasets import Layer, Site, PROFILES

layers = [
    Layer(name="Seawater", thickness=1.0, conductivity=3.33),
    Layer(name="Sediments", thickness=2.5, conductivity=0.2),
]
site = Site(layers=layers, description="Example", name="Demo Site")

print(site.get_resistivities())

# Build from sequences
alt_site = Site.init(
    conductivities=[3.33, 0.2, 0.0003],
    thicknesses=[1.0, 2.5, float("inf")],
    names=["Seawater", "Sediments", "Mantle"],
    description="Three-layer example",
    site_name="Alt Site",
)

# Reuse predefined profiles
print(PROFILES.CS.get_names())

API reference

scubas.datasets.Layer dataclass

Single geophysical layer describing thickness and electrical properties.

Source code in scubas/datasets.py 
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
@dataclass
class Layer:
    """
    Single geophysical layer describing thickness and electrical properties.
    """

    name: str
    thickness: float
    conductivity: float

    @property
    def resistivity(self) -> float:
        """
        Resistivity (Ohm·m) derived as the reciprocal of conductivity.
        """
        if self.conductivity == 0.0:
            raise ValueError("Conductivity must be non-zero to compute resistivity.")
        return 1.0 / self.conductivity

resistivity: float property

Resistivity (Ohm·m) derived as the reciprocal of conductivity.

scubas.datasets.Site dataclass

Collection of layers associated with a location or conductivity profile.

Source code in scubas/datasets.py 
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
@dataclass
class Site:
    """
    Collection of layers associated with a location or conductivity profile.
    """

    layers: Sequence[Layer]
    description: str
    name: str

    def __post_init__(self) -> None:
        self.layers = list(self.layers)

    def get_thicknesses(self, index: int = -1) -> Union[List[float], float]:
        """
        Thickness values (km) for all layers or the specified index.
        """
        values = [layer.thickness for layer in self.layers]
        if index >= 0:
            return values[index]
        return values

    def get_conductivities(self, index: int = -1) -> Union[List[float], float]:
        """
        Conductivity values (S/m) for all layers or the specified index.
        """
        values = [layer.conductivity for layer in self.layers]
        if index >= 0:
            return values[index]
        return values

    def get_resistivities(self, index: int = -1) -> Union[List[float], float]:
        """
        Resistivity values (Ohm·m) derived from the layer conductivities.
        """
        values = [layer.resistivity for layer in self.layers]
        if index >= 0:
            return values[index]
        return values

    def get_names(self, index: int = -1) -> Union[List[str], str]:
        """
        Layer names for all layers or the specified index.
        """
        values = [layer.name for layer in self.layers]
        if index >= 0:
            return values[index]
        return values

    def get(self, index: int = -1) -> pd.DataFrame:
        """
        Return a pandas dataframe with layer metadata.
        """
        return pd.DataFrame(
            {
                "names": np.atleast_1d(self.get_names(index)),
                "conductivities": np.atleast_1d(self.get_conductivities(index)),
                "resistivities": np.atleast_1d(self.get_resistivities(index)),
                "thicknesses": np.atleast_1d(self.get_thicknesses(index)),
            }
        )

    def set_thicknesses(self, value: float, index: int = 0) -> "Site":
        """
        Mutate the thickness of a layer at the specified index.
        """
        self.layers[index].thickness = value
        return self

    def calcZ(
        self,
        freqs: Iterable[float],
        layer: int = 0,
        return_all_layers: bool = False,
    ) -> Union[np.ndarray, Tuple[np.ndarray, np.ndarray]]:
        """
        Compute magnetotelluric impedance for the layered model.
        """
        freqs = np.asarray(freqs, dtype=float)
        resistivities = np.asarray(self.get_resistivities(), dtype=float)
        thicknesses = np.asarray(self.get_thicknesses(), dtype=float)
        n_layers = len(resistivities)
        n_freqs = len(freqs)
        omega = 2 * np.pi * freqs
        complex_factor = 1j * omega * C.mu_0

        k = np.sqrt(1j * omega[np.newaxis, :] * C.mu_0 / resistivities[:, np.newaxis])
        Z = np.zeros(shape=(n_layers, n_freqs), dtype=complex)
        with np.errstate(divide="ignore", invalid="ignore"):
            Z[-1, :] = complex_factor / k[-1, :]
            r = np.zeros(shape=(n_layers, n_freqs), dtype=complex)
            for idx in range(n_layers - 2, -1, -1):
                ratio = k[idx, :] * Z[idx + 1, :] / complex_factor
                r[idx, :] = (1 - ratio) / (1 + ratio)
                exponent = np.exp(-2 * k[idx, :] * thicknesses[idx])
                numerator = complex_factor * (1 - r[idx, :] * exponent)
                denominator = k[idx, :] * (1 + r[idx, :] * exponent)
                Z[idx, :] = numerator / denominator
        if freqs[0] == 0.0:
            Z[:, 0] = 0.0
        Z_output = np.zeros(shape=(4, n_freqs), dtype=complex)
        Z_output[1, :] = Z[layer, :] * (1.0e-3 / C.mu_0)
        Z_output[2, :] = -Z_output[1, :]
        if return_all_layers:
            return Z_output, Z
        return Z_output

    def calcP(
        self,
        freqs: Iterable[float],
        layer: int = 0,
        return_Z: bool = False,
    ) -> Union[np.ndarray, Tuple[np.ndarray, np.ndarray, np.ndarray]]:
        """
        Calculate skin depth for the model at the specified frequencies.
        """
        freqs = np.asarray(freqs, dtype=float)
        omega = 2 * np.pi * freqs

        if return_Z:
            Z_output, Z_all = self.calcZ(
                freqs=freqs, layer=layer, return_all_layers=True
            )
        else:
            Z_output = self.calcZ(freqs=freqs, layer=layer, return_all_layers=False)

        skin_depth = Z_output[1, :] / (1j * omega * C.mu_0) / (1.0e-3 / C.mu_0)
        if return_Z:
            return skin_depth, Z_output, Z_all
        return skin_depth

    @staticmethod
    def init(
        conductivities: Sequence[float],
        thicknesses: Sequence[float],
        names: Sequence[str],
        description: str,
        site_name: str,
    ) -> "Site":
        """
        Initialise a site from parallel sequences of layer properties.
        """
        layers = [
            Layer(name, thickness, conductivity)
            for conductivity, thickness, name in zip(conductivities, thicknesses, names)
        ]
        return Site(layers=layers, description=description, name=site_name)

get_thicknesses(index=-1)

Thickness values (km) for all layers or the specified index.

Source code in scubas/datasets.py 
59
60
61
62
63
64
65
66
def get_thicknesses(self, index: int = -1) -> Union[List[float], float]:
    """
    Thickness values (km) for all layers or the specified index.
    """
    values = [layer.thickness for layer in self.layers]
    if index >= 0:
        return values[index]
    return values

get_conductivities(index=-1)

Conductivity values (S/m) for all layers or the specified index.

Source code in scubas/datasets.py 
68
69
70
71
72
73
74
75
def get_conductivities(self, index: int = -1) -> Union[List[float], float]:
    """
    Conductivity values (S/m) for all layers or the specified index.
    """
    values = [layer.conductivity for layer in self.layers]
    if index >= 0:
        return values[index]
    return values

get_resistivities(index=-1)

Resistivity values (Ohm·m) derived from the layer conductivities.

Source code in scubas/datasets.py 
77
78
79
80
81
82
83
84
def get_resistivities(self, index: int = -1) -> Union[List[float], float]:
    """
    Resistivity values (Ohm·m) derived from the layer conductivities.
    """
    values = [layer.resistivity for layer in self.layers]
    if index >= 0:
        return values[index]
    return values

get_names(index=-1)

Layer names for all layers or the specified index.

Source code in scubas/datasets.py 
86
87
88
89
90
91
92
93
def get_names(self, index: int = -1) -> Union[List[str], str]:
    """
    Layer names for all layers or the specified index.
    """
    values = [layer.name for layer in self.layers]
    if index >= 0:
        return values[index]
    return values

get(index=-1)

Return a pandas dataframe with layer metadata.

Source code in scubas/datasets.py 
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
def get(self, index: int = -1) -> pd.DataFrame:
    """
    Return a pandas dataframe with layer metadata.
    """
    return pd.DataFrame(
        {
            "names": np.atleast_1d(self.get_names(index)),
            "conductivities": np.atleast_1d(self.get_conductivities(index)),
            "resistivities": np.atleast_1d(self.get_resistivities(index)),
            "thicknesses": np.atleast_1d(self.get_thicknesses(index)),
        }
    )

set_thicknesses(value, index=0)

Mutate the thickness of a layer at the specified index.

Source code in scubas/datasets.py 
108
109
110
111
112
113
def set_thicknesses(self, value: float, index: int = 0) -> "Site":
    """
    Mutate the thickness of a layer at the specified index.
    """
    self.layers[index].thickness = value
    return self

calcZ(freqs, layer=0, return_all_layers=False)

Compute magnetotelluric impedance for the layered model.

Source code in scubas/datasets.py 
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
def calcZ(
    self,
    freqs: Iterable[float],
    layer: int = 0,
    return_all_layers: bool = False,
) -> Union[np.ndarray, Tuple[np.ndarray, np.ndarray]]:
    """
    Compute magnetotelluric impedance for the layered model.
    """
    freqs = np.asarray(freqs, dtype=float)
    resistivities = np.asarray(self.get_resistivities(), dtype=float)
    thicknesses = np.asarray(self.get_thicknesses(), dtype=float)
    n_layers = len(resistivities)
    n_freqs = len(freqs)
    omega = 2 * np.pi * freqs
    complex_factor = 1j * omega * C.mu_0

    k = np.sqrt(1j * omega[np.newaxis, :] * C.mu_0 / resistivities[:, np.newaxis])
    Z = np.zeros(shape=(n_layers, n_freqs), dtype=complex)
    with np.errstate(divide="ignore", invalid="ignore"):
        Z[-1, :] = complex_factor / k[-1, :]
        r = np.zeros(shape=(n_layers, n_freqs), dtype=complex)
        for idx in range(n_layers - 2, -1, -1):
            ratio = k[idx, :] * Z[idx + 1, :] / complex_factor
            r[idx, :] = (1 - ratio) / (1 + ratio)
            exponent = np.exp(-2 * k[idx, :] * thicknesses[idx])
            numerator = complex_factor * (1 - r[idx, :] * exponent)
            denominator = k[idx, :] * (1 + r[idx, :] * exponent)
            Z[idx, :] = numerator / denominator
    if freqs[0] == 0.0:
        Z[:, 0] = 0.0
    Z_output = np.zeros(shape=(4, n_freqs), dtype=complex)
    Z_output[1, :] = Z[layer, :] * (1.0e-3 / C.mu_0)
    Z_output[2, :] = -Z_output[1, :]
    if return_all_layers:
        return Z_output, Z
    return Z_output

calcP(freqs, layer=0, return_Z=False)

Calculate skin depth for the model at the specified frequencies.

Source code in scubas/datasets.py 
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
def calcP(
    self,
    freqs: Iterable[float],
    layer: int = 0,
    return_Z: bool = False,
) -> Union[np.ndarray, Tuple[np.ndarray, np.ndarray, np.ndarray]]:
    """
    Calculate skin depth for the model at the specified frequencies.
    """
    freqs = np.asarray(freqs, dtype=float)
    omega = 2 * np.pi * freqs

    if return_Z:
        Z_output, Z_all = self.calcZ(
            freqs=freqs, layer=layer, return_all_layers=True
        )
    else:
        Z_output = self.calcZ(freqs=freqs, layer=layer, return_all_layers=False)

    skin_depth = Z_output[1, :] / (1j * omega * C.mu_0) / (1.0e-3 / C.mu_0)
    if return_Z:
        return skin_depth, Z_output, Z_all
    return skin_depth

init(conductivities, thicknesses, names, description, site_name) staticmethod

Initialise a site from parallel sequences of layer properties.

Source code in scubas/datasets.py 
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
@staticmethod
def init(
    conductivities: Sequence[float],
    thicknesses: Sequence[float],
    names: Sequence[str],
    description: str,
    site_name: str,
) -> "Site":
    """
    Initialise a site from parallel sequences of layer properties.
    """
    layers = [
        Layer(name, thickness, conductivity)
        for conductivity, thickness, name in zip(conductivities, thicknesses, names)
    ]
    return Site(layers=layers, description=description, name=site_name)

scubas.datasets.PROFILES = SimpleNamespace(None={key: Site.init(None=definition) for (key, definition) in {'BM': dict(conductivities=[0.2, 0.0003333, 0.02, 0.1, 1.12201], thicknesses=[2000, 75000, 332000, 250000, np.inf], names=['Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Lower Mantle'], description="Ben's benchmark conductivity model.", site_name="Ben's Model"), 'OM': dict(conductivities=[3.3333333, 0.2, 0.0003333, 0.02, 0.1, 1.12201], thicknesses=[1000, 2000, 75000, 332000, 250000, np.inf], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Lower Mantle'], description="Ben's model with a 1 km oceanic layer.", site_name='Ocean Model'), 'DB': dict(conductivities=[5e-05, 0.005, 0.001, 0.01, 0.3333333], thicknesses=[15000, 10000, 125000, 200000, np.inf], names=['Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Lower Mantle'], description="Conductivity profile used in David's study.", site_name="David's Model"), 'Quebec': dict(conductivities=[5e-05, 0.005, 0.001, 0.01, 0.3333333], thicknesses=[15000, 10000, 125000, 200000, np.inf], names=['Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Lower Mantle'], description='Conductivity model representing Quebec.', site_name="David's Model"), 'UN': dict(conductivities=[0.2, 0.2, 0.2, 0.2, 1.12201], thicknesses=[2000, 75000, 332000, 250000, np.inf], names=['Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Lower Mantle'], description="Uniform conductivity model inspired by David's work.", site_name='Uniform Model'), 'CS': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[100, 3000, 20000, 140000, 246900, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='IGS continental shelf reference profile.', site_name='Continental Shelf'), 'SO': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[1000, 2000, 10000, 70000, 327000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Shallow ocean conductivity profile used in IGS analyses.', site_name='Shallow Ocean'), 'DO': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4000, 2000, 10000, 70000, 327000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Deep ocean reference profile for IGS studies.', site_name='Deep Ocean'), 'LD0': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[0, 1000, 20000, 140000, 249000, 250000, 340000], names=['0-Water', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Land-based conductivity profile from IGS study.', site_name='Land-0m-depth'), 'LD': dict(conductivities=[0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[1000, 20000, 140000, 249000, 250000, 340000], names=['Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Land-based conductivity profile from IGS study.', site_name='Land'), 'CS_W': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[100, 8000, 15000, 150000, 236900, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989: Continental shelf (western) section.', site_name='Continental Shelf West'), 'DO_1': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4000, 4000, 10000, 145000, 247000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 1.', site_name='Deep Ocean'), 'DO_2': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[5200, 2000, 10000, 140000, 252800, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 2.', site_name='Deep Ocean'), 'DO_3': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4000, 2000, 10000, 140000, 254000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 3.', site_name='Deep Ocean'), 'DO_4': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4800, 1000, 10000, 70000, 324200, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 4.', site_name='Deep Ocean'), 'DO_5': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4000, 1000, 10000, 70000, 324200, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 5.', site_name='Deep Ocean'), 'MAR': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[3000, 0, 10000, 25000, 372000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: Mid-Atlantic Ridge profile.', site_name='Mid-Atlantic Ridge'), 'DO_6': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[4500, 1500, 10000, 70000, 324000, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989 study: deep ocean variant 6.', site_name='Deep Ocean'), 'CS_E': dict(conductivities=[3.3333333, 0.3333333, 0.00033333333, 0.001, 0.01, 0.1, 1], thicknesses=[100, 3000, 20000, 120000, 266900, 250000, 340000], names=['Seawater', 'Sediments', 'Crust', 'Lithosphere', 'Upper Mantle', 'Transition Zone', 'Lower Mantle'], description='Space Weather 1989: Continental shelf (eastern) section.', site_name='Continental Shelf East')}.items()}) module-attribute