init_elements()
local xC, nc = ...
local xB = xC / tkz.phi
local xD = (xC * xC) / xB
local xJ = (xC + xD) / 2
local r = xD - xJ
z.A = point(0, 0)
z.B = point(xB, 0)
z.C = point(xC, 0)
L.AC = line(z.A, z.C)
-- z.B = L.AC:gold_ratio()
z.D = point(xD, 0)
z.c = z.C:rotation(-math.pi / 2, z.B)
z.a = z.A:rotation(math.pi / 2, z.C)
z.d = z.D:rotation(-math.pi / 2, z.C)
L.AC = line(z.A, z.C)
z.mAC = L.AC.mid
L.AB = line(z.A, z.B)
z.mAB = L.AB.mid
L.BC = line(z.B, z.C)
z.mBC = L.BC.mid
L.CD = line(z.C, z.D)
z.mCD = L.CD.mid
C.AC = circle(z.A, z.C)
for i = 1, nc do
	z["J" .. i] = point(xJ, 2 * r * i)
	z["H" .. i] = point(xJ, 2 * r * i - r)
	z["J" .. i .. "p"], z["H" .. i .. "p"] = C.AC:inversion(z["J" .. i], z["H" .. i])
	L.AJ = line(z.A, z["J" .. i])
	C.JH = circle(z["J" .. i], z["H" .. i])
	z["S" .. i], z["T" .. i] = intersection(L.AJ, C.JH)
	z["S" .. i .. "p"], z["T" .. i .. "p"] = C.AC:inversion(z["S" .. i], z["T" .. i])
	L.SpTp = line(z["S" .. i .. "p"], z["T" .. i .. "p"])
	z["I" .. i] = L.SpTp.mid
end