extends SkeletonModifier3D
class_name SpringBones

## Secondary motion for cloth, hair and accessories — a position-based cloth
## solver over the rig's own cosmetic bones.
##
## The animation clips drive the BODY and nothing else — tools/retarget.py
## deliberately exports no tracks at all for skirt, hair or accessory bones (see
## its `export_optimize_animation_keep_anim_armature=False`). This is what moves
## them, and it is the half of the pipeline that makes clothes read as clothes.
##
## WHY A PARTICLE SOLVER RATHER THAN A SPRING PER BONE
##
## The previous version integrated one damped spring per bone and then tried to
## repair the result with four more passes stacked behind it: resolve the
## collision against the spring's target, resolve it again against the answer,
## relax the cross-panel links and rebuild every pose from the corrected tips,
## then walk a separate cyclic-descent "lift" up each chain. Each stage wrote
## bone poses the next stage read back and partly undid, and the lift wrote poses
## that were never fed back into the spring state at all — so the spring spent
## the next frame pulling against a pose it did not know about. That is a
## feedback loop, and it is what the skirt "glitching out when it is moved"
## actually was. It also could not stop the thigh coming through, because every
## one of those stages could only ROTATE a bone about its own head, and the
## contacts that matter sit near that head.
##
## This solves the whole garment the way cloth is actually solved — the same
## shape as Magica Cloth 2's BoneCloth, which is what the Hoyoverse-style
## character pipelines use:
##
##   1. Every JOINT of every chain is a particle. Bone i spans particle i to
##      particle i+1, so a bone's head is a particle too and can move.
##   2. Predict: integrate each particle with inertia, gravity and wind, in the
##      chain anchor's frame so simply travelling costs nothing.
##   3. Relax ALL the constraints together, several Gauss-Seidel iterations:
##      bone length, a cone limit against the animated pose, a backstop on how
##      far the garment may stray, the horizontal links that hold neighbouring
##      panels together, and the leg/body colliders.
##   4. Convert the settled particles into bone rotations ONCE, at the end.
##
## Because collision is a constraint inside that loop rather than a repair pass
## behind it, nothing undoes it and it converges with everything else. And
## because a bone's HEAD is a particle, a contact with no rotational leverage is
## resolved by the whole panel moving — which is what a real skirt does when a
## thigh comes up under it, and what three separate bolt-on stages (a bodily
## chain push, an ancestor lift, a drape weight) were each approximating.
##
## Bone lengths, hull samples and the collider radii are MEASURED from the
## model's own geometry at build time and read from <model>.rig.json — a glTF
## skeleton carries no bone tails at all, so there is nothing in the skeleton
## itself that says which way a panel hangs or how thick a thigh is.

## Per class: w = stiffness pulling back toward the animated pose (rad/s), zeta =
## damping ratio (1.0 is critical, lower overshoots), gravity in m/s^2, wind =
## how far the class streams in the airflow, as a fraction of its own length.
##
## Hair is stiffer and lighter than cloth so it settles quickly instead of
## wobbling; a skirt is slacker and heavier so it lags and swings.
##
## Gravity is small because the AUTHORED rest pose already has the garment
## hanging; this only biases the droop while the body accelerates. It is not
## free — a constant force offsets the resting particle by g/w^2, so a large
## value pulls the hem below where it was modelled, straight into the thigh it
## then has to be pushed back out of.
##
## `stray` is the backstop — see MAX_STRAY. It is a per-class number because the
## two classes want opposite things from it. Hair has nothing to get out of the
## way of, so a tight leash is pure benefit and is what stops a long strand
## flailing. A garment has to be able to lie along a thigh that has swung
## horizontal in a slide, which is most of a panel's own length of travel, so the
## same leash on cloth is the clipping.
## `hinge` / `bend` are the swing limits — see MAX_SWING_ROOT. They are per class
## because the two ends of a chain mean different things to different garments. A
## SKIRT hangs from a waistband and hinges there freely: a knee coming up to hip
## height puts a thigh horizontally through where the front panel hangs, and the
## panel has to ride up onto it, which is most of a right angle. HAIR grows out
## of a scalp — a strand that hinges 70 degrees at the root has come off the
## head. Same constraint, opposite settings, and one shared number was wrong for
## both.
const TUNING := {
	"hair": {"w": 15.0, "zeta": 0.45, "gravity": 1.6, "wind": 0.30, "stray": 0.45,
		"hinge": 26.0, "bend": 52.0},
	"skirt": {"w": 12.0, "zeta": 0.48, "gravity": 2.0, "wind": 0.18, "stray": 1.25,
		"hinge": 78.0, "bend": 55.0},
	"cloth": {"w": 13.0, "zeta": 0.55, "gravity": 1.9, "wind": 0.18, "stray": 1.10,
		"hinge": 62.0, "bend": 52.0},
	"cape": {"w": 10.0, "zeta": 0.50, "gravity": 2.4, "wind": 0.35, "stray": 1.10,
		"hinge": 55.0, "bend": 55.0},
	"coat": {"w": 12.0, "zeta": 0.52, "gravity": 2.1, "wind": 0.22, "stray": 1.10,
		"hinge": 62.0, "bend": 52.0},
	"scarf": {"w": 13.0, "zeta": 0.55, "gravity": 1.8, "wind": 0.35, "stray": 0.70,
		"hinge": 40.0, "bend": 58.0},
	"ribbon": {"w": 15.0, "zeta": 0.52, "gravity": 1.6, "wind": 0.35, "stray": 0.70,
		"hinge": 35.0, "bend": 58.0},
	"tail": {"w": 14.0, "zeta": 0.55, "gravity": 1.5, "wind": 0.20, "stray": 0.60,
		"hinge": 30.0, "bend": 50.0},
	"sleeve": {"w": 15.0, "zeta": 0.60, "gravity": 1.6, "wind": 0.20, "stray": 0.60,
		"hinge": 30.0, "bend": 45.0},
	"breast": {"w": 21.0, "zeta": 0.65, "gravity": 1.0, "wind": 0.0, "stray": 0.30,
		"hinge": 14.0, "bend": 14.0},
	"bust": {"w": 21.0, "zeta": 0.65, "gravity": 1.0, "wind": 0.0, "stray": 0.30,
		"hinge": 14.0, "bend": 14.0},
}
const DEFAULT_TUNING := {"w": 14.0, "zeta": 0.55, "gravity": 1.8, "wind": 0.2,
	"stray": 0.55, "hinge": 40.0, "bend": 52.0}

## The most a bone may turn away from the direction the animated pose gives it —
## GRADED along the chain, tighter where it is anchored and looser at the free
## end. Angles still compound down a chain, so a strand curves well past these
## overall, and how far the garment may travel in total is bounded separately by
## MAX_STRAY.
##
## The ROOT limit is not the 20 degrees Magica Cloth 2's skirt guide suggests,
## and the difference is geometry rather than taste. Their figure is for a chain
## whose first bone is a large share of the skirt; the first segment of one of
## Taila's panels is 49 mm of a 288 mm panel, so 20 degrees there moves the panel
## below it by SEVENTEEN MILLIMETRES. Measured over a movement sweep, the solver
## saw 104 mm of thigh inside the front of the skirt, pushed on it every frame
## and every iteration, and left 95 mm of it — not because it was blind or
## mistuned, but because the joint it had to turn was allowed 17 mm of travel.
##
## The two limits also mean different things. Between segments this is BENDING
## stiffness — how sharply a sheet of cloth may crease — and it belongs tight.
## At the root it is a HINGE at the waistband, and a skirt swinging clear of a
## knee coming up under it really does hinge 40-odd degrees there.
##
## Collision is applied after this within each iteration and may exceed it — the
## same precedence Magica uses ("collision detection has a higher priority than
## this limit").
const MAX_SWING_ROOT := deg_to_rad(42.0)
const MAX_SWING_TIP := deg_to_rad(58.0)

