Change samping frequency

Sometimes we require resample survey data at either a heigher or lower sampling frequency. Here we show how to do both.

# %load_ext autoreload
# %autoreload 2

import cmocean
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import verde as vd

import airbornegeo

data_df = pd.read_csv("data/AGAP_gravity_survey_processed.csv")
data_df = data_df[
    [
        "easting",
        "northing",
        "line",
        "unixtime",
        "distance_along_line",
        "grav_disturbance_filt",
    ]
][::5]  # every 5th point
data_df.head()

survey = airbornegeo.Survey(data_df, line_column="line")
data_df = survey.data

Reduce sampling frequency

For this we will use block-reduction, where we define spatial windows (1D), or blocks (2D), and for all data within the block, retain only the mean or median values.

# extract a single line from the survey
line_df = data_df[data_df.line == 8]
line_df

	easting	northing	line	unixtime	distance_along_line	grav_disturbance_filt
21815	1.064934e+06	327364.065720	8	1.230574e+09	202.954030	56.97
21820	1.065268e+06	327411.697676	8	1.230574e+09	540.476913	55.86
21825	1.065601e+06	327464.248325	8	1.230574e+09	878.033672	54.83
21830	1.065934e+06	327524.876461	8	1.230574e+09	1216.298627	53.83
21835	1.066267e+06	327589.259632	8	1.230574e+09	1555.357945	52.86
...	...	...	...	...	...	...
25255	1.321821e+06	373028.825792	8	1.230578e+09	261161.216753	-0.68
25260	1.322151e+06	373081.094428	8	1.230578e+09	261495.865268	-0.63
25265	1.322480e+06	373141.185876	8	1.230578e+09	261829.523736	-0.47
25270	1.322807e+06	373203.581372	8	1.230578e+09	262162.742950	-0.22
25275	1.323135e+06	373265.681016	8	1.230578e+09	262496.663302	0.07

693 rows × 6 columns

Block-reduce by distance

By setting reduce_by to ‘distance_along_line’ and spacing to 2000, we can block reduce the data to have 1 point every 2000 meters.

line_survey = airbornegeo.Survey(line_df, line_column="line")
blocked_line = line_survey.block_reduce(
    np.median,
    spacing=2000,
    reduce_by="distance_along_line",
).data

ax = line_df.plot.line(
    "distance_along_line",
    "grav_disturbance_filt",
    style="bp",
    ms=1,
    label="Original data",
)
ax = blocked_line.plot.line(
    "distance_along_line",
    "grav_disturbance_filt",
    style="rx",
    ms=3,
    ax=ax,
    label="Block-reduced data",
)
ax.set_ylabel("Gravity disturbance (mGal)")

Text(0, 0.5, 'Gravity disturbance (mGal)')

Block-reduce by time

By setting reduce_by to ‘unixtime’ and spacing to 60, we can block reduce the data to have 1 point every minute along the flight.

line_survey = airbornegeo.Survey(line_df, line_column="line")
blocked_line = line_survey.block_reduce(
    np.median,
    spacing=60,
    reduce_by="unixtime",
).data

ax = line_df.plot.line(
    "unixtime", "grav_disturbance_filt", style="bp", ms=1, label="Original data"
)
ax = blocked_line.plot.line(
    "unixtime",
    "grav_disturbance_filt",
    style="rx",
    ms=3,
    ax=ax,
    label="Block-reduced data",
)
ax.set_ylabel("Gravity disturbance (mGal)")

Text(0, 0.5, 'Gravity disturbance (mGal)')

Block-reduce all lines in a survey

By supplying ‘line’ to the groupby_column, the block-reduce occurs only on 1 line at a time. This means for lines that are closer together than the spacing, or where lines cross, the values from other lines within the same block are not included in the block-reduction.

blocked_survey = survey.block_reduce(
    np.median,
    spacing=5000,
    reduce_by="distance_along_line",
)
blocked_survey.data

max_abs = vd.maxabs(survey.data.grav_disturbance_filt, percentile=95)

ax = survey.data.plot.scatter(
    "easting",
    "northing",
    c="grav_disturbance_filt",
    s=0.1,
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Original data",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.6)
plt.show()

ax = blocked_survey.data.plot.scatter(
    "easting",
    "northing",
    c="grav_disturbance_filt",
    s=0.1,
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Block-reduced data",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.6)
plt.show()

Increase or change sampling frequency

Below we show how to increase the sampling frequency, or resample the data at a specified frequency. This can be done based on any column, but typically you would resample based on a distance along line column, or a time column.

As we can see below, the line data has a value every ~25 seconds.

# extract a single line from the survey
line_df = data_df[data_df.line == 8]
# for demonstration, only retain every 5th point
line_df = line_df[::5]

line_df.head()

	easting	northing	line	unixtime	distance_along_line	grav_disturbance_filt
21815	1.064934e+06	327364.065720	8	1.230574e+09	202.954030	56.97
21840	1.066602e+06	327650.536827	8	1.230575e+09	1895.906697	51.91
21865	1.068273e+06	327963.921140	8	1.230575e+09	3596.305820	48.32
21890	1.069939e+06	328274.457926	8	1.230575e+09	5291.379412	47.31
21915	1.071605e+06	328613.652292	8	1.230575e+09	6991.055147	47.73

line_survey = airbornegeo.Survey(line_df, line_column="line")
resampled_line = line_survey.resample(
    spacing=1,
    resample_by="unixtime",
    maxdist=10,  # only retain data within 10 seconds of original data
).data

ax = resampled_line.plot.line(
    "unixtime", "grav_disturbance_filt", style="bp", ms=0.6, label="Resampled data"
)
ax = line_df.plot.line(
    "unixtime", "grav_disturbance_filt", style="rx", ms=3, ax=ax, label="Original data"
)
ax.set_ylabel("Gravity disturbance (mGal)")

Text(0, 0.5, 'Gravity disturbance (mGal)')

Resample all lines in a survey

# for demonstration, only retain every 20th point
survey.data = survey.data[::20]
data_df = survey.data
data_df

resampled_survey = survey.resample(
    spacing=1,
    resample_by="unixtime",
    maxdist=60,  # only retain data within 60 seconds of original data
)
resampled_survey.data

max_abs = vd.maxabs(survey.data.grav_disturbance_filt, percentile=95)

ax = survey.data.plot.scatter(
    "easting",
    "northing",
    c="grav_disturbance_filt",
    s=0.1,
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    colorbar=False,
    title="Original data",
)
ax.set_aspect("equal")
plt.colorbar(ax.collections[0], ax=ax, shrink=0.6)
plt.show()

ax = resampled_survey.data.plot.scatter(
    "easting",
    "northing",
    c="grav_disturbance_filt",
    s=0.1,
    cmap=cmocean.cm.balance,
    vmin=-max_abs,
    vmax=max_abs,
    title="Resampled data",
)
ax.set_aspect("equal")
plt.show()