Starting point are real world measurements containing one unique sinusoidal signal.
For signal processing purposes it might be beneficial to reduce the offset around zero before evaluating the signal either with fft() or with a parameter fit.
This was the reason to write a function, that eliminates / reduces the offset point by point, by means of a moving window.
The function works fine, but the time needed for processing is quite high, if the sampling rate goes up.
Now I wonder, if a similar function already exists that does the job in a more efficient manner, e.g. by spreading tasks on several cores.
In case others are looking for something similar or have suggestions, I put my example code below:
using PlotlyJS, Printf
using RobustModels
using Debugger
b_plt = false
b_trunc_signal = false
num_periods = 5
sampling_rate = 100000
frequ1 = .1
ampl1 = 2.1
phase1 = 10.0 # in degree
offset_sp = 0.1
grad_offset = 0.0001 # gradient of the offset
# ---
ppp_frequ1 = trunc(Int, sampling_rate / frequ1)
window_length = trunc(Int, num_periods * ppp_frequ1)
# --- time and data vector:
vec_t = collect(range(0, step = 1.0 / sampling_rate, length = window_length))
td_data_ = offset_sp .+ (grad_offset * ampl1) .* vec_t .+ ( ampl1 .* cos.( (2 * pi * frequ1 .* vec_t) .+ deg2rad(phase1)) )
# --- Functions ---------------------------------------------------------------------------------------------
function _MyLibCutOffset(_signal::Vector{<:Number}, _sizeMV::Int, _b_truncSignal::Bool=true, _b_verbose::Bool=true)
b_test_without_mean_command = true
b_DBG = true
_signal = Float64.(_signal)
_n_pts_signal = length(_signal)
if _n_pts_signal <= _sizeMV
error("Signal too short!")
end
# max_periods = trunc(Int, _n_pts_signal / _sizeMV)
pts_MV_correted = _n_pts_signal - _sizeMV
if _b_truncSignal
signal_symetric = Vector{Float64}(undef, pts_MV_correted)
idx_shift = 0
else
signal_symetric = Vector{Float64}(undef, _n_pts_signal)
if b_test_without_mean_command
i_mean = 0.0
else
i_mean = RobustModels.mean(_signal[1:_sizeMV])
end
signal_symetric[1:floor(Int,_sizeMV/2)] = _signal[1:floor(Int,_sizeMV/2)] .- i_mean
idx_shift = floor(Int,_sizeMV/2)
end
# ---
if _b_verbose
@info(string("MyLibCutOffset || _b_truncSignal: ", _b_truncSignal, ",\t_n_pts_signal: ", _n_pts_signal, ",\t_sizeMV: ", _sizeMV, ",\tpts_MV_correted: ", pts_MV_correted, ",\tidx_shift: ", idx_shift))
end
for i_pts = 1:pts_MV_correted
_range = range(i_pts, length = _sizeMV)
try
if b_test_without_mean_command
i_mean = 0.0
else
i_mean = RobustModels.mean(_signal[_range])
end
catch e
println(string("i_pts: ", i_pts,", \ti_pts + idx_shift: ", i_pts + idx_shift, ", _range[end]: ", _range[end], ", length(_range): ", length(_range)))
end
signal_symetric[i_pts + idx_shift] = _signal[i_pts + floor(Int,_sizeMV/2)] - i_mean
end
# --- DBG ---
if b_DBG
println(string("pts_MV_correted - ceil(Int,_sizeMV/2): ", pts_MV_correted - ceil(Int,_sizeMV/2)))
end
if ~_b_truncSignal
diff_last_corr = signal_symetric[pts_MV_correted + idx_shift] - _signal[pts_MV_correted + floor(Int,_sizeMV/2)]
signal_symetric[pts_MV_correted+idx_shift+1:end] = _signal[pts_MV_correted+idx_shift+1:end] .- (i_mean + diff_last_corr)
if b_DBG
println(string("pts_MV_correted - ceil(Int,_sizeMV/2)+1: ", pts_MV_correted - ceil(Int,_sizeMV/2)+1,
", idx_shift+floor(Int,_sizeMV/2): ", idx_shift+floor(Int,_sizeMV/2)))
println(string("diff_last_corr: ", diff_last_corr))
end
end
if b_test_without_mean_command
if _b_truncSignal
idx_start = ceil(Int, pts_MV_correted/2)
idx_end = _n_pts_signal - floor(Int, pts_MV_correted/2)
signal_symetric = _signal[idx_start : idx_end]
else
signal_symetric = _signal
end
end
# ---
return signal_symetric
end
# --- Main --------------------------------------------------------------------------------------------------
signal_symetric_ = _MyLibCutOffset(td_data_, ppp_frequ1, b_trunc_signal)
n_pts_signal = length(signal_symetric_)
vec_t_symetric = vec_t[floor(Int, (window_length - n_pts_signal)/2) + 1 : window_length-ceil(Int, (window_length - n_pts_signal)/2)]
# --- info --------------------------------------------------------------------------------------------------
@info(string("window_length: ", window_length, ", n_pts_signal: ", n_pts_signal, ", vec_t_symetric: ", length(vec_t_symetric)))
@info(string("Mean input: ", RobustModels.mean(td_data_), ", Mean output: ", RobustModels.mean(signal_symetric_)))
# --- Plot --------------------------------------------------------------------------------------------------
if b_plt
if window_length > 1e6
error("Too many points!")
end
PlotlyJS.PlotlyJSDisplay()
hdl1 = PlotlyJS.scatter(; x = vec_t, y = td_data_, name = "raw")
hdl2 = PlotlyJS.scatter(; x = vec_t_symetric, y = signal_symetric_, name = "corr")
hdl_plt = PlotlyJS.Plot([hdl1, hdl2])
end
I am not sure about the option of a lowpass-filter, I have to check this idea.