AlgebraicShapeCombination.js (6322B)
1 import * as XYShape from "../XYShape.js"; 2 import { convolutionOperationToFn } from "./PointSet.js"; 3 import * as E_A_Floats from "../utility/E_A_Floats.js"; 4 const toDiscretePointMassesFromTriangulars = (s, { inverse } = { inverse: false }) => { 5 let n = XYShape.T.length(s); 6 const { xs, ys } = s; 7 xs.unshift(xs[0]); 8 ys.unshift(ys[0]); 9 xs.push(xs[n - 1]); 10 ys.push(ys[n - 1]); 11 n = xs.length; 12 const xsSq = new Array(n); 13 const xsProdN1 = new Array(n - 1); 14 const xsProdN2 = new Array(n - 2); 15 for (let i = 0; i <= n - 1; i++) { 16 xsSq[i] = xs[i] * xs[i]; 17 } 18 for (let i = 0; i <= n - 2; i++) { 19 xsProdN1[i] = xs[i] * xs[i + 1]; 20 } 21 for (let i = 0; i <= n - 3; i++) { 22 xsProdN2[i] = xs[i] * xs[i + 2]; 23 } 24 const masses = new Array(n - 2); 25 const means = new Array(n - 2); 26 const variances = new Array(n - 2); 27 if (inverse) { 28 for (let i = 1; i <= n - 2; i++) { 29 masses[i - 1] = ((xs[i + 1] - xs[i - 1]) * ys[i]) / 2; 30 const a = xs[i - 1]; 31 const c = xs[i]; 32 const b = xs[i + 1]; 33 const inverseMean = (2 * 34 ((a * Math.log(a / c)) / (a - c) + (b * Math.log(c / b)) / (b - c))) / 35 (a - b); 36 const inverseVar = (2 * (Math.log(c / a) / (a - c) + (b * Math.log(b / c)) / (b - c))) / 37 (a - b) - 38 inverseMean ** 2; 39 means[i - 1] = inverseMean; 40 variances[i - 1] = inverseVar; 41 } 42 return { n: n - 2, masses, means, variances }; 43 } 44 else { 45 for (let i = 1; i <= n - 2; i++) { 46 masses[i - 1] = ((xs[i + 1] - xs[i - 1]) * ys[i]) / 2; 47 means[i - 1] = (xs[i - 1] + xs[i] + xs[i + 1]) / 3; 48 variances[i - 1] = 49 (xsSq[i - 1] + 50 xsSq[i] + 51 xsSq[i + 1] - 52 xsProdN1[i - 1] - 53 xsProdN1[i] - 54 xsProdN2[i - 1]) / 55 18; 56 } 57 return { n: n - 2, masses, means, variances }; 58 } 59 }; 60 export const combineShapesContinuousContinuous = (op, s1, s2) => { 61 const t1m = toDiscretePointMassesFromTriangulars(s1); 62 const t2m = toDiscretePointMassesFromTriangulars(s2, { inverse: false }); 63 const combineMeansFn = { 64 Add: (m1, m2) => m1 + m2, 65 Subtract: (m1, m2) => m1 - m2, 66 Multiply: (m1, m2) => m1 * m2, 67 }[op]; 68 const combineVariancesFn = { 69 Add: (v1, v2) => v1 + v2, 70 Subtract: (v1, v2) => v1 + v2, 71 Multiply: (v1, v2, m1, m2) => v1 * v2 + v1 * m2 ** 2 + v2 * m1 ** 2, 72 }[op]; 73 let outputMinX = Infinity; 74 let outputMaxX = -Infinity; 75 const masses = new Array(t1m.n * t2m.n); 76 const means = new Array(t1m.n * t2m.n); 77 const variances = new Array(t1m.n * t2m.n); 78 for (let i = 0; i < t1m.n; i++) { 79 for (let j = 0; j < t2m.n; j++) { 80 const k = i * t2m.n + j; 81 masses[k] = t1m.masses[i] * t2m.masses[j]; 82 const mean = combineMeansFn(t1m.means[i], t2m.means[j]); 83 const variance = combineVariancesFn(t1m.variances[i], t2m.variances[j], t1m.means[i], t2m.means[j]); 84 means[k] = mean; 85 variances[k] = variance; 86 const minX = mean - 2 * Math.sqrt(variance) * 1.644854; 87 const maxX = mean + 2 * Math.sqrt(variance) * 1.644854; 88 if (minX < outputMinX) { 89 outputMinX = minX; 90 } 91 if (maxX > outputMaxX) { 92 outputMaxX = maxX; 93 } 94 } 95 } 96 const nOut = 300; 97 const outputXs = E_A_Floats.range(outputMinX, outputMaxX, nOut); 98 const outputYs = new Array(nOut).fill(0); 99 for (let j = 0; j < masses.length; j++) { 100 if (variances[j] > 0 && 101 masses[j] > 0) { 102 for (let i = 0; i < outputXs.length; i++) { 103 const dx = outputXs[i] - means[j]; 104 const contribution = (masses[j] * Math.exp(-(dx ** 2) / (2 * variances[j]))) / 105 Math.sqrt(2 * 3.14159276 * variances[j]); 106 outputYs[i] = outputYs[i] + contribution; 107 } 108 } 109 } 110 return { xs: outputXs, ys: outputYs }; 111 }; 112 export const combineShapesContinuousDiscrete = (op, continuousShape, discreteShape, opts) => { 113 const t1n = XYShape.T.length(continuousShape); 114 const t2n = XYShape.T.length(discreteShape); 115 const opFunc = convolutionOperationToFn(op); 116 const fn = opts.discretePosition === "First" 117 ? (a, b) => opFunc(b, a) 118 : opFunc; 119 const outXYShapes = new Array(t2n); 120 switch (op) { 121 case "Add": 122 case "Subtract": 123 for (let j = 0; j <= t2n - 1; j++) { 124 const dxyShape = new Array(t1n); 125 for (let i = 0; i <= t1n - 1; i++) { 126 const index = opts.discretePosition === "First" && op === "Subtract" 127 ? t1n - 1 - i 128 : i; 129 dxyShape[index] = [ 130 fn(continuousShape.xs[i], discreteShape.xs[j]), 131 continuousShape.ys[i] * discreteShape.ys[j], 132 ]; 133 } 134 outXYShapes[j] = dxyShape; 135 } 136 break; 137 case "Multiply": 138 for (let j = 0; j <= t2n - 1; j++) { 139 const dxyShape = new Array(t1n); 140 for (let i = 0; i <= t1n - 1; i++) { 141 const index = discreteShape.xs[j] > 0 ? i : t1n - 1 - i; 142 dxyShape[index] = [ 143 fn(continuousShape.xs[i], discreteShape.xs[j]), 144 (continuousShape.ys[i] * discreteShape.ys[j]) / 145 Math.abs(discreteShape.xs[j]), 146 ]; 147 } 148 outXYShapes[j] = dxyShape; 149 } 150 break; 151 default: 152 throw new Error(`Unknown operation ${op}`); 153 } 154 return outXYShapes 155 .map(XYShape.T.fromZippedArray) 156 .reduce((acc, x) => XYShape.PointwiseCombination.addCombine(XYShape.XtoY.continuousInterpolator("Linear", "UseZero"), acc, x), XYShape.T.empty); 157 }; 158 export const isOrdered = (a) => { 159 return E_A_Floats.isSorted(a.xs); 160 }; 161 //# sourceMappingURL=AlgebraicShapeCombination.js.map