Here is C++ example I came to, which can work on any signal, including nondeterministic not differentiable signal.
int main()
{
// user limits
const constexpr size_t min_period = 2;
const constexpr size_t max_period = 10;
assert(max_period >= min_period);
struct Corr
{
int period;
float corr;
};
// non normalized example with periodicity
float in_arr[] = { 3, 4, 8.9, 8, 3, 9, 7, 4, 9.3, 6, 8, 9.1, 2 };
const constexpr size_t in_arr_size = sizeof(in_arr) / sizeof(in_arr[0]);
assert(min_period < in_arr_size);
assert(max_period < in_arr_size);
Corr out_arr[in_arr_size - 1];
float in_arr_square[in_arr_size];
// calibration
auto num_offset_shifts = max_period + 1;
num_offset_shifts = (std::max)(num_offset_shifts, min_period + 1);
if (in_arr_size < min_period + num_offset_shifts) {
// recalculate
num_offset_shifts = in_arr_size - min_period;
}
const auto num_autocorr_values = min_period + num_offset_shifts;
for (size_t i = 0; i < in_arr_size; i++) {
in_arr_square[i] = in_arr[i] * in_arr[i];
}
for (size_t i = 0, offset_shift = min_period; max_period >= offset_shift && num_offset_shifts >= min_period; i++, offset_shift++, num_offset_shifts--) {
auto & autocorr = out_arr[i];
autocorr.period = offset_shift;
float corr_numerator_value = 0;
float corr_denominator_first_value = 0;
float corr_denominator_second_value = 0;
for (size_t j = 0; j < num_offset_shifts; j++) {
const float corr_value = in_arr[j] * in_arr[j + offset_shift];
corr_denominator_first_value += in_arr_square[j];
corr_denominator_second_value += in_arr_square[j + offset_shift];
corr_numerator_value += corr_value * corr_value;
}
autocorr.corr = std::sqrt(corr_numerator_value * num_autocorr_values * num_autocorr_values / (std::max)(corr_denominator_first_value, corr_denominator_second_value));
}
auto begin_it = &out_arr[0];
auto end_it = &out_arr[in_arr_size - 1];
std::sort(begin_it, end_it, [&](const Corr & l, const Corr & r) -> bool
{
return l.corr > r.corr;
});
return 0;
}
The sorted result:
+ [0] {period=3 corr=97.1897888 }
+ [1] {period=6 corr=93.0486755 }
+ [2] {period=9 corr=85.9366302 }
+ [3] {period=8 corr=85.2749481 }
+ [4] {period=5 corr=85.0467377 }
+ [5] {period=2 corr=82.3933029 }
+ [6] {period=4 corr=76.3764877 }
+ [7] {period=7 corr=74.6706696 }
+ [8] {period=10 corr=49.8527908 }
+ [9] {period=-858993460 corr=-107374176. }
+ [10]{period=-858993460 corr=-107374176. }
+ [11]{period=-858993460 corr=-107374176. }
It sorted by correlation value which related to highest periodicity.
It only works if several conditions are met:
- The input signal is clean enough without noise. If it isn't, then you have to filter it first. For example, if the input values in range [0; 1.0] and a highest value does repeat, then you can cutoff lowest values below some constant like
0.5
to 0.5
: 0, 0.3, 0.6, 0.9
-> 0.5, 0.5, 0.6, 0.9
.
- The input signal has enough dispersion in some extent. If it isn't, then you have to preprocess values to increase the dispersion.
Note: There could be more conditions to met to make it work.
Properties:
- The bigger values with a period has higher correlation results.
The formula:
autocorr.corr = std::sqrt(corr_numerator_value * num_autocorr_values * num_autocorr_values / (std::max)(corr_denominator_first_value, corr_denominator_second_value));
Based on functions multiplication: corr(x) = fa * fb / max(fa * fa, fb * fb)
square root - conversion back to linear
num_autocorr_values
- weight factor
Note: You still have to test all periods on quality. The correlation just arranges the period by some factor, but the best result still can be not with greater correlation value.