Qu

Signal processing

Interpolate, apply & cut

signal-aware interpolation, elementwise transforms, and range extraction (commit 48ea648).

catalog/qu_interpolate_apply_cut.qu

Three small, composable Signal operations: interpolate_at(sig, t_query, [method]) -- interp1 against the signal's OWN implicit time axis (sig.t = [0, 1/Fs, 2/Fs, ...]), so the caller never has to reconstruct that axis by hand. Returns a plain Num/Vec (never re-wrapped as a Signal -- arbitrary query times have no honest Fs to attach). sig.apply(fnName) -- runs a user-defined (:= or `function`) scalar function elementwise over a Signal's samples, preserving it AS a Signal (same Fs in, same Fs out) -- unlike interpolate_at, this one genuinely stays a Signal because the sample rate is unchanged. cut(sig, t_start, t_end) -- slices out the range [t_start, t_end], INCLUSIVE of both ends, by rounding to the nearest sample index at the signal's own Fs.

What it exercises

How it works

Build a real signal: 5 Hz sine, Fs=100 Hz, 0.5 s (50 samples)

Fs = 100
N  = 50
n  = 0 to N - 1
t_axis = n / Fs
values = sin(2 * pi * 5 * t_axis)
sig = signal(values, Fs)
print("sig: N={length(values)} samples, Fs={sig.Fs} Hz, duration={N / Fs} s")

interpolate_at vs interp1 on the signal's own time axis -- must agree

exactly (interpolate_at just reuses interp1_eval with sig.t as x).

queries = [0.055, 0.123, 0.301]
via_interpolate_at = interpolate_at(sig, queries)
via_interp1 = interp1(sig.t, sig, queries)
print("query times        = {queries}")
print("interpolate_at(sig)= {via_interpolate_at}")
print("interp1(sig.t, sig)= {via_interp1}")
print("max abs difference = {max(abs(via_interpolate_at - via_interp1)):.2e} (expect 0 -- same algorithm)")

# Sanity: value at t=0.055 sits between samples 5 (t=0.05) and 6 (t=0.06).
y5 = sin(2 * pi * 5 * 0.05)
y6 = sin(2 * pi * 5 * 0.06)
hand_linear = y5 + (0.055 - 0.05) / (0.06 - 0.05) * (y6 - y5)
print("hand-computed linear interp at t=0.055 = {hand_linear:.6f} vs interpolate_at = {via_interpolate_at[0]:.6f}")

sig.apply("double") -- elementwise transform, still a Signal, same Fs

double(x) := x * 2
doubled = sig.apply("double")
print("doubled[0:4]  = {doubled[0:4]}")
print("original[0:4] = {sig[0:4]}")
print("doubled.Fs == sig.Fs? {doubled.Fs == sig.Fs} (Fs preserved: {doubled.Fs})")
print("doubled really is 2x: max abs diff = {max(abs(doubled - 2 * sig)):.2e}")

cut(sig, t_start, t_end) -- extract a range, verify length/duration

segment = cut(sig, 0.1, 0.3)
# Both ends are kept, so the count is the span in samples PLUS ONE -- the
# fencepost. A window from t=0.1 to t=0.3 at 100 Hz holds samples 10
# through 30, which is 21 of them, not 20.
expected_len = round((0.3 - 0.1) * Fs) + 1
print("cut(sig, 0.1, 0.3): length={length(segment)} (expected {expected_len}), Fs={segment.Fs}")
# ... and the duration spanned is one sample interval longer than the
# difference of the two times, for the same reason.
expected_dur = (0.3 - 0.1) + 1 / Fs
print("segment duration = {length(segment) / segment.Fs} s (expected {expected_dur} s)")
print("segment[0] should equal original sample at t=0.1 (index 10): segment[0]={segment[0]:.6f}, sig[10]={sig[10]:.6f}")

Results

printed
sig: N=50 samples, Fs=100 Hz, duration=0.5 s
query times        = [0.055, 0.123, 0.301]
interpolate_at(sig)= [0.975528, -0.654155, -0.030902]
interp1(sig.t, sig)= [0.975528, -0.654155, -0.030902]
max abs difference = 0.00e0 (expect 0 -- same algorithm)
hand-computed linear interp at t=0.055 = 0.975528 vs interpolate_at = 0.975528
doubled[0:4]  = signal(Fs=100) [0, 0.618034, 1.175571, 1.618034, 1.902113]
original[0:4] = signal(Fs=100) [0, 0.309017, 0.587785, 0.809017, 0.951057]
doubled.Fs == sig.Fs? true (Fs preserved: 100)
doubled really is 2x: max abs diff = 0.00e0
cut(sig, 0.1, 0.3): length=21 (expected 21), Fs=100
segment duration = 0.21 s (expected 0.21 s)
segment[0] should equal original sample at t=0.1 (index 10): segment[0]=0.000000, sig[10]=0.000000

