Qu

Signal processing

Multisine excitation

multisine synthesis with a per-tone figure breakdown.

catalog/qu_multisine.qu

The broadcast-synthesis idiom -- a frequency column against a time row -- with a wide log-spaced tone set snapped to DFT bins, plus a `figure()` loop that opens one figure per tone so each sine can be inspected on its own.

How it works

Setup

fo = logspace(-2, 3, 64)
Fs = 5e3
N = 1000
df = Fs / N
dt = 1 / Fs

f = round(fo / df) * df
f = unique(f)
K = length(f)
print("freqs: {f} -> df={df}")

tmeas = N * dt
t = 0 to (N - 1) * dt step dt

A = ones(K, 1)
phi = randn(K, 1)

f as vector(K, 1)
A as vector(K, 1)
phi as vector(K, 1)
t as vector(1, N)

Synthesize

Mt := t as matrix(K, N)
Mx := A .* cos(2*pi*f .* Mt + phi)
x = sum(Mx, axis=0)

Plot

plot(t, x)
xlabel("Time (s)")
ylabel("Amplitude")
grid on
box on
show plot

Per-tone figures

for k = 0 to K - 1
    figure()
    plot(t, Mx[k, :])
    title("Sine {k}")
end for

Results

printed
freqs: [0, 5, 10, 15, 20, 25, 30, 35, ... (26 elements)] -> df=5
-10 -5 0 5 10 15 0 0.02 0.04 0.06 0.08 0.10 0.12 0.14 0.16 0.18 0.20 Time (s) Amplitude
Drawn by running the program above.

The whole program

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

# qu_multisine.qu — multisine synthesis with a per-tone figure breakdown.

#

# The broadcast-synthesis idiom -- a frequency column against a time row --

# with a wide log-spaced tone set snapped to DFT bins, plus a `figure()` loop that opens one figure per

# tone so each sine can be inspected on its own.


backend auto

#%% Setup

fo = logspace(-2, 3, 64)
Fs = 5e3
N = 1000
df = Fs / N
dt = 1 / Fs

f = round(fo / df) * df
f = unique(f)
K = length(f)
print("freqs: {f} -> df={df}")

tmeas = N * dt
t = 0 to (N - 1) * dt step dt

A = ones(K, 1)
phi = randn(K, 1)

f as vector(K, 1)
A as vector(K, 1)
phi as vector(K, 1)
t as vector(1, N)

#%% Synthesize

Mt := t as matrix(K, N)
Mx := A .* cos(2*pi*f .* Mt + phi)
x = sum(Mx, axis=0)

#%% Plot

plot(t, x)
xlabel("Time (s)")
ylabel("Amplitude")
grid on
box on
show plot

#%% Per-tone figures

for k = 0 to K - 1
    figure()
    plot(t, Mx[k, :])
    title("Sine {k}")
end for