## How far a particle may end up from where the animation alone would have put
## it, as a fraction of the WHOLE chain's length.
##
## The cone limit above is measured against the PARENT, so a four-bone chain can
## legitimately curl to 4x its per-joint limit and the far end of a long strand
## gets thrown a long way by errors that are individually reasonable. This bounds
## the total. It is a hard positional clamp and it is the single thing that turns
## "hair spazzing about" into hair — Magica calls it Max Distance and it is the
## first parameter their setup guide reaches for.
##
## Measured against the whole garment, NOT against the chain hanging above each
## joint. Scaling it per joint sounds right and throttles exactly the wrong
## place: the second segment of one of Taila's skirt panels has only 98 mm of
## chain above it, which allowed it 44 mm of travel while the thigh under it
## needed it to move 100 mm — so the leg went through the front of the skirt and
## every collider measurement said the solver had fired. The garment is one piece
## of cloth and the whole of it has the same room to move.
##
## Overridden per class in TUNING; this is only the fallback for a class with no
## entry of its own.
const MAX_STRAY := 0.55

## Solver iterations. Gauss-Seidel over length, cone, backstop, collision and the
## cross-panel links. Four is enough for the chains here (four bones deep); the
## extra two buy the collision convergence that stops a thigh showing through on
## the frames where a leg sweeps hard.
const ITERATIONS := 3
## Collision-and-length-only passes after the main loop. See the solve.
const SETTLE_PASSES := 1

## Detail level, set by the owner from how much the viewer can actually see —
## 0 is the full solve, 3 is the cheapest.
##
## This is not optional polish. The solver relaxes 127 joints against five
## capsules several times a frame, in GDScript, and collision is three quarters
## of that: a full-quality character costs about 2.6 ms, which is a sixth of a
## 60 fps frame for ONE of them. A shooter has several on screen and only the
## nearest is being looked at closely enough for a hem to matter.
##
## Each level drops a pass or the collision, in the order that costs the least
## visually — the shape constraints hold the silhouette on their own, so a
## distant character still has a skirt that swings; it just stops being carefully
## kept off a thigh nobody can resolve at that range.
const LOD_ITERS := [3, 2, 2, 1]
const LOD_SETTLE := [1, 1, 0, 0]
const LOD_COLLIDE := [true, true, true, false]
var lod: int = 0: set = _set_lod
var _iters: int = ITERATIONS
var _settles: int = SETTLE_PASSES
var _collide: bool = true

func _set_lod(v: int) -> void:
	lod = clampi(v, 0, LOD_ITERS.size() - 1)
	_iters = LOD_ITERS[lod]
	_settles = LOD_SETTLE[lod]
	_collide = LOD_COLLIDE[lod]
## Inverse mass of a chain's ROOT joint, relative to the joints below it.
##
## Not zero. A contact against the top of a thigh has almost no lever on the bone
## that owns it, so no rotation can clear it — which is what defeated every
## earlier round of tuning. Letting the root carry a little of the correction
## lets the whole panel ride up over the leg instead, rigidly, so the waistband
## seam cannot be stretched open by it. Small, and hard-clamped by ROOT_STRAY, so
## the garment can never detach from the body.
const ROOT_INV_MASS := 0.25
const ROOT_STRAY := 0.045
## A frame delta longer than this is a hitch or a load spike. Integrating it
## launches every chain across the map, so it is clamped instead.
const MAX_STEP := 1.0 / 30.0
## Origin jump (metres in one frame) that means a teleport — respawn, or the
## model being reparented — rather than movement. Chains snap instead of whip.
const TELEPORT := 1.5
## Speed at which the wind term reaches full strength.
const WIND_SPEED := 9.0
## WHY THERE IS NO "DRAPE" TERM HERE ANY MORE.
##
## There used to be one: each cloth bone measured how much of it rested on each
## leg, and took a share of that leg's motion before the solver ran, so a thigh
## arrived with the skirt already moving out of its way. Hoyoverse-style rigs do
## carry a partial constraint from the leg onto the upper skirt bones, and the
## idea is sound — but it was here to paper over a solver that could not push a
## panel off a thigh at all, and it cost far more than it bought once that was
## fixed. Measured over the movement sweep in debug/cloth_clip_check.gd, with it
## against without:
##
##     run    101 -> 92 mm      fall     82 -> 49 mm
##     air     96 -> 88 mm      dash    136 -> 95 mm
##     slide   97 -> 96 mm      idle after a dash   103 -> 20 mm
##
## It was worse in every state but a walk, and that last row is the one that
## mattered: WITH the drape the garment never came home after a hard move, it
## just sat displaced. It also cost an order of magnitude in stability — 0.48
## deg/frame of skirt movement at a dead-still idle against 0.05 without, because
## the target it aimed at sat inside the leg the collision was pushing out of, so
## the two ran against each other forever. That is what the skirt "glitching out"
## was, and it is gone with the term that caused it.
##
## The mesh-neighbour data it used is still read: the horizontal links below need
## it, and that is what it was always most useful for.
## HORIZONTAL links between neighbouring cloth chains.
##
## A skirt's bones are linked vertically by the parent-child hierarchy and not at
## all sideways, so every panel solves as if the ones sewn to it did not exist.
## Magica Cloth 2's BoneCloth skirt guide calls connecting them "the most
## important work when expressing a skirt", and notes that without it "the
## accuracy of collision detection will be significantly reduced" — which is
## exactly what was happening here: each panel individually satisfied its
## constraints and the garment as a whole came apart.
##
## Relaxed inside the same iteration loop as everything else, so a panel pushed
## off a thigh takes the ones sewn to it with it.
const LINK_STIFF := 0.45
## The same link resisting COMPRESSION. Much weaker: see the solve loop.
const LINK_SQUASH := 0.12
## Most links per bone, strongest first. A skirt vertex is shared with a handful
## of neighbours; linking every one of them just stiffens the garment into a
## bell.
const LINK_MAX := 4
## Classes a leg may carry. Cloth that hangs over the hips — not hair.
## Ceiling on how fast a particle may travel relative to its anchor, as a
## multiple of the bone length it hangs from. Bounds the tip speed to something
## proportional to the bone instead of a flat number that means nothing to a
## 20 mm hair segment and everything to a 140 mm skirt hem.
const MAX_RATE := 14.0
## Ceiling on how fast a chain's anchor is believed to be travelling, and how
## quickly that estimate may change. Backstops against a bad frame delta.
const MAX_CARRY := 20.0
const CARRY_SMOOTH := 0.15
const MAX_SUBSTEPS := 4
## Overlap a contact is allowed to keep before the solver acts on it.
##
## The rest-clearance table is measured in the AUTHORED rest pose, and a
## character never stands in it — the idle clip alone moves the legs enough that
## cloth hanging against a thigh grazes its capsule every frame. Without slop
## those grazes are real contacts and the idle buzzes. Six millimetres is well
## under anything the eye can find on a hem and it makes the idle silent.
const CONTACT_SLOP := 0.006
## Fractions along a bone tested against the colliders when the sidecar carries
## no hull for it. A skirt panel is a sheet, so testing only its tip lets the
## middle of the panel clip.
const SAMPLES := [0.35, 0.7, 1.0]
## How much of a collision correction the bone's HEAD takes when the contact sits
## at the head end. 1.0 would let a panel shear off its own waistband.
const COLLIDE_HEAD_SHARE := 0.85

