Levelling with equivalent sources

If we don’t have crossing lines, and therfore don’t have cross-over errors, we can instead use the data of nearby lines to level each line in a survey.

# %load_ext autoreload
# %autoreload 2

import logging

import cmocean
import matplotlib.pyplot as plt
import pandas as pd
import plotly.io as pio
import verde as vd

import airbornegeo

logging.getLogger("airbornegeo").setLevel("WARN")
logging.basicConfig()
pio.renderers.default = "notebook"

Load data

This is a subset of the BAS AGAP survey over Antarctica’s Gamburtsev Subglacial Mountains. The file is download and subset in the notebook AGAP_magnetic_survey, and the BAS processing steps are repeated in the notebook processing_AGAP_magnetic_survey.

data_df = pd.read_csv("data/AGAP_magnetic_survey_processed_blocked.csv")
data_df = data_df[
    [
        "easting",
        "northing",
        "height",
        "line",
        "distance_along_line",
        "mag",
    ]
]
data_df = data_df.rename(columns={"mag": "mag_unlevelled"})
data_df = data_df.dropna(subset=["mag_unlevelled"], how="any")

# pick a subset of lines
data_df = data_df[data_df.line >= 42]
data_df = data_df[data_df.line < 90]
data_df = data_df[data_df.line >= 75]
print(data_df.line.unique())
data_df.head()

[75 76 77 78 79 80 81 82 83 84 85 86 87 88 89]

	easting	northing	height	line	distance_along_line	mag_unlevelled
150229	1.101638e+06	204725.949316	4292.5	75	0.000000	-69.320
150230	1.101353e+06	204631.030460	4293.7	75	300.826802	-47.195
150231	1.101191e+06	204576.698238	4294.0	75	471.722344	-44.700
150232	1.100996e+06	204510.658613	4293.9	75	677.065902	-41.600
150233	1.100802e+06	204443.572754	4293.5	75	882.564878	-38.310

survey = airbornegeo.Survey(
    data_df,
    line_column="line",
    distance_column="distance_along_line",
    height_column="height",
)

ax = data_df[::10].plot.scatter(
    "easting",
    "northing",
    c="line",
    s=0.6,
    cmap="rainbow",
)
ax.set_aspect("equal")

max_abs = vd.maxabs(data_df.mag_unlevelled, percentile=95)

ax = data_df[::10].plot.scatter(
    "easting",
    "northing",
    c="mag_unlevelled",
    s=0.6,
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Unlevelled",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.6)
plt.show()

Perform iterative equivalent source levelling

This function will go through each line and perform the following steps:

1. find the data within max_dist of the line data (exluding the line data itself)

2. fit equivalent sources to this nearby data (with the provided damping, depth, and block_size parameters)

3. predict the forward effect of these equivalent sources onto the line

4. compute a misfit between the observed and predicted values

5. fit a trend ff the specified degree (0 for DC-shift, 1 for tilt, 2 for curve etc.)

6. add this trend to the observed line data.

This results in a levelled line. The algorithm then processed to the next line. The order it goes through the lines is randomized, as set by the seed parameter for reproducibility. by default, all lines with be levelled. If you only want to level a subset of lines, you can pass a list of line names to parameter lines_to_level. Once all (or subset of) lines are levelled, this concludes the first iteration. The entire processed can be repeated by setting parameter max_iterations to greater than 1.

The iterations will terminate on one of 4 stopping criteria:

1. the max_iterations is reached

2. the root mean square (RMS) of the levelling correction values for the iteration is below the set rms_tolerance

3. the RMS of the levelling correction values has increased more than the set rms_percent_increase_tolerance relative to the minimum RMS value for past iterations; this helps to stop run-away iterations

4. the RMS of the levelling correction values has not changed by more than the set rms_percent_change_tolerance between 2 subsequent iterations; this help save time by ending iterations when they aren’t offering much improvement.

survey.equivalent_source_levelling(
    data_column="mag_unlevelled",
    result_column="mag_levelled",
    max_dist=20e3,
    degree=2,
    damping=None,
    depth="default",
    data_block_size=1e3,
    source_block_size=10e3,
    max_iterations=10,
)
data_df = survey.data
data_df.head()

fig, axs = plt.subplots(2, 1, figsize=(10, 6))

max_abs = vd.maxabs(data_df.mag_unlevelled, percentile=95)

ax = data_df[::10].plot.scatter(
    "easting",
    "northing",
    c="mag_unlevelled",
    s=1,
    ax=axs[0],
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Unlevelled",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.8)

ax = data_df[::10].plot.scatter(
    "easting",
    "northing",
    c="mag_levelled",
    s=1,
    ax=axs[1],
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Levelled",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.8)
plt.tight_layout()
plt.show()