Changed line-offset algorithm to account for spherical coordinates
Since the current algorithm is only correct at the equator i'm suggesting the following change:
- Use spherical coordinate and dot product between the two points to get 'real' distance
- Transform a local cartesian coordiante system from the midpoint that accounts for the fact that longitudes are squashed closer to the poles.
Hi @PereUbu7. Understand it's been a while since you proposed this. Are you still interested in getting this merged with some assistance from us?
@smallsaucepan Hi there! I'd love to. Just got to refresh and recapitulate what this was all about :)
@smallsaucepan Hi again. I've written some unit tests for the process segment method in vanilla js. Could you guys help putting this in the test framework you're actually using and maybe remove any old tests if there are any? :)
` "use strict"
const IsWithin = (value1, value2, variance) => { return value1 > (value2 - variance) && value1 < (value2 + variance); }
const MeterInKm = 0.001; const KmInDegrees = 0.00899;
const TestAbsoluteDistance = (point1, point2) => { console.log("Testing absolut distance with points:"+ point1 + ";" + point2);
/* Offset points by one km */
const offsetPoints = turf.processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
/* Distances from original points to new points */
const distance1 = turf.distance(point1, offsetPoint1);
const distance2 = turf.distance(point2, offsetPoint2);
/* Distances should be one km within half a meter (0.05% error) */
console.assert(IsWithin(distance1, 1, 0.5*MeterInKm), distance1);
console.assert(IsWithin(distance2, 1, 0.5*MeterInKm), distance2);
};
const TestDistanceFromPoints = (point1, point2) => { console.log("Testing distance with points:"+ point1 + ";" + point2);
/* Offset points by one km */
const offsetPoints = turf.processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
/* Distances from original points to new points */
const distance1 = turf.distance(point1, offsetPoint1);
const distance2 = turf.distance(point2, offsetPoint2);
/* Distances should be same within half a meter (0.05% error) */
console.assert(IsWithin(distance1, distance2, 0.5*MeterInKm), distance1, distance2);
};
const TestBearing = (point1, point2) => { console.log("Testing bearing with points:" + point1 + ";" + point2);
/* Offset points by one km */
const offsetPoints = turf.processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
const bearing1 = turf.bearing(point1, point2);
const bearing2 = turf.bearing(offsetPoint1, offsetPoint2);
/* Bearing between new points should be same as between old points within 0.01 degree */
console.assert(IsWithin(bearing1, bearing2, 0.01));
};
TestAbsoluteDistance([0.0, 0.0], [0.0, 0.0 + KmInDegrees]); TestAbsoluteDistance([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]); TestAbsoluteDistance([0.0, 0.0], [0.0 + KmInDegrees, 0.0]); TestAbsoluteDistance([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestAbsoluteDistance([0.0, 45.0], [0.0, 45.0 + KmInDegrees]); TestAbsoluteDistance([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]); TestAbsoluteDistance([0.0, 45.0], [0.0 + KmInDegrees, 45.0]); TestAbsoluteDistance([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestAbsoluteDistance([0.0, 65.0], [0.0, 65.0 + KmInDegrees]); TestAbsoluteDistance([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]); TestAbsoluteDistance([0.0, 65.0], [0.0 + KmInDegrees, 65.0]); TestAbsoluteDistance([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestAbsoluteDistance([30.0, 65.0], [30.0, 65.0 + KmInDegrees]); TestAbsoluteDistance([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]); TestAbsoluteDistance([30.0, 65.0], [30.0 + KmInDegrees, 65.0]); TestAbsoluteDistance([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([0.0, 0.0], [0.0, 0.0 + KmInDegrees]); TestDistanceFromPoints([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]); TestDistanceFromPoints([0.0, 0.0], [0.0 + KmInDegrees, 0.0]); TestDistanceFromPoints([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestDistanceFromPoints([0.0, 45.0], [0.0, 45.0 + KmInDegrees]); TestDistanceFromPoints([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]); TestDistanceFromPoints([0.0, 45.0], [0.0 + KmInDegrees, 45.0]); TestDistanceFromPoints([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestDistanceFromPoints([0.0, 65.0], [0.0, 65.0 + KmInDegrees]); TestDistanceFromPoints([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]); TestDistanceFromPoints([0.0, 65.0], [0.0 + KmInDegrees, 65.0]); TestDistanceFromPoints([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([30.0, 65.0], [30.0, 65.0 + KmInDegrees]); TestDistanceFromPoints([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]); TestDistanceFromPoints([30.0, 65.0], [30.0 + KmInDegrees, 65.0]); TestDistanceFromPoints([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]);
TestBearing([0.0, 0.0], [0.0, 0.0 + KmInDegrees]); TestBearing([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]); TestBearing([0.0, 0.0], [0.0 + KmInDegrees, 0.0]); TestBearing([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestBearing([0.0, 45.0], [0.0, 45.0 + KmInDegrees]); TestBearing([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]); TestBearing([0.0, 45.0], [0.0 + KmInDegrees, 45.0]); TestBearing([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestBearing([0.0, 65.0], [0.0, 65.0 + KmInDegrees]); TestBearing([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]); TestBearing([0.0, 65.0], [0.0 + KmInDegrees, 65.0]); TestBearing([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestBearing([30.0, 65.0], [30.0, 65.0 + KmInDegrees]); TestBearing([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]); TestBearing([30.0, 65.0], [30.0 + KmInDegrees, 65.0]); TestBearing([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]); `
So, the TestAbsoluteDistance checks that the actual offset is correct, and the TestDistanceFromPoints checks that the two new points are at the same distance from their orignal points (on a local scale), and the TestBearing checks that the line that crosses the two new points has the same bearing as the old line (on a local scale again). Do you guys think we need some more tests?
Is is possible to add some comment to the lineOffset function that is is only correct locally (up to some 10s kms)?
Greetings, Martin
I think it's a good idea, but in that case we should mark the old method as deprecated and then try to implement the geodesic version correctly.
This PR contains the planar version. The old (current) method seems to be planar too, but seems only to be correct close to the equator.
@PereUbu7 this should get you started with the tests, it should copy/paste into packages/turf-line-offset/test.ts cleanly. You'll have to add @turf/distance and @turf/bearing as devDependencies to turf-line-offset. The only question I have is whether we want to move processSegment into its own file so that we can export it for testing, or if you can just rework the tests to use the public API of lineOffset instead of hitting processSegment directly. I'd prefer to just test the public API if we can get away with it.
const IsWithin = (value1: number, value2: number, variance: number) => {
return value1 > value2 - variance && value1 < value2 + variance;
};
const MeterInKm = 0.001;
const KmInDegrees = 0.00899;
const TestAbsoluteDistance = (
point1: [number, number],
point2: [number, number]
) => {
console.log("Testing absolut distance with points:" + point1 + ";" + point2);
/* Offset points by one km */
const offsetPoints = processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
/* Distances from original points to new points */
const distance1 = distance(point1, offsetPoint1);
const distance2 = distance(point2, offsetPoint2);
test(`TestAbsoluteDistance - ${point1}; ${point2}`, (t) => {
/* Distances should be one km within half a meter (0.05% error) */
t.ok(IsWithin(distance1, 1, 0.5 * MeterInKm));
t.ok(IsWithin(distance2, 1, 0.5 * MeterInKm));
t.end();
});
};
const TestDistanceFromPoints = (
point1: [number, number],
point2: [number, number]
) => {
/* Offset points by one km */
const offsetPoints = processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
/* Distances from original points to new points */
const distance1 = distance(point1, offsetPoint1);
const distance2 = distance(point2, offsetPoint2);
test(`TestDistanceFromPoints ${point1};${point2}`, (t) => {
/* Distances should be same within half a meter (0.05% error) */
t.ok(IsWithin(distance1, distance2, 0.5 * MeterInKm));
t.end();
});
};
const TestBearing = (point1: [number, number], point2: [number, number]) => {
console.log("Testing bearing with points:" + point1 + ";" + point2);
/* Offset points by one km */
const offsetPoints = processSegment(point1, point2, KmInDegrees);
const offsetPoint1 = offsetPoints[0];
const offsetPoint2 = offsetPoints[1];
const bearing1 = bearing(point1, point2);
const bearing2 = bearing(offsetPoint1, offsetPoint2);
test(`TestBearing ${point1};${point2}`, (t) => {
/* Bearing between new points should be same as between old points within 0.01 degree */
t.ok(IsWithin(bearing1, bearing2, 0.01));
t.end();
});
};
TestAbsoluteDistance([0.0, 0.0], [0.0, 0.0 + KmInDegrees]);
TestAbsoluteDistance([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]);
TestAbsoluteDistance([0.0, 0.0], [0.0 + KmInDegrees, 0.0]);
TestAbsoluteDistance([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestAbsoluteDistance([0.0, 45.0], [0.0, 45.0 + KmInDegrees]);
TestAbsoluteDistance([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]);
TestAbsoluteDistance([0.0, 45.0], [0.0 + KmInDegrees, 45.0]);
TestAbsoluteDistance([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestAbsoluteDistance([0.0, 65.0], [0.0, 65.0 + KmInDegrees]);
TestAbsoluteDistance([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]);
TestAbsoluteDistance([0.0, 65.0], [0.0 + KmInDegrees, 65.0]);
TestAbsoluteDistance([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestAbsoluteDistance([30.0, 65.0], [30.0, 65.0 + KmInDegrees]);
TestAbsoluteDistance([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]);
TestAbsoluteDistance([30.0, 65.0], [30.0 + KmInDegrees, 65.0]);
TestAbsoluteDistance([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([0.0, 0.0], [0.0, 0.0 + KmInDegrees]);
TestDistanceFromPoints([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]);
TestDistanceFromPoints([0.0, 0.0], [0.0 + KmInDegrees, 0.0]);
TestDistanceFromPoints([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestDistanceFromPoints([0.0, 45.0], [0.0, 45.0 + KmInDegrees]);
TestDistanceFromPoints([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]);
TestDistanceFromPoints([0.0, 45.0], [0.0 + KmInDegrees, 45.0]);
TestDistanceFromPoints([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestDistanceFromPoints([0.0, 65.0], [0.0, 65.0 + KmInDegrees]);
TestDistanceFromPoints([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]);
TestDistanceFromPoints([0.0, 65.0], [0.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([30.0, 65.0], [30.0, 65.0 + KmInDegrees]);
TestDistanceFromPoints([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]);
TestDistanceFromPoints([30.0, 65.0], [30.0 + KmInDegrees, 65.0]);
TestDistanceFromPoints([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]);
TestBearing([0.0, 0.0], [0.0, 0.0 + KmInDegrees]);
TestBearing([0.0 + KmInDegrees, 0.0], [0.0, 0.0 + KmInDegrees]);
TestBearing([0.0, 0.0], [0.0 + KmInDegrees, 0.0]);
TestBearing([0.0, 0.0 + KmInDegrees], [0.0 + KmInDegrees, 0.0]);
TestBearing([0.0, 45.0], [0.0, 45.0 + KmInDegrees]);
TestBearing([0.0 + KmInDegrees, 45.0], [0.0, 45.0 + KmInDegrees]);
TestBearing([0.0, 45.0], [0.0 + KmInDegrees, 45.0]);
TestBearing([0.0, 45.0 + KmInDegrees], [0.0 + KmInDegrees, 45.0]);
TestBearing([0.0, 65.0], [0.0, 65.0 + KmInDegrees]);
TestBearing([0.0 + KmInDegrees, 65.0], [0.0, 65.0 + KmInDegrees]);
TestBearing([0.0, 65.0], [0.0 + KmInDegrees, 65.0]);
TestBearing([0.0, 65.0 + KmInDegrees], [0.0 + KmInDegrees, 65.0]);
TestBearing([30.0, 65.0], [30.0, 65.0 + KmInDegrees]);
TestBearing([30.0 + KmInDegrees, 65.0], [30.0, 65.0 + KmInDegrees]);
TestBearing([30.0, 65.0], [30.0 + KmInDegrees, 65.0]);
TestBearing([30.0, 65.0 + KmInDegrees], [30.0 + KmInDegrees, 65.0]);