var _chains: Array = []
var _colliders: Array = []
## [particle a, particle b, rest distance] for every horizontal link. Indices are
## into the flat particle arrays below.
var _links: Array = []
## Joint particles: one per bone plus a final tip, per chain. Kept in the chain
## ANCHOR's frame, so simply travelling through the world cancels exactly and the
## solver only ever deals with the motion that should actually move cloth.
var _q: PackedVector3Array = PackedVector3Array()
var _v: PackedVector3Array = PackedVector3Array()
## Scratch, world space, rebuilt every frame. Kept as members so the per-frame
## solve does not reallocate.
var _pw: PackedVector3Array = PackedVector3Array()       # working positions
var _ref: PackedVector3Array = PackedVector3Array()      # animated reference
var _rdir: PackedVector3Array = PackedVector3Array()     # reference bone dirs
var _rbasis: Array = []                                  # animated bone bases
var _imass: PackedFloat32Array = PackedFloat32Array()
## The collider capsules in world space, rebuilt ONCE per frame.
##
## They were being rebuilt inside the innermost loop — for every bone, against
## every collider, on every relaxation pass — which is two skeleton queries and
## two transform multiplies about five thousand times a frame for a result that
## cannot change while the solve is running. That alone was most of the 10.9 ms
## per character this used to cost.
var _cap_a: PackedVector3Array = PackedVector3Array()
var _cap_b: PackedVector3Array = PackedVector3Array()
var _cap_rh: PackedFloat32Array = PackedFloat32Array()
var _cap_rt: PackedFloat32Array = PackedFloat32Array()
## Scratch for one bone's posed hull points, reused instead of reallocated.
var _pts: PackedVector3Array = PackedVector3Array()
## Where each chain's anchor was last frame, to measure how fast it is carrying
## the cloth through the world.
var _anchor_prev: PackedVector3Array = PackedVector3Array()
var _carry_prev: PackedVector3Array = PackedVector3Array()
var _time: float = 0.0
var _settled: bool = false
## Force a fixed timestep instead of the engine clock. 0 = normal. Set by
## debug/cloth_settle_check.gd so settling can be measured at a REAL frame rate:
## headless runs uncapped, where the clamped sub-millisecond delta makes every
## chain look motionless whatever the tuning.
var fixed_delta: float = 0.0
## Diagnostics for the debug scripts — how many collider pushes fired, and how
## deep the worst contact on each bone got.
var _hits: int = 0
var _hit_frames: int = 0
var _hit_depth: Dictionary = {}
var _applied: Dictionary = {}
var _left: Dictionary = {}
var debug_bone: int = -1
var debug_swing: float = 0.0
var debug_len: float = 0.0
var debug_origin: Vector3 = Vector3.ZERO
var _calls: int = 0
var _frames_seen: Dictionary = {}
var _dt_min: float = 999.0
var _dt_max: float = 0.0


## How many times the solver ran per rendered frame, and the delta spread.
func debug_rate() -> String:
	var n: int = maxi(_frames_seen.size(), 1)
	var out := "calls/frame %.2f  dt %.5f..%.5f" % [float(_calls) / n, _dt_min, _dt_max]
	_calls = 0
	_frames_seen.clear()
	_dt_min = 999.0
	_dt_max = 0.0
	return out


## bone index -> [total correction applied (mm), deepest overlap seen (mm)].
func debug_effort_report() -> Dictionary:
	var out := {}
	for b in _applied:
		out[b] = [_applied[b] * 1000.0, _left.get(b, 0.0) * 1000.0]
	_applied.clear()
	_left.clear()
	return out


## bone index -> worst penetration (metres) seen since the last call.
func debug_hit_report() -> Dictionary:
	var d := _hit_depth.duplicate()
	_hit_depth.clear()
	return d


## Average collider pushes per frame since the last call. A chain that is
## settling should report ~0; a steady non-zero count means the solver and the
## colliders are fighting, which no amount of damping will settle.
func debug_collisions_per_frame() -> float:
	var v := float(_hits) / maxf(_hit_frames, 1.0)
	_hits = 0
	_hit_frames = 0
	return v


