Created
July 11, 2018 01:15
-
-
Save qgustavor/e9c09635cbf8e64944936fedf75d2983 to your computer and use it in GitHub Desktop.
Time Series Synchronization (not working, I don't know why)
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
// Source: http://lexfridman.com/fast-cross-correlation-and-time-series-synchronization-in-python/ | |
// Using: https://github.com/nicolaspanel/numjs | |
function cross_correlation_using_fft (x, y) { | |
let f1 = nj.fft(x) | |
let f2 = nj.fft(nj.flip(y, 1)) | |
let cc = nj.ifft(f1.multiply(f2)).slice(0, [0, 1]) | |
return fftshift(cc) | |
} | |
// shift < 0 means that y starts 'shift' time steps before x | |
// shift > 0 means that y starts 'shift' time steps after x | |
function compute_shift (x, y) { | |
console.assert(x.shape[0] === y.shape[0], x.shape[0] + ' !== ' + y.shape[0]) | |
const c = cross_correlation_using_fft(x, y) | |
console.assert(c.shape[0] === x.shape[0], c.shape[0] + ' !== ' + x.shape[0]) | |
const zero_index = Math.floor(x.shape[0] / 2) - 1 | |
const shift = zero_index - c.tolist().indexOf(c.max()) | |
return shift | |
} | |
// fftshift from https://github.com/oramics/fftshift/blob/master/index.js | |
function rotate (src, n) { | |
var len = src.length | |
reverse(src, 0, len) | |
reverse(src, 0, n) | |
reverse(src, n, len) | |
return src | |
} | |
function reverse (src, from, to) { | |
--from | |
while (++from < --to) { | |
var tmp = src[from] | |
src[from] = src[to] | |
src[to] = tmp | |
} | |
} | |
function fftshift (src) { | |
// Convert ndArray to Array then back to ndArray again | |
const arr = src.tolist() | |
const len = arr.length | |
return nj.array(rotate(src, Math.floor(len / 2))).reshape(src.shape) | |
} | |
// Test function | |
compute_shift( | |
nj.concatenate(nj.arange(10).reshape(10, 1), nj.arange(100, 110, 1).reshape(10, 1)), | |
nj.concatenate(nj.arange(10).reshape(10, 1), nj.arange(200, 210, 1).reshape(10, 1)) | |
) | |
// Expected: 100 | |
// Returns: 5 |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment