diff --git a/src/frechetDistance.ts b/src/frechetDistance.ts index 6c08690..6684b82 100644 --- a/src/frechetDistance.ts +++ b/src/frechetDistance.ts @@ -3,48 +3,53 @@ import { Curve, pointDistance } from './geometry'; /** * Discrete Frechet distance between 2 curves * based on http://www.kr.tuwien.ac.at/staff/eiter/et-archive/cdtr9464.pdf + * modified to be iterative and have better memory usage * @param curve1 * @param curve2 */ const frechetDist = (curve1: Curve, curve2: Curve) => { - const results: number[][] = []; - for (let i = 0; i < curve1.length; i++) { - results.push([]); - for (let j = 0; j < curve2.length; j++) { - results[i].push(-1); - } - } - - const recursiveCalc = (i: number, j: number) => { - if (results[i][j] > -1) return results[i][j]; + const longCurve = curve1.length >= curve2.length ? curve1 : curve2; + const shortCurve = curve1.length >= curve2.length ? curve2 : curve1; + const calcVal = ( + i: number, + j: number, + prevResultsCol: number[], + curResultsCol: number[] + ): number => { if (i === 0 && j === 0) { - results[i][j] = pointDistance(curve1[0], curve2[0]); - } else if (i > 0 && j === 0) { - results[i][j] = Math.max( - recursiveCalc(i - 1, 0), - pointDistance(curve1[i], curve2[0]) - ); - } else if (i === 0 && j > 0) { - results[i][j] = Math.max( - recursiveCalc(0, j - 1), - pointDistance(curve1[0], curve2[j]) - ); - } else if (i > 0 && j > 0) { - results[i][j] = Math.max( - Math.min( - recursiveCalc(i - 1, j), - recursiveCalc(i - 1, j - 1), - recursiveCalc(i, j - 1) - ), - pointDistance(curve1[i], curve2[j]) + return pointDistance(longCurve[0], shortCurve[0]); + } + if (i > 0 && j === 0) { + return Math.max( + prevResultsCol[0], + pointDistance(longCurve[i], shortCurve[0]) ); - } else { - results[i][j] = Infinity; } - return results[i][j]; + const lastResult = curResultsCol[curResultsCol.length - 1]; + if (i === 0 && j > 0) { + return Math.max(lastResult, pointDistance(longCurve[0], shortCurve[j])); + } + + return Math.max( + Math.min(prevResultsCol[j], prevResultsCol[j - 1], lastResult), + pointDistance(longCurve[i], shortCurve[j]) + ); }; - return recursiveCalc(curve1.length - 1, curve2.length - 1); + let prevResultsCol: number[] = []; + for (let i = 0; i < longCurve.length; i++) { + const curResultsCol: number[] = []; + for (let j = 0; j < shortCurve.length; j++) { + // we only need the results from i - 1 and j - 1 to continue the calculation + // so we only need to hold onto the last column of calculated results + // prevResultsCol is results[i-1][:] in the original algorithm + // curResultsCol is results[i][:j-1] in the original algorithm + curResultsCol.push(calcVal(i, j, prevResultsCol, curResultsCol)); + } + prevResultsCol = curResultsCol; + } + + return prevResultsCol[shortCurve.length - 1]; }; export default frechetDist; diff --git a/src/geometry.ts b/src/geometry.ts index 8392b1a..9431b91 100644 --- a/src/geometry.ts +++ b/src/geometry.ts @@ -13,8 +13,7 @@ export const subtract = (v1: Point, v2: Point): Point => ({ y: v1.y - v2.y }); -const magnitude = (vector: Point) => - Math.sqrt(Math.pow(vector.x, 2) + Math.pow(vector.y, 2)); +const magnitude = ({ x, y }: Point) => Math.sqrt(x * x + y * y); /** * Calculate the distance between 2 points diff --git a/test/frechetDistance.test.ts b/test/frechetDistance.test.ts index 6c5d943..2cfae67 100644 --- a/test/frechetDistance.test.ts +++ b/test/frechetDistance.test.ts @@ -1,5 +1,5 @@ import frechetDistance from '../src/frechetDistance'; -import { subdivideCurve } from '../src/geometry'; +import { rebalanceCurve, subdivideCurve } from '../src/geometry'; describe('frechetDist', () => { it('is 0 if the curves are the same', () => { @@ -41,4 +41,65 @@ describe('frechetDist', () => { expect(frechetDistance(curve1, curve2)).toBe(1); expect(frechetDistance(curve2, curve1)).toBe(1); }); + + it('gives correct results 1', () => { + const curve1 = [ + { x: 1, y: 0 }, + { x: 2.4, y: 43 }, + { x: -1, y: 4.3 }, + { x: 4, y: 4 } + ]; + const curve2 = [{ x: 0, y: 0 }, { x: 14, y: 2.4 }, { x: 4, y: 4 }]; + + expect(frechetDistance(curve1, curve2)).toBeCloseTo(39.0328); + }); + + it('gives correct results 2', () => { + const curve1 = [ + { x: 63.44852183813086, y: 24.420192387119634 }, + { x: 19.472881275654252, y: 77.306125067647 }, + { x: 22.0150089075698, y: 5.115699052924483 }, + { x: 90.85925658487311, y: 80.37914225209231 }, + { x: 96.81784894898642, y: 81.33960258698878 }, + { x: 75.45756084113779, y: 96.87017085629488 }, + { x: 87.77706429291412, y: 15.70163068744641 }, + { x: 37.36893642596093, y: 44.86136460914203 }, + { x: 37.35720453846581, y: 90.65479959420186 }, + { x: 41.28185352889147, y: 34.02195976325355 }, + { x: 27.65820587389076, y: 12.382281496757997 }, + { x: 42.43674529129338, y: 33.38959395979349 }, + { x: 3.377463737709774, y: 52.387593489371966 }, + { x: 50.93481600582428, y: 16.868378936261696 }, + { x: 68.46675900966153, y: 52.04265123799294 }, + { x: 1.9235036598383326, y: 55.87935516876048 }, + { x: 28.02334783421687, y: 98.08317663407114 }, + { x: 53.74539146366855, y: 33.27918237496243 }, + { x: 49.39670128874036, y: 47.59663728140997 }, + { x: 47.51990428391566, y: 11.23339071630216 }, + { x: 53.31256301680558, y: 55.4279696833061 }, + { x: 38.797168750480026, y: 26.172634107810833 }, + { x: 45.604650160570515, y: 71.69212699940685 }, + { x: 36.83931368726911, y: 38.74324014933978 }, + { x: 68.76987877419623, y: 1.2518741233677577 }, + { x: 91.27606575268427, y: 96.2141050404784 }, + { x: 24.407614843135406, y: 76.20115332073458 }, + { x: 8.764170623754097, y: 37.003392529458104 }, + { x: 52.97112238152346, y: 9.76631343977752 }, + { x: 88.85357966283867, y: 60.767524033054144 } + ]; + const curve2 = [{ x: 0, y: 0 }, { x: 14, y: 2.4 }, { x: 4, y: 4 }]; + + expect(frechetDistance(curve1, curve2)).toBeCloseTo(121.5429); + }); + + it("doesn't overflow the node stack if the curves are very long", () => { + const curve1 = rebalanceCurve([{ x: 1, y: 0 }, { x: 4, y: 4 }], { + numPoints: 5000 + }); + const curve2 = rebalanceCurve([{ x: 0, y: 0 }, { x: 4, y: 4 }], { + numPoints: 5000 + }); + + expect(frechetDistance(curve1, curve2)).toBe(1); + }); });