## Build from the sidecar written by tools/retarget.py. Returns how many bones
## are being driven, so the caller can log or disable itself when there are none.
##
## `hull_override` is bone name -> the points, in that bone's own rest space,
## that the collision should test. It comes from the MESH ITSELF at load time
## (see SkinnedPlayerModel._cloth_hulls) and it replaces the ten farthest-point
## samples the sidecar carries. Farthest-point sampling spans a panel's corners
## and hem and leaves its MIDDLE unsampled, which is precisely where a thigh
## comes up through a skirt: measured over a movement sweep the solver reported
## every contact resolved while 158 vertices were 95 mm inside a leg, because
## not one of them was a point it was looking at.
func setup(skel: Skeleton3D, info: Dictionary, hull_override: Dictionary = {}) -> int:
	_chains.clear()
	_links.clear()
	_colliders.clear()
	if skel == null or info.is_empty():
		return 0

	# Colliders first: each chain's per-bone clearance is measured against them.
	for c in info.get("colliders", []):
		var ca := skel.find_bone(String(c.get("bone", "")))
		var cb := skel.find_bone(String(c.get("child", "")))
		if ca >= 0 and cb >= 0:
			var tail := float(c.get("radius_tail", c.get("radius", 0.1)))
			_colliders.append({
				"a": ca, "b": cb,
				# How far down the bone the capsule starts. The top of a thigh is
				# hip, not limb — see tools/retarget.py::_leg_colliders.
				"from": float(c.get("from", 0.0)),
				"lid": bool(c.get("lid", false)),
				"rh": float(c.get("radius_head", tail)),
				"rt": tail,
			})

	var total := 0
	for entry in info.get("chains", []):
		var parent_name: String = entry.get("root_parent", "")
		var parent := skel.find_bone(parent_name) if parent_name != "" else -1
		if parent < 0:
			continue
		var bones := PackedInt32Array()
		var tips := PackedVector3Array()
		var hulls: Array = []
		var names: Array = entry.get("bones", [])
		var raw_tips: Array = entry.get("tips", [])
		var raw_hulls: Array = entry.get("hulls", [])
		for i in names.size():
			var idx := skel.find_bone(String(names[i]))
			if idx < 0:
				continue
			var tip := Vector3.ZERO
			if i < raw_tips.size():
				var t: Array = raw_tips[i]
				if t.size() == 3:
					tip = Vector3(t[0], t[1], t[2])
			if tip.length() < 0.001:
				continue  # no measurable extent — nothing to swing
			# Chains must stay a strict parent->child run: the particle chain
			# below assumes bone i+1 hangs off bone i, and a gap would put a
			# distance constraint across a joint that is not there.
			if bones.size() > 0 and skel.get_bone_parent(idx) != bones[bones.size() - 1]:
				break
			bones.append(idx)
			tips.append(tip)
			var hull: PackedVector3Array = hull_override.get(
				skel.get_bone_name(idx), PackedVector3Array())
			if hull.is_empty() and i < raw_hulls.size():
				for h in raw_hulls[i]:
					if h.size() == 3:
						hull.append(Vector3(h[0], h[1], h[2]))
			hulls.append(hull)
		if bones.is_empty():
			continue
		var cls := String(entry.get("class", ""))
		var tune: Dictionary = TUNING.get(cls, DEFAULT_TUNING)
		# Rest LINK geometry: particle i sits at bone i's head, particle n at the
		# last bone's tip. So the link below bone i reaches to the next bone's
		# head, which is NOT the same vector as that bone's mesh tip — the last
		# segment of a panel runs on past its own child joint to the hem.
		var seg := PackedFloat32Array()
		var segdir := PackedVector3Array()
		var span := 0.0
		for i in bones.size():
			var v: Vector3 = tips[i] if i == bones.size() - 1 \
				else skel.get_bone_rest(bones[i + 1]).origin
			if v.length() < 0.0005:
				v = tips[i]
			seg.append(v.length())
			segdir.append(v.normalized())
			span += v.length()
		_chains.append({
			"parent": parent,
			"bones": bones,
			"tips": tips,
			"seg": seg,
			"segdir": segdir,
			"span": span,
			"hulls": hulls,
			"class": cls,
			# Everything collides now, including hair. It used to be garments
			# only, because a collision push happened AFTER the integrator and so
			# was deaf to every spring parameter — long back hair reaching past
			# the hips got shoved out of a thigh and hauled back every frame at
			# stride frequency, which is the blur. Inside the relaxation there is
			# no such fight, and hair that does not collide simply passes through
			# the body: measured 373 mm of leg inside the hair over a movement
			# sweep, against 183 mm for the skirt that did collide.
			"radii": _rest_clearances(skel, bones, tips, hulls),
			"names": names,
			"neighbours": entry.get("neighbours", []),
			"w": float(tune["w"]),
			"zeta": float(tune["zeta"]),
			"gravity": float(tune["gravity"]),
			"wind": float(tune.get("wind", 0.0)),
			"stray": float(tune.get("stray", MAX_STRAY)),
			"hinge": deg_to_rad(float(tune.get("hinge",
				rad_to_deg(MAX_SWING_ROOT)))),
			"bend": deg_to_rad(float(tune.get("bend", rad_to_deg(MAX_SWING_TIP)))),
		})
		total += bones.size()

	# Flat particle layout: n+1 particles per chain.
	var k := 0
	for ci in _chains.size():
		_chains[ci]["base"] = k
		k += (_chains[ci]["bones"] as PackedInt32Array).size() + 1
	_build_links(skel)

	print("SpringBones: %d chains, %d particles, %d collision hull points" % [
		_chains.size(), k, _hull_points()])

	_q.resize(k)
	_v.resize(k)
	_pw.resize(k)
	_ref.resize(k)
	_rdir.resize(k)
	_imass.resize(k)
	_rbasis.resize(k)
	_cap_a.resize(_colliders.size())
	_cap_b.resize(_colliders.size())
	_cap_rh.resize(_colliders.size())
	_cap_rt.resize(_colliders.size())
	for c in _colliders.size():
		_cap_rh[c] = float(_colliders[c]["rh"])
		_cap_rt[c] = float(_colliders[c]["rt"])
	_anchor_prev.resize(_chains.size())
	_carry_prev.resize(_chains.size())
	# Inverse mass: heavier (stiffer) near the anchor, lighter toward the free
	# end, so a correction prefers to move the hem rather than the waistband.
	for ci in _chains.size():
		var base: int = _chains[ci]["base"]
		var n: int = (_chains[ci]["bones"] as PackedInt32Array).size()
		for i in n + 1:
			var along := 0.0 if n < 1 else float(i) / float(n)
			_imass[base + i] = ROOT_INV_MASS if i == 0 else lerpf(0.55, 1.0, along)
	_settled = false
	return total


## Total collision sample points across every chain — the solver's real cost.
func _hull_points() -> int:
	var n := 0
	for c in _chains:
		for h: PackedVector3Array in c["hulls"]:
			n += h.size() if h.size() > 0 else SAMPLES.size()
	return n


## The two ends of a collider capsule, with its `from` offset applied — so every
## place that tests against a limb agrees on where that limb starts.
func _capsule(skel: Skeleton3D, to_world: Transform3D, col: Dictionary) -> Array:
	var a: Vector3 = to_world * skel.get_bone_global_pose(col["a"]).origin
	var b: Vector3 = to_world * skel.get_bone_global_pose(col["b"]).origin
	return [a.lerp(b, float(col["from"])), b]


## Same, in the skeleton's REST pose.
func _capsule_rest(skel: Skeleton3D, col: Dictionary) -> Array:
	var a := skel.get_bone_global_rest(col["a"]).origin
	var b := skel.get_bone_global_rest(col["b"]).origin
	return [a.lerp(b, float(col["from"])), b]


## Per (bone, collider) capsule radius, capped so the REST pose never collides.
##
## The authored rest pose is by definition not clipping — the artist modelled the
## skirt over these legs. But a capsule sized from the thigh's own vertices is
## fat enough (0.11 m at Taila's hip) to swallow the skirt bones that hang right
## against it, and those bones then get shoved out and pulled straight back in
## every single frame: a limit cycle that never decays. That is cloth and hair
## "never settling", and it also leaves the collider saturated and useless
## against the clipping it exists to prevent.
##
## Capping to just inside the rest clearance makes the rest pose a valid state,
## so idle is quiet and a push only ever means the leg has genuinely swung into
## the cloth.
##
## The cap is PER POINT, not per bone. Scaling a whole bone's radius by its worst
## point switches collision off entirely for every skirt panel whose top
## naturally hangs against the thigh — which is all the ones that matter.
##
## And it is measured against the capsule radius WHERE THE POINT RESTS, not
## against the widest radius anywhere on the limb. A leg tapers 2:1 from hip to
## knee, so comparing every point to the hip radius declares the whole skirt
## already inside the capsule and hands all of it the reduced allowance meant for
## the waistband alone.
func _rest_clearances(skel: Skeleton3D, bones: PackedInt32Array,
		tips: PackedVector3Array, hulls: Array) -> Array:
	var out: Array = []
	for i in bones.size():
		var rest := skel.get_bone_global_rest(bones[i])
		var origin := rest.origin
		var tip: Vector3 = rest * tips[i]
		var pts := _sample_points(rest, origin, tip, hulls[i])
		var n := _colliders.size() * pts.size()
		var caps := PackedFloat32Array()
		caps.resize(n)
		# How far this bone's geometry reaches from its own head, for the
		# broad-phase skip in _collide_bone.
		var span := 0.0
		for p: Vector3 in pts:
			span = maxf(span, p.distance_to(origin))
		for c in _colliders.size():
			var col: Dictionary = _colliders[c]
			var ends := _capsule_rest(skel, col)
			var a: Vector3 = ends[0]
			var b: Vector3 = ends[1]
			var ab := b - a
			var d2 := ab.length_squared()
			for j in pts.size():
				var t: float = 0.0 if d2 < 0.000001 \
					else clampf((pts[j] - a).dot(ab) / d2, 0.0, 1.0)
				var d: float = (pts[j] - (a + ab * t)).length()
				# A limb tapers 2:1 hip to knee, so measure against the radius
				# WHERE THE POINT IS.
				var local: float = lerpf(float(col["rh"]), float(col["rt"]), t)
				caps[c * pts.size() + j] = maxf(local, d * 0.9) if d >= local else d * 0.9
		out.append({"cap": caps, "span": span})
	return out


