From b5093c94cef0c77aa62b85f0a2e0518e7e7a7c20 Mon Sep 17 00:00:00 2001 From: Leo Fang <leofang@bnl.gov> Date: Mon, 18 May 2020 16:10:39 -0400 Subject: [PATCH 1/6] fix NEP-18 --- autocorr/multitau.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/autocorr/multitau.py b/autocorr/multitau.py index 1b0bd1cd..17594edb 100644 --- a/autocorr/multitau.py +++ b/autocorr/multitau.py @@ -44,8 +44,8 @@ def even(x): # calculate levels levels = np.int_(np.log2(N / m)) + 1 dims = (a.shape[0], (levels + 1) * (m // 2)) - g2 = np.zeros(dims, dtype=np.float32) - tau = np.zeros(dims[1], dtype=np.float32) + g2 = np.zeros_like(a, shape=dims, dtype=np.float32) + tau = np.zeros_like(a, shape=(dims[1],), dtype=np.float32) # zero level delta_t = 1 From 6da3e063ba15c0b83fe9c5e7b1c06d3138d021f3 Mon Sep 17 00:00:00 2001 From: Leo Fang <leofang@bnl.gov> Date: Mon, 18 May 2020 23:41:22 -0400 Subject: [PATCH 2/6] Use recent NumPy to support NEP-18 --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index 163fdd06..1b4707e5 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,2 +1,2 @@ -numpy +numpy >=1.17 pybind11 From 3181e20c338bbb0a172438321d9e4694171a6336 Mon Sep 17 00:00:00 2001 From: tankonst <tatiana.v.konst@gmail.com> Date: Sun, 16 Aug 2020 12:12:31 -0400 Subject: [PATCH 3/6] adding function for calculating the entire ROI --- autocorr/multitau_roi.py | 78 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 78 insertions(+) create mode 100755 autocorr/multitau_roi.py diff --git a/autocorr/multitau_roi.py b/autocorr/multitau_roi.py new file mode 100755 index 00000000..d3318a95 --- /dev/null +++ b/autocorr/multitau_roi.py @@ -0,0 +1,78 @@ +import numpy as np + + +def multitau(signal, lags_per_level=16): + """Autocorrelation of a signal using multi-tau method. + + For details, please refer to D. Magatti and F. Ferri doi:10.1364/ao.40.004011 + + Parameters + ---------- + signal: 2-D array + input signal, autocorrlation is calulated along `1-axis` + lags_per_level: integer, optional + number of lag-times per level, 16 is a as good number as any other + + Returns: + -------- + autocorrelations: numpy.ndarray + should be self-explanatory + lag-times: numpy.ndarray + lag times in log2, corresponding to autocorrelation values + """ + # if N is add subtract 1 + def even(x): + return x if x % 2 == 0 else x - 1 + + # 1-D data hack + if signal.ndim == 1: + N = even(signal.shape[0]) + a = signal[np.newaxis, :N] + elif signal.ndim == 2: + # copy data a local array + N = even(signal.shape[1]) + a = signal[:, :N] + elif signal.ndim > 2: + raise ValueError('Flatten the [2,3,..] dimensions before passing to autocorrelate.') + + if N < lags_per_level: + raise ValueError('Lag times per level must be greater than length of signal.') + + # shorthand for long names + m = lags_per_level + + # calculate levels + levels = np.int_(np.log2(N / m)) + 1 + dims = (a.shape[0], (levels + 1) * (m // 2)) + g2 = np.zeros_like(a, shape=(dims[1],), dtype=np.float32) + tau = np.zeros_like(a, shape=(dims[1],), dtype=np.float32) + + # zero level + delta_t = 1 + for i in range(m): + tau[i] = i * delta_t + d = a.shape[0] + t1 = np.mean(a[:, :N - i]) + t2 = np.mean(a[:, i:]) + cov = np.diagonal(a[:, :-i].T @ a[:, i:]) + g2[i] = np.mean(cov) / t1 / t2/ d + a = (a[:, :N:2] + a[:, 1:N:2]) / 2 + N = even(N // 2) + + for level in range(1, levels): + delta_t *= 2 + for n in range(m // 2): + idx = m + (level - 1) * (m // 2) + n + shift = m // 2 + n + tau[idx] = tau[idx - 1] + delta_t + d = a.shape[0] + t1 = np.mean(a[:, :-shift]) + t2 = np.mean(a[:, shift:]) + cov = np.diagonal(a[:, :-shift].T @ a[:, shift:]) + g2[idx] = np.mean(cov) / t1 / t2/ d + a = (a[:, :N:2] + a[:, 1:N:2]) / 2 + N = even(N // 2) + if N < lags_per_level: + break + + return g2, tau From a59fca9cc36d7cd160a94691db41cc626c7be39f Mon Sep 17 00:00:00 2001 From: tankonst <tatiana.v.konst@gmail.com> Date: Sun, 16 Aug 2020 12:34:40 -0400 Subject: [PATCH 4/6] fixing __init__.py --- autocorr/__init__.py | 1 + 1 file changed, 1 insertion(+) diff --git a/autocorr/__init__.py b/autocorr/__init__.py index 36d6c5a7..4911847f 100644 --- a/autocorr/__init__.py +++ b/autocorr/__init__.py @@ -3,6 +3,7 @@ from .multitau import multitau +from .multitau_roi import multitau_roi from .fftautocorr import fftautocorr # noqa from .cAutocorr import fftautocorr as fftautocorr_mt # noqa try: From ce452918279e04cbce88ef7aadafd76c357761d4 Mon Sep 17 00:00:00 2001 From: tankonst <tatiana.v.konst@gmail.com> Date: Sun, 16 Aug 2020 12:46:51 -0400 Subject: [PATCH 5/6] fixing __init__.py and multitau_roi.py --- autocorr/__init__.py | 2 +- autocorr/multitau_roi.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/autocorr/__init__.py b/autocorr/__init__.py index 4911847f..745636c3 100644 --- a/autocorr/__init__.py +++ b/autocorr/__init__.py @@ -5,7 +5,7 @@ from .multitau import multitau from .multitau_roi import multitau_roi from .fftautocorr import fftautocorr # noqa -from .cAutocorr import fftautocorr as fftautocorr_mt # noqa +# from .cAutocorr import fftautocorr as fftautocorr_mt # noqa try: from .cAutocorr import multitau_mt except ImportError: diff --git a/autocorr/multitau_roi.py b/autocorr/multitau_roi.py index d3318a95..09dc332e 100755 --- a/autocorr/multitau_roi.py +++ b/autocorr/multitau_roi.py @@ -1,7 +1,7 @@ import numpy as np -def multitau(signal, lags_per_level=16): +def multitau_roi(signal, lags_per_level=16): """Autocorrelation of a signal using multi-tau method. For details, please refer to D. Magatti and F. Ferri doi:10.1364/ao.40.004011 From aac84eb5bb8eb637d8af106515c34bf94185390f Mon Sep 17 00:00:00 2001 From: tankonst <tatiana.v.konst@gmail.com> Date: Sun, 16 Aug 2020 14:38:28 -0400 Subject: [PATCH 6/6] Corrected zero-sized arrays warning --- autocorr/multitau_roi.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/autocorr/multitau_roi.py b/autocorr/multitau_roi.py index 09dc332e..4ab6245c 100755 --- a/autocorr/multitau_roi.py +++ b/autocorr/multitau_roi.py @@ -52,10 +52,11 @@ def even(x): for i in range(m): tau[i] = i * delta_t d = a.shape[0] - t1 = np.mean(a[:, :N - i]) + t1 = np.mean(a[:, : N - i]) t2 = np.mean(a[:, i:]) - cov = np.diagonal(a[:, :-i].T @ a[:, i:]) + cov = np.diagonal(a[:, :N-i].T @ a[:, i:]) g2[i] = np.mean(cov) / t1 / t2/ d + a = (a[:, :N:2] + a[:, 1:N:2]) / 2 N = even(N // 2) @@ -66,10 +67,11 @@ def even(x): shift = m // 2 + n tau[idx] = tau[idx - 1] + delta_t d = a.shape[0] - t1 = np.mean(a[:, :-shift]) + t1 = np.mean(a[:, : N -shift]) t2 = np.mean(a[:, shift:]) - cov = np.diagonal(a[:, :-shift].T @ a[:, shift:]) + cov = np.diagonal(a[:, : N-shift].T @ a[:, shift:]) g2[idx] = np.mean(cov) / t1 / t2/ d + a = (a[:, :N:2] + a[:, 1:N:2]) / 2 N = even(N // 2) if N < lags_per_level: