To ensure the accuracy of MVA, only when the ratio of the middle to the minimum eigenvalue (labeled QMVA for simplicity) is larger than 3 are the results used for further analysis.
The following method implicitly assumes the data is evenly sampled, otherwise, data resampling is needed.
This program computes the principal variance directions and variances of a vector quantity as well as the associated eigenvalues.*
Type
Details
data
ndarray
Returns
vrot:
an array of (npoints, ndim) containing the rotated data in the new coordinate system, ijk. Vi(maximum direction)=vrot[0,:] Vj(intermediate direction)=vrot[1,:] Vk(minimum variance direction)=Vrot[2,:]
Exported source
def minvar(data: np.ndarray):""" see `pyspedas.cotrans.minvar` This program computes the principal variance directions and variances of a vector quantity as well as the associated eigenvalues. Parameters ----------- data: Vxyz, an (npoints, ndim) array of data(ie Nx3) Returns ------- vrot: an array of (npoints, ndim) containing the rotated data in the new coordinate system, ijk. Vi(maximum direction)=vrot[0,:] Vj(intermediate direction)=vrot[1,:] Vk(minimum variance direction)=Vrot[2,:] v: an (ndim,ndim) array containing the principal axes vectors Maximum variance direction eigenvector, Vi=v[*,0] Intermediate variance direction, Vj=v[*,1] (descending order) w: the eigenvalues of the computation """# Min var starts here# data must be Nx3 vecavg = np.nanmean(np.nan_to_num(data, nan=0.0), axis=0) mvamat = np.zeros((3, 3))for i inrange(3):for j inrange(3): mvamat[i, j] = ( np.nanmean(np.nan_to_num(data[:, i] * data[:, j], nan=0.0))- vecavg[i] * vecavg[j] )# Calculate eigenvalues and eigenvectors w, v = np.linalg.eigh(mvamat, UPLO="U")# Sorting to ensure descending order w = np.abs(w) idx = np.flip(np.argsort(w))# IDL compatabilityifTrue:if np.sum(w) ==0.0: idx = [0, 2, 1] w = w[idx] v = v[:, idx]# Rotate intermediate var direction if system is not Right Handed YcrossZdotX = v[0, 0] * (v[1, 1] * v[2, 2] - v[2, 1] * v[1, 2])if YcrossZdotX <0: v[:, 1] =-v[:, 1]# v[:, 2] = -v[:, 2] # Should not it is being flipped at Z-axis?# Ensure minvar direction is along +Z (for FAC system)if v[2, 2] <0: v[:, 2] =-v[:, 2] v[:, 1] =-v[:, 1] vrot = np.array([np.dot(row, v) for row in data])return vrot, v, w
Note: - see datetime_to_numeric in xarray.core.duck_array_ops for more details about converting datetime to numeric - Xarray uses the numpy dtypes datetime64[ns] and timedelta64[ns] to represent datetime data.*
time = pd.date_range("2000-01-01", periods=10)x = np.linspace(0, np.pi, 10)# generate data circular in three dimensions, so the biggest distance is between the first and the last pointdata = np.zeros((10, 3))data[:, 0] = np.cos(x)data[:, 1] = np.sin(x)data = xr.DataArray(data, coords={"time": time}, dims=["time", "space"])fit_maxiumum_variance_direction(data[:, 0])["fit.best_fit"]