## Distance links between bones of DIFFERENT chains that share mesh vertices.
##
## The garment's own topology decides what is sewn to what — `neighbours` in the
## sidecar is the weight each pair of cloth bones shares over the same vertices,
## which is the artist's answer and not a guess from names or rest distance.
##
## Only cross-chain pairs: within a chain the particle links already hold the
## bones together, and re-stating that here only fights them.
func _build_links(skel: Skeleton3D) -> void:
	_links.clear()
	# bone -> [chain, index]; the particle at that bone's TIP is base+index+1.
	var where := {}
	for ci in _chains.size():
		var bones: PackedInt32Array = _chains[ci]["bones"]
		for i in bones.size():
			where[bones[i]] = [ci, i]

	var seen := {}
	for ci in _chains.size():
		var chain: Dictionary = _chains[ci]
		if String(chain["class"]) == "hair":
			continue   # hair strands hang free; linking them stiffens them into rope
		var bones: PackedInt32Array = chain["bones"]
		var names: Array = chain["names"]
		var nbrs: Array = chain["neighbours"]
		for i in bones.size():
			var at: int = names.find(skel.get_bone_name(bones[i]))
			if at < 0 or at >= nbrs.size():
				continue
			var nb: Dictionary = nbrs[at]
			# Strongest few only — see LINK_MAX.
			var ranked: Array = nb.keys()
			ranked.sort_custom(func(a, b): return float(nb[a]) > float(nb[b]))
			var made := 0
			for other_name in ranked:
				if made >= LINK_MAX:
					break
				var oi := skel.find_bone(String(other_name))
				if oi < 0 or not where.has(oi) or where[oi][0] == ci:
					continue
				var pa: int = int(chain["base"]) + i + 1
				var ob: Array = where[oi]
				var pb: int = int(_chains[ob[0]]["base"]) + int(ob[1]) + 1
				var key := "%d_%d" % [mini(pa, pb), maxi(pa, pb)]
				if seen.has(key):
					continue
				seen[key] = true
				made += 1
				# Rest separation of the two TIPS, which is what the links hold.
				var oc: Dictionary = _chains[ob[0]]
				var wa: Vector3 = skel.get_bone_global_rest(bones[i]) \
					* (chain["tips"] as PackedVector3Array)[i]
				var wb: Vector3 = skel.get_bone_global_rest(oi) \
					* (oc["tips"] as PackedVector3Array)[ob[1]]
				_links.append([pa, pb, wa.distance_to(wb)])
	print("SpringBones: %d horizontal links between cloth chains" % _links.size())


static func _sample_points(xform: Transform3D, origin: Vector3, tip: Vector3,
		hull: PackedVector3Array) -> Array:
	var out: Array = []
	if hull.size() > 0:
		for h in hull:
			out.append(xform * h)
		return out
	for s in SAMPLES:
		out.append(origin.lerp(tip, s))
	return out


static func _seg_point_distance(p: Vector3, a: Vector3, b: Vector3) -> float:
	var ab := b - a
	var d2 := ab.length_squared()
	var t := 0.0 if d2 < 0.000001 else clampf((p - a).dot(ab) / d2, 0.0, 1.0)
	return p.distance_to(a + ab * t)


## The ENGINE's frame delta, never a wall clock.
##
## `Time.get_ticks_usec()` deltas measured 0.002 s to 0.033 s inside a single
## second — a 16x spread. That alone is survivable, but the anchor velocity
## divides a position change by it, and that position change came from the
## animation advancing by the ENGINE's delta. When the two disagree the computed
## velocity is wrong by that ratio, so the wind and the damping reference are
## wrong by 16x frame to frame.
func _delta() -> float:
	if fixed_delta > 0.0:
		return fixed_delta
	var dt := get_physics_process_delta_time() if Engine.is_in_physics_frame() \
		else get_process_delta_time()
	return clampf(dt, 1.0 / 480.0, MAX_STEP)


# ── The solve ─────────────────────────────────────────────────────────────────

func _process_modification() -> void:
	var skel := get_skeleton()
	if skel == null or _chains.is_empty():
		return
	var dt := _delta()
	_calls += 1
	_frames_seen[Engine.get_process_frames()] = true
	_dt_min = minf(_dt_min, dt)
	_dt_max = maxf(_dt_max, dt)
	_time += dt
	_hit_frames += 1
	var to_world := skel.global_transform
	for c in _colliders.size():
		var ends := _capsule(skel, to_world, _colliders[c])
		_cap_a[c] = ends[0]
		_cap_b[c] = ends[1]

	# 1. Where the animation alone would put every joint this frame, including
	#    the share of each leg's swing the cloth resting on it is carried by.
	var anchors: Array = []
	for ci in _chains.size():
		anchors.append(_build_reference(skel, to_world, ci))

	# 2. Integrate. The state lives in the anchor's frame, so travelling at a
	#    steady speed excites nothing and only real acceleration, rotation, wind
	#    and gravity move the cloth.
	for ci in _chains.size():
		_predict(skel, ci, anchors[ci], dt)

	# 3. Relax everything together.
	#
	# LINKS FIRST, chains second, so the last thing to touch any particle in an
	# iteration is its collision. Order is not cosmetic here. With the chains
	# first, every iteration ended by pulling neighbouring panels back toward
	# their rest separation — straight back into the leg that had just been
	# cleared. Measured on a slide: 93 mm of thigh inside a front panel went in,
	# 95 mm came out, and with the links switched off entirely the same frame
	# solved down to 27 mm. Collision outranks the garment's own shape, which is
	# also the precedence Magica Cloth 2 documents.
	for _it in _iters:
		for L in _links:
			var a: int = L[0]
			var b: int = L[1]
			var d: Vector3 = _pw[b] - _pw[a]
			var cur := d.length()
			if cur < 0.00001:
				continue
			var wa: float = _imass[a]
			var wb: float = _imass[b]
			var sum := wa + wb
			if sum <= 0.0001:
				continue
			# Cloth is stiff in TENSION and nearly free in COMPRESSION — a
			# garment gathers and folds happily but does not tear. Pulling
			# equally hard on both is what made these links fight the legs: a
			# panel riding up over a thigh opens away from the ones sewn beside
			# it, and a symmetric link reads that as the garment coming apart.
			var rest: float = float(L[2])
			var k: float = LINK_STIFF if cur > rest else LINK_SQUASH
			var move: Vector3 = d * ((cur - rest) / cur * k / sum)
			_pw[a] += move * wa
			_pw[b] -= move * wb
		for ci in _chains.size():
			_solve_chain(skel, to_world, ci, true)
	# A short tail of COLLISION AND BONE LENGTH ONLY.
	#
	# The cone limit and the backstop are shape constraints — they say what the
	# garment should look like — and re-imposing them after each collision means
	# the two argue instead of converging: quadrupling the main iteration count
	# moved a 92 mm overlap to 76 mm and stalled, which is a standing fight
	# rather than slow relaxation. These last passes let the collision finish,
	# holding only the constraint that cannot be given up (a bone has a length).
	# Shape has already been decided above and a few millimetres of overlap is
	# all that is left to clear, so nothing visible is given away by it.
	for _t in _settles:
		for ci in _chains.size():
			_solve_chain(skel, to_world, ci, false)

	# What the relaxation actually converged to, measured but not corrected. Read
	# alongside the same measurement taken off the finished skeleton, this is the
	# only way to tell "the constraints could not agree" from "they agreed and
	# the answer was lost on the way into the pose".
	if _want_residual:
		for ci in _chains.size():
			_measure_residual(skel, to_world, ci)

	# 4. Settled particles -> bone poses, and back into the anchor's frame.
	for ci in _chains.size():
		_write_poses(skel, anchors[ci], ci)
	_settled = true


