gregoryschwingmdphd/spinesurg-ct-annotator
0
1"""Draw both readers' marks on the same film, so the disagreement is visible.
2
3The pilot numbers say readers disagree by a median 0.0236 of image width and that most of
4it is head COUNT rather than precision. That is a claim about pictures, and it should be
5checked as one: six films, both readers overlaid, chosen to span the range rather than to
6flatter it -- the two they agreed on most, two where they disagreed about how many heads
7there are, and the two worst.
8
9Colour carries reader identity and nothing else: slot 1 blue, slot 2 orange, in fixed
10order, and every mark is also drawn with the reader's number beside it so identity is
11never colour alone. Strokes get a dark halo because a radiograph is black ground AND
12bright bone, and a single flat colour is illegible over one or the other.
13
14 python annot/render_pilot.py [--out DIR] [--n 6]
15"""
16from __future__ import annotations
17
18import argparse
19import json
20import math
21import os
22import statistics
23import sys
24from pathlib import Path
25
26from PIL import Image, ImageDraw, ImageFont
27
28PILOT = "gregoryschwingmdphd/xrsp-femhead-asp-pilot"
29IMAGES = "gregoryschwingmdphd/xrsp-femhead-images"
30TOL = 0.005
31
32# categorical slots 1 and 2, fixed order, never cycled: reader A is always blue
33COL = [(57, 135, 229), (235, 104, 52)]
34HALO = (10, 10, 12)
35DERIVED = (245, 200, 60) # the hip point each reader implies
36INK = (255, 255, 255)
37ROLES = ("A", "S", "P")
38
39
40def font(sz, bold=True):
41 for f in (("segoeuib.ttf", "arialbd.ttf", "DejaVuSans-Bold.ttf") if bold
42 else ("segoeui.ttf", "arial.ttf", "DejaVuSans.ttf")):
43 try:
44 return ImageFont.truetype(f, sz)
45 except OSError:
46 continue
47 return ImageFont.load_default()
48
49
50def halo_line(d, xy, fill, w):
51 d.line(xy, fill=HALO, width=w + 3)
52 d.line(xy, fill=fill, width=w)
53
54
55def halo_ellipse(d, box, fill, w):
56 d.ellipse(box, outline=HALO, width=w + 3)
57 d.ellipse(box, outline=fill, width=w)
58
59
60def label(d, xy, text, fill, f, anchor="mm"):
61 x, y = xy
62 for dx in (-2, -1, 0, 1, 2):
63 for dy in (-2, -1, 0, 1, 2):
64 if dx or dy:
65 d.text((x + dx, y + dy), text, font=f, fill=HALO, anchor=anchor)
66 d.text((x, y), text, font=f, fill=fill, anchor=anchor)
67
68
69def reads(c):
70 s = c.get("slots") or {}
71 return [s[k] for k in ("1", "2")
72 if (s.get(k) or {}).get("done") and (s.get(k) or {}).get("points")]
73
74
75def heads(p):
76 return [q for q in (p.get("heads") or []) if q]
77
78
79def pair_gap(A, B, asp):
80 d = lambda u, v: math.hypot(u[0] - v[0], (u[1] - v[1]) * asp) # noqa: E731
81 if len(A) == 2 and len(B) == 2:
82 return min(max(d(A[0], B[0]), d(A[1], B[1])),
83 max(d(A[0], B[1]), d(A[1], B[0])))
84 return min(d(u, v) for u in A for v in B)
85
86
87def mid(P):
88 return [sum(q[0] for q in P) / len(P), sum(q[1] for q in P) / len(P)]
89
90
91def draw_reader(d, p, W, H, ox, oy, k, sc, R_px):
92 """One reader's marks, in their slot colour, with their number on every head."""
93 col = COL[k % 2]
94 f_small = font(max(11, int(R_px * 0.30)))
95 lw = max(2, int(R_px * 0.045))
96 rr = max(3, int(R_px * 0.075))
97 hs = heads(p)
98 rad = p.get("radii") or []
99 lms = p.get("landmarks") or []
100 ex = p.get("extra") or []
101 for i, hc in enumerate(hs):
102 cx = (hc[0] * W - ox) * sc
103 cy = (hc[1] * H - oy) * sc
104 R = ((rad[i] if i < len(rad) and rad[i] else 0.085) * W) * sc
105 halo_ellipse(d, [cx - R, cy - R, cx + R, cy + R], col, lw)
106 t = R * 0.30
107 halo_line(d, [cx - t, cy, cx + t, cy], col, lw)
108 halo_line(d, [cx, cy - t, cx, cy + t], col, lw)
109 # the reader's number rides on the circle, so identity is never colour alone
110 label(d, (cx, cy - R - rr * 3), f"R{k+1}", col, f_small)
111 lm = lms[i] if i < len(lms) else {}
112 for role in ROLES:
113 q = (lm.get(role) or {})
114 if q.get("src") == "obs" and q.get("xy"):
115 x = (q["xy"][0] * W - ox) * sc
116 y = (q["xy"][1] * H - oy) * sc
117 d.rectangle([x - rr, y - rr, x + rr, y + rr], outline=HALO, width=lw + 2)
118 d.rectangle([x - rr, y - rr, x + rr, y + rr], outline=col, width=lw)
119 # role letters are deliberately NOT drawn here. With two readers on one
120 # small joint the six letters collide with each other and with the R1/R2
121 # tags, and the story this figure tells -- how the two circles differ --
122 # is carried by the circles. The squares still mark the landmarks.
123 for q in (ex[i] if i < len(ex) else []):
124 x = (q[0] * W - ox) * sc
125 y = (q[1] * H - oy) * sc
126 halo_ellipse(d, [x - rr * 0.7, y - rr * 0.7, x + rr * 0.7, y + rr * 0.7],
127 col, max(1, lw - 1))
128 if hs:
129 m = mid(hs)
130 mx, my = (m[0] * W - ox) * sc, (m[1] * H - oy) * sc
131 r2 = max(4, int(R_px * 0.10))
132 d.ellipse([mx - r2, my - r2, mx + r2, my + r2], fill=DERIVED, outline=HALO, width=2)
133 return (mx, my)
134 return None
135
136
137def main():
138 ap = argparse.ArgumentParser()
139 ap.add_argument("--out", default="annot/pilot_review")
140 ap.add_argument("--n", type=int, default=6)
141 a = ap.parse_args()
142
143 from huggingface_hub import snapshot_download, hf_hub_download
144 tok = os.environ.get("HF_TOKEN")
145 root = Path(snapshot_download(PILOT, repo_type="dataset",
146 allow_patterns="cases/*.json", max_workers=16,
147 token=tok))
148 cases = [json.loads(p.read_text()) for p in sorted(root.glob("cases/*.json"))]
149 rows = []
150 for c in cases:
151 R = reads(c)
152 if len(R) != 2:
153 continue
154 P = [r["points"] for r in R]
155 A, B = heads(P[0]), heads(P[1])
156 if not A or not B:
157 continue
158 W, H = float(P[0].get("w") or 1), float(P[0].get("h") or 1)
159 gap = pair_gap(A, B, H / W)
160 hipgap = math.hypot(mid(A)[0] - mid(B)[0], (mid(A)[1] - mid(B)[1]) * H / W)
161 rows.append(dict(case=c["case_id"], P=P, who=[r.get("annotator", "?") for r in R],
162 gap=gap, hip=hipgap, na=len(A), nb=len(B), W=W, H=H))
163 rows.sort(key=lambda r: r["gap"])
164 print(f" {len(rows)} films with two reads")
165
166 pick, seen = [], set()
167
168 def take(cand, tag):
169 for r in cand:
170 if r["case"] in seen:
171 continue
172 seen.add(r["case"])
173 r["tag"] = tag
174 pick.append(r)
175 return
176
177 take(rows, "closest agreement")
178 take(rows[1:], "close agreement")
179 take([r for r in rows if r["na"] != r["nb"]], "one saw two heads, one saw one")
180 take([r for r in rows if r["na"] != r["nb"]], "one saw two heads, one saw one")
181 take(rows[::-1], "worst disagreement")
182 take(rows[::-1], "second worst")
183 pick = pick[:a.n]
184
185 out = Path(a.out)
186 out.mkdir(parents=True, exist_ok=True)
187 panels = []
188 for r in pick:
189 cid = r["case"]
190 fp = None
191 for cand in (f"images/{cid}.jpg", f"{cid}.jpg"):
192 try:
193 fp = hf_hub_download(IMAGES, cand, repo_type="dataset", token=tok)
194 break
195 except Exception: # noqa: BLE001
196 continue
197 if not fp:
198 print(f" ! no film for {cid}")
199 continue
200 im = Image.open(fp).convert("RGB")
201 W, H = im.size
202 pts = [(q[0] * W, q[1] * H) for p in r["P"] for q in heads(p)]
203 rad = max([(p.get("radii") or [0.085])[0] * W for p in r["P"]] or [0.085 * W])
204 pad = rad * 2.6
205 x0 = max(0, int(min(x for x, _ in pts) - pad))
206 y0 = max(0, int(min(y for _, y in pts) - pad))
207 x1 = min(W, int(max(x for x, _ in pts) + pad))
208 y1 = min(H, int(max(y for _, y in pts) + pad))
209 crop = im.crop((x0, y0, x1, y1))
210 TARGET = 560
211 sc = TARGET / max(crop.width, crop.height)
212 crop = crop.resize((max(1, int(crop.width * sc)), max(1, int(crop.height * sc))),
213 Image.LANCZOS)
214 d = ImageDraw.Draw(crop)
215 hips = []
216 for k, p in enumerate(r["P"]):
217 hips.append(draw_reader(d, p, W, H, x0, y0, k, sc, rad * sc))
218 if all(hips):
219 halo_line(d, [hips[0][0], hips[0][1], hips[1][0], hips[1][1]], DERIVED, 2)
220
221 # caption strip under the film
222 fh = font(15)
223 fs = font(13, bold=False)
224 strip = 76
225 panel = Image.new("RGB", (crop.width, crop.height + strip), (16, 16, 18))
226 panel.paste(crop, (0, 0))
227 pd = ImageDraw.Draw(panel)
228 y = crop.height + 8
229 pd.text((10, y), f"{cid} {r['tag']}", font=fh, fill=INK)
230 pd.text((10, y + 21),
231 f"hip points {r['hip']:.4f} W apart = {r['hip']/TOL:.1f}x the tolerance",
232 font=fs, fill=(200, 200, 205))
233 for k in range(2):
234 cx = 14 + k * 190
235 pd.rectangle([cx - 4, y + 44, cx + 6, y + 54], fill=COL[k])
236 pd.text((cx + 14, y + 42),
237 f"R{k+1} {r['who'][k][:14]} · {len(heads(r['P'][k]))} head"
238 f"{'s' if len(heads(r['P'][k])) != 1 else ''}",
239 font=fs, fill=(210, 210, 215))
240 panels.append(panel)
241 panel.save(out / f"{cid}.png", optimize=True)
242 print(f" {cid:16s} {r['tag']:32s} hip gap {r['hip']:.4f} W "
243 f"({r['na']}/{r['nb']} heads)")
244
245 if panels:
246 cols = 3 if len(panels) >= 3 else len(panels)
247 rowsn = math.ceil(len(panels) / cols)
248 cw = max(p.width for p in panels)
249 ch = max(p.height for p in panels)
250 gap = 10
251 head_h = 58
252 grid = Image.new("RGB", (cols * cw + (cols + 1) * gap,
253 head_h + rowsn * ch + (rowsn + 1) * gap), (16, 16, 18))
254 gd = ImageDraw.Draw(grid)
255 gd.text((gap, 12), "Preston pilot — both readers on the same film",
256 font=font(20), fill=INK)
257 gd.text((gap, 36),
258 "blue = reader 1, orange = reader 2 · squares A/S/P · small rings = extra "
259 "rim points · yellow dot = the hip point that read implies",
260 font=font(13, bold=False), fill=(190, 190, 196))
261 for i, p in enumerate(panels):
262 cx = gap + (i % cols) * (cw + gap)
263 cy = head_h + gap + (i // cols) * (ch + gap)
264 grid.paste(p, (cx, cy))
265 grid.save(out / "_grid.png", optimize=True)
266 print(f"\n wrote {out}/_grid.png ({grid.width}x{grid.height})")
267 return 0
268
269
270if __name__ == "__main__":
271 sys.exit(main())
272 