statistics.py 950 Bytes
 Dion Haefner committed Aug 31, 2016 1 ``````import numpy as np `````` Dion Haefner committed Dec 14, 2016 2 ``````from scipy import integrate, signal `````` Dion Haefner committed Aug 31, 2016 3 `````` `````` Dion Haefner committed Oct 20, 2016 4 5 6 7 8 9 10 11 12 ``````def norm_L2(a,x): a_squared = np.squeeze(a**2) n = len(np.shape(a_squared)) for i in range(n): a_squared = integrate.simps(a_squared,x=x[i],axis=0) return np.sqrt(a_squared) def autocorrelation(H, pad=None): `````` Dion Haefner committed Aug 31, 2016 13 14 15 16 17 18 `````` """ Calculates the autocorrelation of an array via ffts. .. math:: \\text{autocovariance}(H,x) = \\text{fft}^{-1}(|\\text{fft}(H - \\overline{H},k)|^2,x) `````` Dion Haefner committed Oct 20, 2016 19 `````` `````` Dion Haefner committed Aug 31, 2016 20 21 22 `````` \\text{autocorrelation}(H,x) = \\text{autocovariance}(x) / \\text{autocovariance}(0) :param H: input array `````` Dion Haefner committed Oct 20, 2016 23 `````` :param pad: Pads each axis of the input with pad[i] zeros before and after each FFT. `````` Dion Haefner committed Aug 31, 2016 24 25 `````` """ `````` Dion Haefner committed Dec 14, 2016 26 27 28 29 `````` #pad = np.array(H.shape) if pad is None else pad #H = H - np.mean(H) #ft_signal = np.fft.fftn(H,pad) #B = np.fft.ifftn(np.abs(ft_signal)**2,pad) `````` Dion Haefner committed Dec 15, 2016 30 `````` reverse = [slice(n,0,-1) for n in H.shape] `````` Dion Haefner committed Dec 14, 2016 31 32 `````` autocovariance = signal.fftconvolve(H,H[reverse],mode="full")[reverse] return autocovariance / autocovariance.flat[0]``````