## Deepest remaining overlap per bone, from the settled particles. Diagnostic
## only — nothing is corrected here.
func _measure_residual(skel: Skeleton3D, to_world: Transform3D, ci: int) -> void:
	var chain: Dictionary = _chains[ci]
	var bones: PackedInt32Array = chain["bones"]
	var seg: PackedFloat32Array = chain["seg"]
	var hulls: Array = chain["hulls"]
	var radii: Array = chain["radii"]
	var base: int = chain["base"]
	for i in bones.size():
		if i >= radii.size():
			continue
		var caps: PackedFloat32Array = radii[i].get("cap", PackedFloat32Array())
		if caps.is_empty():
			continue
		var a := base + i
		var d: Vector3 = _pw[a + 1] - _pw[a]
		if d.length() < 0.00001:
			continue
		var u := d.normalized()
		var posed := Transform3D(Basis(Quaternion(_rdir[a], u)) * (_rbasis[a] as Basis), _pw[a])
		var pts := _sample_points(posed, _pw[a], _pw[a] + u * seg[i], hulls[i])
		var worst := 0.0
		for c in _colliders.size():
			var col: Dictionary = _colliders[c]
			var ends := _capsule(skel, to_world, col)
			var ca: Vector3 = ends[0]
			var ab: Vector3 = ends[1] - ca
			var d2 := ab.length_squared()
			for j in pts.size():
				var p: Vector3 = pts[j]
				var t: float = 0.0 if d2 < 0.000001 \
					else clampf((p - ca).dot(ab) / d2, 0.0, 1.0)
				var kk := c * pts.size() + j
				var allow: float = caps[kk] if kk < caps.size() else 1.0
				var r: float = minf(lerpf(float(col["rh"]), float(col["rt"]), t), allow)
				worst = maxf(worst, (r - CONTACT_SLOP) - p.distance_to(ca + ab * t))
		if worst > 0.0:
			_resid[bones[i]] = maxf(_resid.get(bones[i], 0.0), worst)


## Turn on the residual measurement above. Off by default — it is a second full
## collision sweep and the game does not need it.
var _want_residual: bool = false
var _resid: Dictionary = {}

func debug_residual_report() -> Dictionary:
	_want_residual = true
	var d := _resid.duplicate()
	_resid.clear()
	return d


## The chain's animated pose for this frame: every joint where rigidly following
## the body would put it: every joint where rigidly following the animation
## would leave it.
##
## Returns the anchor transform, and fills _ref / _rdir / _rbasis for the chain.
func _build_reference(skel: Skeleton3D, to_world: Transform3D, ci: int) -> Transform3D:
	var chain: Dictionary = _chains[ci]
	var bones: PackedInt32Array = chain["bones"]
	var segdir: PackedVector3Array = chain["segdir"]
	var seg: PackedFloat32Array = chain["seg"]
	var base: int = chain["base"]
	var anchor: Transform3D = to_world * skel.get_bone_global_pose(chain["parent"])
	var pw := anchor
	for i in bones.size():
		# Built from the REST offset, not the current pose, or last frame's
		# result would compound into a permanent drift.
		var rigid: Transform3D = pw * skel.get_bone_rest(bones[i])
		_ref[base + i] = rigid.origin
		_rbasis[base + i] = rigid.basis
		_rdir[base + i] = (rigid.basis * segdir[i]).normalized()
		pw = Transform3D(rigid.basis, rigid.origin + _rdir[base + i] * seg[i])
	_ref[base + bones.size()] = pw.origin
	return anchor


## Integrate one chain's particles and seed the working positions.
func _predict(skel: Skeleton3D, ci: int, anchor: Transform3D, dt: float) -> void:
	var chain: Dictionary = _chains[ci]
	var bones: PackedInt32Array = chain["bones"]
	var base: int = chain["base"]
	var n: int = bones.size()
	var w: float = chain["w"]
	var pull: float = w * w
	var damp: float = 2.0 * float(chain["zeta"]) * w
	var gravity := Vector3.DOWN * float(chain["gravity"])
	var seg: PackedFloat32Array = chain["seg"]

	# How fast the chain's anchor is travelling through the world — the reference
	# the wind is measured against. Low-passed: a raw per-call difference
	# alternates violently when the solver is stepped more than once for the same
	# animation frame, and the wind force then flips sign every other step.
	var carry := (anchor.origin - _anchor_prev[ci]) / dt if _settled else Vector3.ZERO
	carry = _carry_prev[ci].lerp(carry.limit_length(MAX_CARRY), CARRY_SMOOTH)
	_carry_prev[ci] = carry
	_anchor_prev[ci] = anchor.origin

	var inv_anchor := anchor.affine_inverse()
	var reset := not _settled
	if not reset:
		# A teleport — respawn, or the model being reparented. Snap, do not whip.
		var got: Vector3 = anchor * _q[base + n]
		if got.distance_to(_ref[base + n]) > TELEPORT:
			reset = true

	# Cloth streams backwards through the air it is moving into. A STEADY wind
	# gives a steady deflection, which is hair that blows back once and then
	# hangs there rigid — travelling at a constant speed excites nothing else now
	# that the solver works in the anchor's frame — so the gust varies. Two
	# incommensurate rates so the pattern does not read as a loop, and a
	# per-chain phase so strands do not move in lockstep.
	var wind_dir := Vector3.ZERO
	var wind_mag := 0.0
	if carry.length() > 0.05:
		var ph: float = float(ci) * 1.7
		var gust := minf(carry.length() / WIND_SPEED, 1.0) * float(chain["wind"])
		var back := -carry.normalized()
		var side := back.cross(Vector3.UP)
		side = side.normalized() if side.length() > 0.001 else Vector3.ZERO
		var t := 1.0 + 0.45 * sin(_time * 6.3 + ph) + 0.25 * sin(_time * 10.7 + ph * 1.6)
		wind_dir = back * t + side * (0.35 * sin(_time * 8.1 + ph))
		wind_mag = pull * gust

	# Substep so the pull can never overshoot, however stiff it is or however
	# long the frame was. A single explicit step is only stable while w*dt stays
	# small; past that it gains energy every frame and the chain flies off.
	var sub := clampi(int(ceil(dt * w / 0.4)), 1, MAX_SUBSTEPS)
	var h := dt / float(sub)
	if reset:
		_q[base] = inv_anchor * _ref[base]
		_v[base] = Vector3.ZERO
	# The root joint is not pinned outright — see ROOT_INV_MASS. It keeps a
	# fraction of last frame's bodily shift and relaxes back toward the body, so
	# a panel that had to ride up over a thigh comes home once the leg has gone
	# by instead of stepping back the instant the contact clears.
	_pw[base] = _ref[base] \
		+ ((anchor * _q[base]) - _ref[base]).limit_length(ROOT_STRAY) * 0.55
	for i in range(1, n + 1):
		var k := base + i
		if reset:
			_q[k] = inv_anchor * _ref[k]
			_v[k] = Vector3.ZERO
		var p: Vector3 = anchor * _q[k]
		var vel: Vector3 = anchor.basis * _v[k]
		var length: float = seg[i - 1]
		for _s in sub:
			var acc := (_ref[k] - p) * pull + gravity + wind_dir * (wind_mag * length)
			vel += acc * h
			vel *= exp(-damp * h)
			vel = vel.limit_length(length * MAX_RATE)
			p += vel * h
		_pw[k] = p
		_v[k] = anchor.basis.inverse() * vel