The whole program

Runnable as it stands: qu run catalog/qu_interpolate_apply_cut.qu

# qu_interpolate_apply_cut.qu -- signal-aware interpolation, elementwise

# transforms, and range extraction (commit 48ea648).

#

# Three small, composable Signal operations:

#   interpolate_at(sig, t_query, [method]) -- interp1 against the signal's

#     OWN implicit time axis (sig.t = [0, 1/Fs, 2/Fs, ...]), so the caller

#     never has to reconstruct that axis by hand. Returns a plain

#     Num/Vec (never re-wrapped as a Signal -- arbitrary query times have

#     no honest Fs to attach).

#   sig.apply(fnName) -- runs a user-defined (:= or `function`) scalar

#     function elementwise over a Signal's samples, preserving it AS a

#     Signal (same Fs in, same Fs out) -- unlike interpolate_at, this one

#     genuinely stays a Signal because the sample rate is unchanged.

#   cut(sig, t_start, t_end) -- slices out the range [t_start, t_end],

#     INCLUSIVE of both ends, by rounding to the nearest sample index at

#     the signal's own Fs.

#

# Constructs exercised (all real, verified-working builtins):

#   signal(data, Fs) / sig.Fs / sig.t     -- Signal construction + fields

#   interp1(x, y, xq)                     -- reference for interpolate_at

#   interpolate_at(sig, t_query, [method])

#   double(x) := x * 2  /  sig.apply("double")

#   cut(sig, t_start, t_end)


backend auto

#%% Build a real signal: 5 Hz sine, Fs=100 Hz, 0.5 s (50 samples)

Fs = 100
N  = 50
n  = 0 to N - 1
t_axis = n / Fs
values = sin(2 * pi * 5 * t_axis)
sig = signal(values, Fs)
print("sig: N={length(values)} samples, Fs={sig.Fs} Hz, duration={N / Fs} s")

#%% interpolate_at vs interp1 on the signal's own time axis -- must agree

# exactly (interpolate_at just reuses interp1_eval with sig.t as x).

queries = [0.055, 0.123, 0.301]
via_interpolate_at = interpolate_at(sig, queries)
via_interp1 = interp1(sig.t, sig, queries)
print("query times        = {queries}")
print("interpolate_at(sig)= {via_interpolate_at}")
print("interp1(sig.t, sig)= {via_interp1}")
print("max abs difference = {max(abs(via_interpolate_at - via_interp1)):.2e} (expect 0 -- same algorithm)")

# Sanity: value at t=0.055 sits between samples 5 (t=0.05) and 6 (t=0.06).

y5 = sin(2 * pi * 5 * 0.05)
y6 = sin(2 * pi * 5 * 0.06)
hand_linear = y5 + (0.055 - 0.05) / (0.06 - 0.05) * (y6 - y5)
print("hand-computed linear interp at t=0.055 = {hand_linear:.6f} vs interpolate_at = {via_interpolate_at[0]:.6f}")

#%% sig.apply("double") -- elementwise transform, still a Signal, same Fs

double(x) := x * 2
doubled = sig.apply("double")
print("doubled[0:4]  = {doubled[0:4]}")
print("original[0:4] = {sig[0:4]}")
print("doubled.Fs == sig.Fs? {doubled.Fs == sig.Fs} (Fs preserved: {doubled.Fs})")
print("doubled really is 2x: max abs diff = {max(abs(doubled - 2 * sig)):.2e}")

#%% cut(sig, t_start, t_end) -- extract a range, verify length/duration

segment = cut(sig, 0.1, 0.3)
# Both ends are kept, so the count is the span in samples PLUS ONE -- the

# fencepost. A window from t=0.1 to t=0.3 at 100 Hz holds samples 10

# through 30, which is 21 of them, not 20.

expected_len = round((0.3 - 0.1) * Fs) + 1
print("cut(sig, 0.1, 0.3): length={length(segment)} (expected {expected_len}), Fs={segment.Fs}")
# ... and the duration spanned is one sample interval longer than the

# difference of the two times, for the same reason.

expected_dur = (0.3 - 0.1) + 1 / Fs
print("segment duration = {length(segment) / segment.Fs} s (expected {expected_dur} s)")
print("segment[0] should equal original sample at t=0.1 (index 10): segment[0]={segment[0]:.6f}, sig[10]={sig[10]:.6f}")