Resources Numerical accuracy Calculation policy
Writes points, surfaces, alignments, profiles, pipe networks and feature lines as LandXML, and reads them back.
The writer formats every number with a fixed number of decimal places: 6 for coordinates, elevations, lengths and radii, 8 for angles and slopes, 4 for the existing-ground profile and the surface areas. Coordinates are written as absolute values, not shifted. The reader parses numbers with Number and keeps them in 64-bit; PIs are rebuilt by crossing tangents.
Open item. Write coordinates, elevations and curve data with enough significant digits to read back the same 64-bit number instead of a fixed 6 (or 4) decimal places.
16 calculations in 1 source file. The code is read from the released app when this page loads. It is the calculation people run, not a copy kept for this page.
Tests bin/test-civilcraft.mjs, bin/test-civilcraft-spirals.mjs, bin/test-civilcraft-arcs.mjs
Formats every number the export writes.
text = v.toFixed(d) with trailing zeros removed; d = 6 unless stated
Rounded to d decimal places, 6 by default. A value below half of the last place is written as 0 and a value that is not finite as 0.
const num = (v, d = 6) => {
if (!Number.isFinite(v)) return '0';
const s = (Math.abs(v) < 0.5 * 10 ** -d ? 0 : v).toFixed(d);
return s.includes('.') ? s.replace(/0+$/, '').replace(/\.$/, '') : s;
};Writes a point as northing, easting and elevation.
Each of the three values goes through the 6-decimal format. The coordinates are absolute project coordinates, so a northing in the millions keeps 13 significant digits.
const ne = (p, d) => num(p.y, d) + ' ' + num(p.x, d);
const nez = (p, d) => ne(p, d) + ' ' + num(p.z, d);Writes a transition curve: length, radii, station, angle, total X and Y, long and short tangents, constant.
theta (degrees) = θ × 180 ÷ π; constant = √(r × L_s)
The angle is converted to degrees and written with 8 decimals; the other values with 6. The clothoid constant is one square root of r × L_s.
L.push(` <Spiral rot="${e.ccw ? 'ccw' : 'cw'}" spiType="${SPI_TYPES[e.tr] || 'clothoid'}" length="${num(e.len)}" radiusStart="${inn ? 'INF' : num(e.r)}" radiusEnd="${inn ? num(e.r) : 'INF'}" staStart="${num(e.sta)}" theta="${num(e.theta * 180 / Math.PI, 8)}" totalX="${num(e.Xc)}" totalY="${num(e.Yc)}" tanLong="${num(e.tanLong)}" tanShort="${num(e.tanShort)}"${!e.tr || e.tr === 'clothoid' ? ` constant="${num(Math.sqrt(e.r * e.Ls))}"` : ''}>` +
`<Start>${ne(e.p0)}</Start><PI>${ne(e.pi)}</PI><End>${ne(e.p1)}</End></Spiral>`);Writes a circular curve: radius, length, station, chord, delta, tangent, external and middle ordinate.
delta (degrees) = Δ × 180 ÷ π
The values come from the alignment element as they are; delta is written with 8 decimals and the rest with 6.
L.push(` <Curve rot="${e.ccw ? 'ccw' : 'cw'}" crvType="arc" radius="${num(e.r)}" length="${num(e.len)}" staStart="${num(e.sta)}" chord="${num(e.chord)}" delta="${num(e.delta * 180 / Math.PI, 8)}" tangent="${num(e.T)}" external="${num(e.E)}" midOrd="${num(e.M)}">` +
`<Start>${ne(e.p0)}</Start><Center>${ne(e.c)}</Center><End>${ne(e.p1)}</End><PI>${ne(e.pi)}</PI></Curve>`);Writes the ground line along an alignment as station and elevation pairs.
Stations and elevations of the existing-ground profile are written with 4 decimals, two fewer than the design profile's PVIs, which use 6.
for (const r of runs) if (r.length >= 2) L.push(` <ProfSurf name="${esc(eg.name)}"><PntList2D>${r.map((p) => num(p.sta, 4) + ' ' + num(p.z, 4)).join(' ')}</PntList2D></ProfSurf>`);Writes a surface's plan area, surface area and highest and lowest elevations, then its points and faces.
The two areas are written with 4 decimals and the elevations with 6. The areas come from the surface's own stats; each point is written through the 6-decimal format.
L.push(` <Definition surfType="TIN" area2DSurf="${num(st.area2d, 4)}" area3DSurf="${num(st.area3d, 4)}" elevMax="${num(st.maxZ)}" elevMin="${num(st.minZ)}">`);
L.push(' <Pnts>');
pts.forEach((p, i) => L.push(` <P id="${i + 1}">${nez(p)}</P>`));Writes each pipe's plan length and slope.
length = hypot(Δx, Δy); slope = (invert up − invert down) ÷ length
The length is the plan distance between structure centers, written with 6 decimals; the slope is written with 8 decimals and is 0 when the length is 0.
const len = Math.hypot(b.x - a.x, b.y - a.y);
L.push(` <Pipe name="${esc(p.name)}" refStart="${esc(sName.get(p.from))}" refEnd="${esc(sName.get(p.to))}" length="${num(len)}" slope="${num(len ? (p.invUp - p.invDn) / len : 0, 8)}"><CircPipe diameter="${num(p.dia)}"/></Pipe>`);The center, radius, sweep and length of an arc segment of a feature line.
k = (1 − b²) ÷ 4b; r = c(1 + b²) ÷ 4|b|; sweep = 4·atan b; L = r × |sweep|
Closed form from the chord and the bulge. A chord shorter than 1e-12 or a bulge of 0 or of 100 and more gives no arc. For a very small bulge the center is far away and k is large, so the center loses digits while r × sweep stays close to the chord.
/** The arc from a to b with bulge a.b (+ = left), or null (featureline.js arcOf, kept here so this module stands alone). */
function arcGeom(a, b) {
const bb = Number.isFinite(a.b) && a.b !== 0 && Math.abs(a.b) < 100 ? a.b : 0, dx = b.x - a.x, dy = b.y - a.y, c = Math.hypot(dx, dy);
if (!bb || !(c > 1e-12)) return null;
const k = (1 - bb * bb) / (4 * bb), r = (c * (1 + bb * bb)) / (4 * Math.abs(bb)), sweep = 4 * Math.atan(bb);
return { cx: (a.x + b.x) / 2 - k * dy, cy: (a.y + b.y) / 2 + k * dx, r, sweep, L: r * Math.abs(sweep) };
}Splits a text list on spaces and commas and parses each number.
Each token is parsed with Number, which gives the nearest 64-bit number to the decimal text. Nothing is rounded on reading.
const nums = (s) => { const t = String(s || '').trim(); return t ? t.split(/[\s,]+/).map(Number) : []; };Decides whether the file is metric or in feet from its linear unit.
The linear unit is looked up in METRES_PER (units.js); a unit longer than 0.5 m counts as metric, anything else as feet. The scaling itself is done by scaleImport with a factor the caller passes.
out.units = METRES_PER[out.linearUnit] !== undefined && METRES_PER[out.linearUnit] > 0.5 ? 'm' : 'ft';Turns lines, curves and spirals into PIs with radius and spiral lengths.
PI = p₁ + d₁·a, a = ((p₂ − p₁) × d₂) ÷ (d₁ × d₂)
For a curve with spirals the PI is the crossing of the first and last tangents; the tangent directions are differences of file points or the radius vector turned 90 degrees. The parallel test scales 1e-12 by the lengths of both direction vectors. A plain curve uses the file's own PI when it has one. For a small deflection the division by a small cross product magnifies rounding in the file's points.
/**
* 0.7 D7.35: the PIs of a CoordGeom's elements (Line, Curve, Spiral). A curve with a spiral either side (or one side) becomes one
* PI with a radius and ls1 / ls2: tangent - spiral - curve - spiral - tangent, the PI being where the first and last tangents meet.
* null for what is not that shape (a spiral between two curves, a spiral that does not start or end in a straight).
*/
function readCoordGeom(els, pt) {
const inf = (v) => String(v).trim().toUpperCase() === 'INF' || !Number.isFinite(Number(v)) || Number(v) === 0;
const isIn = (e) => e.name === 'Spiral' && inf(e.attrs.radiusStart) && !inf(e.attrs.radiusEnd);
const isOut = (e) => e.name === 'Spiral' && !inf(e.attrs.radiusStart) && inf(e.attrs.radiusEnd);
const tangent = (p, c, rot) => { const vx = p.x - c.x, vy = p.y - c.y; return rot === 'cw' ? { x: vy, y: -vx } : { x: -vy, y: vx }; };
const meet = (p1, d1, p2, d2) => {
const den = d1.x * d2.y - d1.y * d2.x;
if (Math.abs(den) < 1e-12 * Math.hypot(d1.x, d1.y) * Math.hypot(d2.x, d2.y)) return null;
const a = ((p2.x - p1.x) * d2.y - (p2.y - p1.y) * d2.x) / den;
return { x: p1.x + d1.x * a, y: p1.y + d1.y * a };
};
const s0 = pt(els[0], 'Start');
if (!s0) return null;
const pis = [{ x: s0.x, y: s0.y, r: 0 }];
let i = 0;
while (i < els.length) {
const e = els[i];
if (e.name === 'Line') {
const en = pt(e, 'End');
if (!en) return null;
if (els[i + 1]?.name === 'Line') pis.push({ x: en.x, y: en.y, r: 0 });
i++; continue;
}
let j = i, sIn = null, cv = null, sOut = null;
if (e.name === 'Spiral') { if (!isIn(e)) return null; sIn = e; j++; }
if (els[j]?.name === 'Curve') { cv = els[j]; j++; }
if (els[j] && isOut(els[j])) { sOut = els[j]; j++; }
if (!cv && !(sIn && sOut)) return null; // a lone spiral, or a spiral that does not end in a curve or another spiral
if (!sIn && !sOut) {
// a plain curve: its own PI, or where the tangents at its ends meet
const s = pt(cv, 'Start'), en = pt(cv, 'End'), c = pt(cv, 'Center');
if (!s || !en) return null;
let pi = pt(cv, 'PI');
const r = Number(cv.attrs.radius) || (c ? Math.hypot(s.x - c.x, s.y - c.y) : 0);
if (!pi && c) pi = tangentsMeet(s, en, c);
if (!pi || !s || !en) return null;
pis.push({ x: pi.x, y: pi.y, r });
} else {
const c = cv ? pt(cv, 'Center') : null, rot = cv ? cv.attrs.rot : null;
let p1, d1, p2, d2;
if (sIn) { const s = pt(sIn, 'Start'), pi = pt(sIn, 'PI'); if (!s || !pi) return null; p1 = s; d1 = { x: pi.x - s.x, y: pi.y - s.y }; }
else { const s = pt(cv, 'Start'); if (!s || !c) return null; p1 = s; d1 = tangent(s, c, rot); }
if (sOut) { const en = pt(sOut, 'End'), pi = pt(sOut, 'PI'); if (!en || !pi) return null; p2 = en; d2 = { x: en.x - pi.x, y: en.y - pi.y }; }
else { const en = pt(cv, 'End'); if (!en || !c) return null; p2 = en; d2 = tangent(en, c, rot); }
const pi = meet(p1, d1, p2, d2);
const r = Number(cv?.attrs.radius) || Number((sIn || sOut).attrs.radiusEnd) || Number(sOut?.attrs.radiusStart);
if (!pi || !(r > 0)) return null;
const rec = { x: pi.x, y: pi.y, r };
// 0.7.1 71.26: the kind (the first spiral's when the two differ: one PI has one kind); a cubic parabola's length is its curve's
const tr = SPI_READ[String((sIn || sOut).attrs.spiType || 'clothoid').toLowerCase()] || 'clothoid';
const lsOf = (e) => { const v = Number(e.attrs.length) || 0; return tr === 'cubic' && v > 0 ? cubicParamFor(r, v) : v; };
if (sIn) rec.ls1 = lsOf(sIn);
if (sOut) rec.ls2 = lsOf(sOut);
if (tr !== 'clothoid') rec.tr = tr;
pis.push(rec);
}
i = j;
}
const last = pt(els[els.length - 1], 'End');
if (!last) return null;
pis.push({ x: last.x, y: last.y, r: 0 });
return pis;
}Where the tangents at a curve's start and end meet, from its center.
The tangents are the radius vectors turned 90 degrees, not normalized, so the cross product has units of length squared and the parallel threshold 1e-12 is absolute, unlike the scaled test used for curves with spirals.
/** The point where the tangents at a curve's start and end meet (from its centre). */
function tangentsMeet(s, e, c) {
const a = { x: -(s.y - c.y), y: s.x - c.x }, b = { x: -(e.y - c.y), y: e.x - c.x };
const den = a.x * b.y - a.y * b.x;
if (Math.abs(den) < 1e-12) return null;
const t = ((e.x - s.x) * b.y - (e.y - s.y) * b.x) / den;
return { x: s.x + a.x * t, y: s.y + a.y * t };
}The bulge of a feature-line arc from its start, center and end.
sweep = (a₁ − a₀) mod 2π; bulge = ± tan(sweep ÷ 4)
Two atan2 angles are subtracted and wrapped into 0…2π with the 64-bit value of 2π. A sweep below 1e-12 is treated as no arc. The radius in the file is not used, so the arc is fixed by the three points.
/** A LandXML Curve (Start, Center, End, rot cw|ccw) → the bulge of the segment from start to end (0 = not an arc). */
function curveBulge(s, c, e, rot) {
const a0 = Math.atan2(s.y - c.y, s.x - c.x), a1 = Math.atan2(e.y - c.y, e.x - c.x), ccw = rot !== 'cw';
let sw = ccw ? a1 - a0 : a0 - a1;
sw = ((sw % (2 * Math.PI)) + 2 * Math.PI) % (2 * Math.PI);
if (sw < 1e-12) return 0;
return (ccw ? 1 : -1) * Math.tan(sw / 4);
}The depth of a structure's sump below its lowest pipe invert.
sump depth = max(0, lowest invert − sump elevation)
One subtraction of two elevations read from the file, clamped at zero.
if (Number.isFinite(st.sumpElev) && v.length) st.sumpDepth = Math.max(0, Math.min(...v) - st.sumpElev);Joins consecutive lines and curves of a plan feature into one line.
A piece that starts more than 1e-6 (absolute, file units) from the end of the one before is counted as a gap and its start is kept as a new point; within 1e-6 the earlier point is kept and the new start is dropped.
const join = (a) => {
const last = pts[pts.length - 1];
if (!last) { pts.push(a); return pts[0]; }
if (Math.hypot(a.x - last.x, a.y - last.y) > 1e-6) { gaps++; pts.push(a); }
return pts[pts.length - 1];
};Multiplies every coordinate, elevation, station, radius and length of an import by one factor.
One 64-bit multiplication per value; a factor of exactly 1 returns the data untouched. Bulges and pipe diameters are not scaled.
export function scaleImport(d, k) {
if (k === 1) return d;
const p = (q) => (q ? { ...q, x: q.x * k, y: q.y * k, ...(q.z !== undefined && q.z !== null ? { z: q.z * k } : {}) } : q);
d.points = d.points.map(p);
for (const s of d.surfaces) { s.pts = s.pts.map(p); s.breaklines = s.breaklines.map((b) => ({ ...b, pts: b.pts.map(p) })); s.boundary = s.boundary && s.boundary.map(p); }
for (const a of d.alignments) {
a.start *= k; a.pis = a.pis.map((q) => ({ ...p(q), r: q.r * k, ...(q.ls1 ? { ls1: q.ls1 * k } : {}), ...(q.ls2 ? { ls2: q.ls2 * k } : {}), ...(q.tr ? { tr: q.tr } : {}) }));
for (const pr of a.profiles) pr.pvis = pr.pvis.map((v) => ({ sta: v.sta * k, z: v.z * k, L: v.L * k }));
for (const ps of a.surfProfiles) ps.pts = ps.pts.map((v) => ({ sta: v.sta * k, z: v.z * k }));
}
for (const f of d.features || []) f.pts = f.pts.map(p); // bulges are ratios: they stay
for (const n of d.networks || []) {
n.structures = n.structures.map((x) => ({ ...x, x: x.x * k, y: x.y * k, rim: x.rim === null ? null : x.rim * k, sumpDepth: x.sumpDepth * k }));
n.pipes = n.pipes.map((p) => ({ ...p, invUp: p.invUp === null ? null : p.invUp * k, invDn: p.invDn === null ? null : p.invDn * k }));
}
return d;
}The arithmetic is 64-bit, but the method or the data sets the accuracy: a solver tolerance, a mesh, a step size, a datum shift, a table. The limit is stated where the result is shown.