## One Gauss-Seidel sweep over a chain: bone length, cone limit, backstop,
## colliders. Runs top-down so each bone sees where its parent has just gone.
func _solve_chain(skel: Skeleton3D, to_world: Transform3D, ci: int,
		shape: bool) -> void:
	var chain: Dictionary = _chains[ci]
	var bones: PackedInt32Array = chain["bones"]
	var seg: PackedFloat32Array = chain["seg"]
	var hulls: Array = chain["hulls"]
	var radii: Array = chain["radii"]
	var base: int = chain["base"]
	var n: int = bones.size()
	var far: float = float(chain["span"]) * float(chain["stray"])
	var hinge: float = chain["hinge"]
	var bend: float = chain["bend"]
	# How far the bone ABOVE has swung out of the animated pose. Each bone's cone
	# limit is measured from its own reference carried through this, so the chain
	# curls instead of every joint being pinned to the animation independently.
	var carry_swing := Quaternion.IDENTITY
	for i in n:
		var a := base + i
		var b := base + i + 1
		var length: float = seg[i]

		# 1. Bone length. Symmetric, so a correction below propagates upward —
		#    which is the whole reason the joints are particles.
		var d: Vector3 = _pw[b] - _pw[a]
		var cur := d.length()
		if cur > 0.00001:
			var wa: float = _imass[a]
			var wb: float = _imass[b]
			var sum := wa + wb
			if sum > 0.0001:
				var fix: Vector3 = d * ((cur - length) / cur / sum)
				_pw[a] += fix * wa
				_pw[b] -= fix * wb
		else:
			_pw[b] = _pw[a] + _rdir[a] * length

		# 2. Cone limit, measured from where this bone would point if it FOLLOWED
		#    ITS PARENT — not from where the animation alone put it.
		#
		#    This is the difference between a chain that can curl and one that
		#    cannot. Measuring every joint against the animated pose pins the
		#    whole strand within one cone of it: the root may swing its 42
		#    degrees, and then its child, which has physically been carried along
		#    by that swing, is told it is already 42 degrees out of place and
		#    must come back. Nothing compounds, so a four-bone panel is no more
		#    mobile than a one-bone one, and a skirt that has to lie along a
		#    thigh in a slide simply cannot get there — measured, the solver
		#    reported the contact every frame and left 99 mm of leg inside the
		#    hem.
		#
		#    Carrying the parent's actual swing onto the child's reference makes
		#    the limit mean what it should: how far this bone may bend RELATIVE
		#    TO THE ONE ABOVE IT, which is bending stiffness, which is a property
		#    of cloth. Total travel stays bounded by the backstop below.
		var ref_dir: Vector3 = (carry_swing * _rdir[a]).normalized()
		if shape:
			var dir: Vector3 = _pw[b] - _pw[a]
			if dir.length() > 0.00001:
				dir = dir.normalized()
				var along := 0.0 if n < 2 else float(i) / float(n - 1)
				var limit: float = lerpf(hinge, bend, minf(along * 2.0, 1.0))
				var ang := ref_dir.angle_to(dir)
				if ang > limit:
					dir = ref_dir.slerp(dir, limit / ang).normalized()
				_pw[b] = _pw[a] + dir * length

			# 3. Backstop: however the cone limits compound, the joint stays
			#    within a fraction of the whole garment of where the animation
			#    put it. See MAX_STRAY.
			var stray: Vector3 = _pw[b] - _ref[b]
			if stray.length() > far:
				_pw[b] = _ref[b] + stray.normalized() * far

		# 4. Colliders, against the REAL GEOMETRY this bone drives.
		if _collide and i < radii.size():
			_collide_bone(skel, to_world, bones[i], a, b, length,
				hulls[i], radii[i])

		# The root's own leash back to the body, re-applied as soon as the bone
		# hanging off it has had its say — the bones below solve against wherever
		# the root ended up, so letting it drift and only clamping at the end
		# would have them chasing a position the garment is not allowed to hold.
		if i == 0:
			var off0: Vector3 = _pw[base] - _ref[base]
			if off0.length() > ROOT_STRAY:
				_pw[base] = _ref[base] + off0.normalized() * ROOT_STRAY

		# What this bone finally settled on, for its child's cone limit.
		var got: Vector3 = _pw[b] - _pw[a]
		carry_swing = Quaternion(_rdir[a], got.normalized()) \
			if got.length() > 0.00001 else Quaternion.IDENTITY


