mirror of
https://gitlab.com/wgp/dougal/software.git
synced 2025-12-06 07:37:08 +00:00
81 lines
2.7 KiB
JavaScript
81 lines
2.7 KiB
JavaScript
const math = require('mathjs');
|
|
|
|
/**
|
|
* Compute PCA (eigenvectors and eigenvalues) for deviation data to assess geometric repeatability.
|
|
* @param {Array<Array<number>>} deviations - Array of [point, line, i_deviation, j_deviation]
|
|
* @returns {Object} - { eigenvalues, eigenvectors, rms, anisotropy, primaryDirection }
|
|
*/
|
|
function computePCA(deviations) {
|
|
// Extract i_deviation and j_deviation
|
|
const deviationMatrix = deviations.map(row => [row[2], row[3]]);
|
|
|
|
// Convert to mathjs matrix
|
|
const D = math.matrix(deviationMatrix);
|
|
|
|
// Compute mean for centering (1 x 2 matrix)
|
|
const mean = math.mean(D, 0);
|
|
|
|
// Explicitly repeat mean to match D's shape (n x 2)
|
|
const n = deviationMatrix.length;
|
|
const meanRepeated = math.repmat(mean, n, 1);
|
|
|
|
// Center the data
|
|
const centered = math.subtract(D, meanRepeated);
|
|
|
|
// Compute covariance matrix: (1/(n-1)) * (D_centered^T * D_centered)
|
|
const covMatrix = math.multiply(
|
|
math.multiply(1 / (n - 1), math.transpose(centered)),
|
|
centered
|
|
);
|
|
|
|
// Perform eigen decomposition
|
|
const result = math.eigs(covMatrix);
|
|
let eigenvalues = result.values;
|
|
const evObjs = result.eigenvectors;
|
|
|
|
// Convert eigenvalues to array if it's a matrix
|
|
eigenvalues = Array.isArray(eigenvalues) ? eigenvalues : eigenvalues.toArray();
|
|
|
|
// Create pairs and convert vector to array if necessary
|
|
const pairs = eigenvalues.map((val, i) => {
|
|
let vec = evObjs[i].vector;
|
|
if (vec.toArray) vec = vec.toArray();
|
|
return { val, vec };
|
|
});
|
|
|
|
// Sort by descending eigenvalues
|
|
pairs.sort((a, b) => b.val - a.val);
|
|
|
|
// Sorted eigenvalues
|
|
const sortedEigenvalues = pairs.map(p => p.val);
|
|
|
|
// Build eigenvector matrix: rows as components, columns as eigenvectors
|
|
const dimension = pairs[0].vec.length; // e.g., 2
|
|
const evecRows = [];
|
|
for (let comp = 0; comp < dimension; comp++) {
|
|
evecRows.push(pairs.map(p => p.vec[comp]));
|
|
}
|
|
const sortedEigenvectors = math.matrix(evecRows);
|
|
|
|
// Compute RMS errors along principal axes
|
|
const rms = sortedEigenvalues.map(val => Math.sqrt(Math.max(val, 0)));
|
|
|
|
// Compute anisotropy (ratio of major to minor axis variance)
|
|
const anisotropy = sortedEigenvalues[0] / (sortedEigenvalues[1] || 1); // Avoid division by zero
|
|
|
|
// Primary direction (angle in degrees of major eigenvector)
|
|
const primaryVector = sortedEigenvectors.subset(math.index([0, 1], 0)).toArray();
|
|
const primaryDirection = Math.atan2(primaryVector[1], primaryVector[0]) * 180 / Math.PI;
|
|
|
|
return {
|
|
eigenvalues: sortedEigenvalues,
|
|
eigenvectors: sortedEigenvectors.toArray(),
|
|
rms: rms, // RMS errors along major/minor axes
|
|
anisotropy: anisotropy, // Ratio of variances
|
|
primaryDirection: primaryDirection // Angle of major axis (degrees)
|
|
};
|
|
}
|
|
|
|
|
|
module.exports = { computePCA };
|