## Push one bone's mesh hull out of every collider by moving its two joints.
##
## The bone stick is not what clips — a skirt panel is a wide sheet, and measured
## over a movement sweep the bones sat ~1 mm clear of the legs while the thigh
## was 85 mm inside the skirt MESH. So the test is against the hull samples from
## the sidecar, which are the actual vertices this bone drives.
##
## The correction is shared between the bone's head and its tip by where the
## contact sits along the bone, which makes it a rigid motion of the panel rather
## than a hinge. A contact at the head end therefore lifts the WHOLE panel — the
## case no amount of rotation could ever reach, because rotating about a point
## moves that point not at all.
func _collide_bone(skel: Skeleton3D, to_world: Transform3D, bone: int,
		a: int, b: int, length: float, hull: PackedVector3Array,
		rec: Dictionary) -> void:
	var caps: PackedFloat32Array = rec.get("cap", PackedFloat32Array())
	if caps.is_empty() or _colliders.is_empty():
		return
	var dir: Vector3 = _pw[b] - _pw[a]
	if dir.length() < 0.00001:
		return
	var wa: float = _imass[a]
	var wb: float = _imass[b]
	if wa + wb <= 0.0001:
		return
	var span: float = maxf(float(rec.get("span", length)), length)
	# Per-BONE broad phase, before the hull is posed at all. Most cloth bones are
	# nowhere near a limb on most frames — Taila's hair is 48 of the 127 — and
	# posing a couple of dozen hull points for them, nine passes a frame, was
	# pure waste.
	var near := false
	for c in _colliders.size():
		if _seg_point_distance(_pw[a], _cap_a[c], _cap_b[c]) 				<= span + maxf(_cap_rh[c], _cap_rt[c]):
			near = true
			break
	if not near:
		return
	var n := _fill_points(a, b, length, hull)
	if n == 0:
		return

	for c in _colliders.size():
		var ca: Vector3 = _cap_a[c]
		var cb: Vector3 = _cap_b[c]
		var ab := cb - ca
		var d2 := ab.length_squared()
		var rh: float = _cap_rh[c]
		var rt: float = _cap_rt[c]
		# Broad phase. Nothing this bone drives can reach further than `span`
		# from its own head, so a limb further away than that plus its own radius
		# cannot be touching it. Skips most of the (bone, collider) pairs.
		if _seg_point_distance(_pw[a], ca, cb) > span + maxf(rh, rt):
			continue
		# ONE correction per limb, for the DEEPEST point on the bone.
		#
		# Not one per point. The hull is a couple of dozen samples spread over a
		# panel, and a thigh inside it violates most of them at once — applying
		# every violation in turn adds up to many times the one overlap that
		# actually exists and flings the panel off the character. The deepest
		# point is the one that has to clear; the rest are the same contact seen
		# from nearby, and the next iteration re-measures whatever is left.
		var deepest := 0.0
		var push := Vector3.ZERO
		var at := Vector3.ZERO
		for j in n:
			var p: Vector3 = _pts[j]
			var t: float = 0.0 if d2 < 0.000001 				else clampf((p - ca).dot(ab) / d2, 0.0, 1.0)
			var kk := c * n + j
			var allow: float = caps[kk] if kk < caps.size() else 1.0
			var r: float = minf(lerpf(rh, rt, t), allow)
			var away := p - (ca + ab * t)
			var dist := away.length()
			if dist < 0.0001:
				continue
			var pen := (r - CONTACT_SLOP) - dist
			if pen > deepest:
				deepest = pen
				push = (away / dist) * pen
				at = p
		if deepest <= 0.0:
			continue
		_hits += 1
		_hit_depth[bone] = maxf(_hit_depth.get(bone, 0.0), deepest)
		_left[bone] = maxf(_left.get(bone, 0.0), deepest)
		_applied[bone] = _applied.get(bone, 0.0) + deepest
		# Where along the bone that point hangs decides which end carries the
		# correction. A contact at the HEAD end therefore lifts the whole panel
		# — the case no rotation can ever reach, because rotating about a point
		# moves that point not at all.
		var u := (_pw[b] - _pw[a]).normalized()
		var s: float = clampf((at - _pw[a]).dot(u) / maxf(length, 0.0001), 0.0, 1.0)
		var share_a := (1.0 - s) * COLLIDE_HEAD_SHARE
		var share_b := s + (1.0 - s) * (1.0 - COLLIDE_HEAD_SHARE)
		# Normalise so THE CONTACT POINT moves by the full overlap — not so the
		# two endpoint moves add up to it.
		#
		# The point rides at fraction `s` along the bone, so it travels
		# (1-s)*head + s*tip. Splitting the push into two shares that sum to it
		# therefore delivers only about half of it to the place that is actually
		# inside the leg, and a contact in the middle of a bone converged at
		# roughly half the rate it should: measured, quadrupling the iteration
		# count moved a 92 mm overlap to 76 mm and no further.
		var eff: float = (1.0 - s) * share_a * wa + s * share_b * wb
		if eff <= 0.0001:
			continue
		_pw[a] += push * (share_a * wa / eff)
		_pw[b] += push * (share_b * wb / eff)
		# Re-read, so the next limb sees where the bone has just gone.
		if _fill_points(a, b, length, hull) == 0:
			return


## Pose one bone's collision hull into the scratch buffer. Returns how many
## points were written. Kept out of the collider loop and off the heap: this runs
## for every cloth bone on every relaxation pass, and building a fresh Array of
## Vector3 each time was a large share of the solver's whole cost.
func _fill_points(a: int, b: int, length: float,
		hull: PackedVector3Array) -> int:
	var dir: Vector3 = _pw[b] - _pw[a]
	if dir.length() < 0.00001:
		return 0
	var u := dir.normalized()
	# The bone's basis under the correction it has taken so far: the animated
	# basis, swung onto the direction the particles settled on.
	var basis: Basis = Basis(Quaternion(_rdir[a], u)) * (_rbasis[a] as Basis)
	var head: Vector3 = _pw[a]
	if hull.size() > 0:
		if _pts.size() < hull.size():
			_pts.resize(hull.size())
		for i in hull.size():
			_pts[i] = head + basis * hull[i]
		return hull.size()
	if _pts.size() < SAMPLES.size():
		_pts.resize(SAMPLES.size())
	for i in SAMPLES.size():
		_pts[i] = head + u * (length * float(SAMPLES[i]))
	return SAMPLES.size()


## Turn the settled particles into bone rotations, and store back the joint
## positions the skeleton ACTUALLY got.
##
## The feedback matters as much as the write. A joint's head is fixed by its
## parent's rotation, so the pose the skeleton ends up in is never quite the pose
## the constraints asked for. Integrating from what the solver WANTED rather than
## from what was rendered leaves a standing gap between the two, and that gap is
## a force nothing in the model agreed to — it is what made the old multi-pass
## version oscillate, because its final ancestor-lift pass wrote poses that were
## never fed back into the spring state at all.
func _write_poses(skel: Skeleton3D, anchor: Transform3D, ci: int) -> void:
	var chain: Dictionary = _chains[ci]
	var bones: PackedInt32Array = chain["bones"]
	var segdir: PackedVector3Array = chain["segdir"]
	var seg: PackedFloat32Array = chain["seg"]
	var base: int = chain["base"]
	var n: int = bones.size()
	var inv_anchor := anchor.affine_inverse()
	# Whatever bodily shift the colliders asked of the chain root. Every segment
	# below it inherits this through the hierarchy, so the panel travels as one
	# piece and the mesh it drives cannot be stretched by it — which is the one
	# way translating cloth is safe where a weight gradient is not.
	var shift: Vector3 = (_pw[base] - _ref[base]).limit_length(ROOT_STRAY)

	var pw := anchor
	for i in n:
		var bone: int = bones[i]
		var rigid: Transform3D = pw * skel.get_bone_rest(bone)
		if i == 0:
			rigid.origin += shift
		# This joint's head is now final; record where it really landed.
		_q[base + i] = inv_anchor * rigid.origin
		# The bone's own rest line, and where the particles want it to point.
		var rest_dir: Vector3 = (rigid.basis * segdir[i]).normalized()
		var want: Vector3 = _pw[base + i + 1] - rigid.origin
		var world := rigid
		if want.length() > 0.00001 and rest_dir.length_squared() > 0.5:
			world.basis = Basis(Quaternion(rest_dir, want.normalized())) * rigid.basis
		var local := pw.affine_inverse() * world
		skel.set_bone_pose_rotation(bone, local.basis.get_rotation_quaternion())
		if i == 0:
			skel.set_bone_pose_position(bone, local.origin)
		pw = world
	# The free tip, from the last bone's achieved orientation.
	var tip_dir: Vector3 = (pw.basis * segdir[n - 1]).normalized()
	_q[base + n] = inv_anchor * (pw.origin + tip_dir * seg[n - 1])