concept-collection / hitandrun-interactive
Non-convex: add local-segment vs. union sampling mode
Add a 'local segment only' checkbox, shown only when non-convex is selected. Default (unchecked) keeps the standard walk that samples uniformly across every in-region segment the line crosses; checked restricts each step to the segment containing the current point — a local walk that can't hop across a concavity and does not sample uniformly. hit_and_run_general gains a local flag; the flag is threaded through the resample/newRegion protocol and also drives the movie overlay (regionSegments keeps only the segment straddling the current point).
Jeremy Magland <jmagland@flatironinstitute.org> committed commit a50c671f17d0 parent fcc5958 Browse files
5 changed files+185−83
CLAUDE.mdmodified+27−18View file
@@ -51,25 +51,34 @@ receive via `HTMLEventReceivedFcn`
5151 4. **(done)** Non-convex regions. A **non-convex region** checkbox switches to a
5252 star polygon (`make_region(false)`) sampled by `hit_and_run_general`; the
5353 movie draws the multiple in-region segments a line makes through it.
54-
55-The figure → script protocol (each request carries `convex`):
56-
57-- **Samples** slider (on release) / **Resample** →
58- `sendToMATLAB('resample', { n, x, y, convex })`. `x, y` are the *current
59- region's* vertices — the script is **stateless**, so it samples the region it
60- is given and replies with `sendEventToHTMLSource(src, 'samples', { x, y, n })`.
54+5. **(done)** Local vs. union sampling. A **local segment only** checkbox (shown
55+ only when non-convex) switches `hit_and_run_general` between sampling the
56+ union of all in-region segments (default; uniform) and only the segment
57+ through the current point (local; non-uniform). It drives both the sampler
58+ and the movie overlay (`regionSegments`'s `localOnly`).
59+
60+The figure → script protocol (each request carries `convex`, plus `local` for
61+non-convex regions):
62+
63+- **Samples** slider (on release) / **Resample** / **local** toggle →
64+ `sendToMATLAB('resample', { n, x, y, convex, local })`. `x, y` are the
65+ *current region's* vertices — the script is **stateless**, so it samples the
66+ region it is given and replies with `sendEventToHTMLSource(src, 'samples', {
67+ x, y, n })`.
6168 - **New region** / **non-convex** toggle → `sendToMATLAB('newRegion', { n,
62- convex })`. The script builds a new region of that type and replies with a full
63- `sendEventToHTMLSource(src, 'data', <payload>)` (same shape as the initial
64- `Data`, plus `convex`).
69+ convex, local })`. The script builds a new region of that type and replies with
70+ a full `sendEventToHTMLSource(src, 'data', <payload>)` (same shape as the
71+ initial `Data`, plus `convex`).
6572
66-`hitandrun_sampler.m`'s `sample_region(vx, vy, n, convex)` dispatches to
67-`hit_and_run` (convex, JIT) or `hit_and_run_general` (non-convex, interpreter).
68-The non-convex path can't JIT (sort + point-in-polygon), so `App.tsx` caps its
69-`N` lower via `SAMPLE_CHOICES_NONCONVEX`.
73+`hitandrun_sampler.m`'s `sample_region(vx, vy, n, convex, local)` dispatches to
74+`hit_and_run` (convex, JIT) or `hit_and_run_general(..., local)` (non-convex,
75+interpreter). The non-convex path can't JIT (sort + point-in-polygon), so
76+`App.tsx` caps its `N` lower via `SAMPLE_CHOICES_NONCONVEX`. `local` only applies
77+to non-convex regions (`App.tsx`: `useLocal = nonConvex && local`); toggling it
78+resamples the same region.
7079
71-Ideas for the next iteration: show the chain path / burn-in, animate the walk,
72-or a "uniform vs. non-uniform target" toggle.
80+Ideas for the next iteration: show the chain path / burn-in, or animate the
81+walk.
7382
7483 ## Client-side vs. server-side interactivity
7584
@@ -82,13 +91,13 @@ or a "uniform vs. non-uniform target" toggle.
8291
8392 | File | Purpose |
8493 |------|---------|
85-| `app/src/App.tsx` | Main React component; reads data, hosts the panel, non-convex toggle + movie `regionSegments` |
94+| `app/src/App.tsx` | Main React component; reads data, hosts the panel, non-convex + local toggles, movie `regionSegments` |
8695 | `app/src/render/RegionView.tsx` | 2D canvas renderer (region outline + sample dots + movie `segments`) |
8796 | `app/src/bridge.ts` | `onData` / `onHostEvent` / `sendToMATLAB` helpers (generic; identical across these projects) |
8897 | `hitandrun_sampler.m` | Opens the figure, sends data, handles `resample`/`newRegion` (figure plumbing: `on_event`, `sample_region`, `pack_data`, `pack_samples`) |
8998 | `helpers/make_region.m` | Random region — convex polygon on an ellipse (no convhull, so the auto-running figure view doesn't race the qhull backend load), or `make_region(false)` for a star polygon |
9099 | `helpers/hit_and_run.m` | Convex hit-and-run sampler (half-plane chord intersection; JIT) |
91-| `helpers/hit_and_run_general.m` | Non-convex hit-and-run sampler (all in-region segments; interpreter) |
100+| `helpers/hit_and_run_general.m` | Non-convex hit-and-run sampler (union of all in-region segments, or `local` = only the one through the current point; interpreter) |
92101 | `hitandrun_demo.m` | User-facing driver: `addpath('helpers')` + seed + call the sampler |
93102
94103 ## Local iteration
README.mdmodified+20−11View file
@@ -25,9 +25,14 @@ Controls:
2525 - **non-convex region** — toggle between a convex region and a non-convex
2626 star-shaped one. The non-convex region uses the slower general sampler, so
2727 `N` is capped lower.
28+- **local segment only** *(non-convex only)* — by default the step samples
29+ across *every* segment where the line crosses the region (the standard walk,
30+ uniform on the whole region); check this to instead sample only within the
31+ segment that contains the current point, a local walk that can't hop across a
32+ concavity and so does not sample uniformly.
2833 - **Play movie** — step through the algorithm: each step draws the chord — one
29- segment for a convex region, possibly several for a non-convex one — and the
30- point that landed on it.
34+ segment for a convex region, possibly several for a non-convex one (or just
35+ the local segment) — and the point that landed on it.
3136
3237 ## How it works
3338
@@ -39,9 +44,9 @@ Controls:
3944 The script and figure talk both ways: the script sends the region + samples via
4045 `uihtml(..., 'Data', ...)`, and the controls call back with
4146 `sendToMATLAB('resample' | 'newRegion', ...)` (each carrying whether the region
42-is `convex`), which re-runs the sampler and returns new points via
43-`sendEventToHTMLSource`. The script is stateless — the figure owns the region
44-and passes it back with each request.
47+is `convex` and, for non-convex, the `local` sampling mode), which re-runs the
48+sampler and returns new points via `sendEventToHTMLSource`. The script is
49+stateless — the figure owns the region and passes it back with each request.
4550
4651 ## Convex vs. non-convex
4752
@@ -55,12 +60,16 @@ and passes it back with each request.
5560 silently falling back). It relies on numbl's scalar-`rand()` JIT support;
5661 `rng(seed)` still controls the shared PRNG.
5762 - **Non-convex** — [`hit_and_run_general.m`](helpers/hit_and_run_general.m)
58- finds *every* crossing along the line, keeps the segments whose midpoint is
59- inside (a point-in-polygon test), and samples uniformly across their union, so
60- concavities are handled correctly. The sort + polygon tests don't JIT, so this
61- runs in the interpreter and `N` is capped lower. `make_region(false)` builds a
62- star polygon that is star-shaped about the origin, so the origin is a valid
63- interior start point.
63+ finds *every* crossing along the line and keeps the in-region segments (a
64+ point-in-polygon test on each interval's midpoint), so concavities are handled
65+ correctly. It then samples in one of two modes: **union** (default) picks a
66+ point uniformly across all those segments, which samples the whole region
67+ uniformly; **local** (the *local segment only* checkbox) restricts to the
68+ single segment straddling the current point — a local walk that can't cross a
69+ concavity, so it does *not* sample uniformly. The sort + polygon tests don't
70+ JIT, so this runs in the interpreter and `N` is capped lower.
71+ `make_region(false)` builds a star polygon that is star-shaped about the
72+ origin, so the origin is a valid interior start point.
6473
6574 ## Deploy
6675
app/src/App.tsxmodified+66−15View file
@@ -53,14 +53,17 @@ interface MovieState {
5353
5454 /** The in-region segment(s) hit-and-run samples along: the line through
5555 * (px,py) with direction (dx,dy), intersected with the polygon. One segment
56- * for a convex region, possibly several for a non-convex one. Mirrors the
57- * sampler's geometry, so it reproduces each step exactly. */
56+ * for a convex region, possibly several for a non-convex one. With `localOnly`
57+ * it keeps just the segment containing the current point (t = 0) — matching the
58+ * sampler's local-segment mode. Mirrors the sampler's geometry, so it
59+ * reproduces each step exactly. */
5860 function regionSegments(
5961 region: Points,
6062 px: number,
6163 py: number,
6264 dx: number,
63- dy: number
65+ dy: number,
66+ localOnly: boolean
6467 ): Segment[] {
6568 if (Math.hypot(dx, dy) < 1e-12) return [];
6669 const m = region.x.length;
@@ -80,14 +83,15 @@ function regionSegments(
8083 const segs: Segment[] = [];
8184 for (let k = 0; k < ts.length - 1; k++) {
8285 const tm = (ts[k] + ts[k + 1]) / 2;
83- if (pointInPolygon(region, px + tm * dx, py + tm * dy)) {
84- segs.push({
85- x0: px + ts[k] * dx,
86- y0: py + ts[k] * dy,
87- x1: px + ts[k + 1] * dx,
88- y1: py + ts[k + 1] * dy,
89- });
90- }
86+ if (!pointInPolygon(region, px + tm * dx, py + tm * dy)) continue;
87+ // Local mode: keep only the in-region interval straddling t = 0.
88+ if (localOnly && !(ts[k] <= 0 && ts[k + 1] >= 0)) continue;
89+ segs.push({
90+ x0: px + ts[k] * dx,
91+ y0: py + ts[k] * dy,
92+ x1: px + ts[k + 1] * dx,
93+ y1: py + ts[k + 1] * dy,
94+ });
9195 }
9296 return segs;
9397 }
@@ -117,6 +121,9 @@ export function App() {
117121 const [n, setN] = useState(DEFAULT_N);
118122 const [busy, setBusy] = useState(false);
119123 const [movie, setMovie] = useState<MovieState | null>(null);
124+ // Non-convex sampling mode: false = sample the union of all in-region
125+ // segments the line makes; true = only the segment through the current point.
126+ const [local, setLocal] = useState(false);
120127
121128 useEffect(() => {
122129 // Initial region + samples (script → page via `Data`).
@@ -170,9 +177,13 @@ export function App() {
170177
171178 const nonConvex = !!data && data.convex === false;
172179 const choices = nonConvex ? SAMPLE_CHOICES_NONCONVEX : SAMPLE_CHOICES;
180+ // `local` only applies to non-convex regions.
181+ const useLocal = nonConvex && local;
173182
174- // Re-draw `count` samples in the current region (script round-trip).
175- const resample = (count: number) => {
183+ // Re-draw `count` samples in the current region (script round-trip). The
184+ // sampling mode (`localMode`) can be overridden — the local checkbox passes
185+ // its new value directly, since the `local` state hasn't committed yet.
186+ const resample = (count: number, localMode: boolean = local) => {
176187 if (!data) return;
177188 setBusy(true);
178189 sendToMATLAB("resample", {
@@ -180,13 +191,14 @@ export function App() {
180191 x: data.region.x,
181192 y: data.region.y,
182193 convex: data.convex !== false,
194+ local: localMode,
183195 });
184196 };
185197
186198 // Generate a brand new region (convex or not) and sample it.
187199 const newRegion = (count: number, convex: boolean) => {
188200 setBusy(true);
189- sendToMATLAB("newRegion", { n: count, convex });
201+ sendToMATLAB("newRegion", { n: count, convex, local });
190202 };
191203
192204 // Checkbox: switch region type. Clamp N to the active set's max first.
@@ -197,6 +209,12 @@ export function App() {
197209 newRegion(clamped, !makeNonConvex);
198210 };
199211
212+ // Checkbox: switch the non-convex sampling mode; resample the same region.
213+ const setLocalMode = (value: boolean) => {
214+ setLocal(value);
215+ resample(n, value);
216+ };
217+
200218 const toggleMovie = () => {
201219 if (movie) {
202220 setMovie(null);
@@ -216,7 +234,14 @@ export function App() {
216234 const { x, y } = data.samples;
217235 cloud = prefixPoints(data.samples, i); // settled points 0..i-1
218236 from = { x: x[f], y: y[f] };
219- segments = regionSegments(data.region, x[f], y[f], x[i] - x[f], y[i] - y[f]);
237+ segments = regionSegments(
238+ data.region,
239+ x[f],
240+ y[f],
241+ x[i] - x[f],
242+ y[i] - y[f],
243+ useLocal
244+ );
220245 newPoint = { x: x[i], y: y[i] };
221246 }
222247
@@ -288,6 +313,21 @@ export function App() {
288313 non-convex region
289314 </label>
290315
316+ {nonConvex && (
317+ <label
318+ style={subCheckLabelStyle}
319+ title="Sample only the segment through the current point instead of every segment the line crosses"
320+ >
321+ <input
322+ type="checkbox"
323+ checked={local}
324+ disabled={!data || busy || !!movie}
325+ onChange={e => setLocalMode(e.target.checked)}
326+ />
327+ local segment only
328+ </label>
329+ )}
330+
291331 <button
292332 style={playBtnStyle}
293333 disabled={!canPlay || busy}
@@ -353,6 +393,17 @@ const checkLabelStyle: CSSProperties = {
353393 cursor: "pointer",
354394 };
355395
396+const subCheckLabelStyle: CSSProperties = {
397+ display: "flex",
398+ alignItems: "center",
399+ gap: 5,
400+ fontSize: 10,
401+ marginTop: 4,
402+ marginLeft: 14,
403+ color: "#475569",
404+ cursor: "pointer",
405+};
406+
356407 const sliderStyle: CSSProperties = {
357408 width: "100%",
358409 marginTop: 2,
helpers/hit_and_run_general.mmodified+54−27View file
@@ -1,12 +1,22 @@
1-function [sx, sy] = hit_and_run_general(vx, vy, N, nBurn)
1+function [sx, sy] = hit_and_run_general(vx, vy, N, nBurn, local)
22 %HIT_AND_RUN_GENERAL Hit-and-run for an arbitrary simple polygon (convex or
3-% non-convex). Along each random line it finds *every* in-region segment and
4-% samples uniformly across their union, so concavities are handled correctly
5-% (unlike the convex-only chord in hit_and_run.m). Starts at the origin —
6-% make_region's non-convex regions are star-shaped about it.
3+% non-convex). Along each random line it finds *every* crossing of the polygon
4+% boundary, so concavities are handled correctly (unlike the convex-only chord
5+% in hit_and_run.m). Starts at the origin — make_region's non-convex regions
6+% are star-shaped about it.
7+%
8+% LOCAL (default false) selects how the step samples along that line:
9+% false — sample uniformly across the *union* of all in-region segments the
10+% line makes (standard hit-and-run; uniform on the whole region).
11+% true — sample only within the single in-region segment that contains the
12+% current point (a local walk that can't jump across a concavity;
13+% it does not sample the region uniformly).
714 %
815 % This runs in the numbl interpreter (sort + point-in-polygon tests don't
916 % JIT), so it's used only for the non-convex demo at modest N.
17+if nargin < 5 || isempty(local)
18+ local = false;
19+end
1020 nv = numel(vx);
1121 px = 0;
1222 py = 0;
@@ -43,30 +53,47 @@ for step = 1:total
4353 end
4454 ts = sort(ts(1:m));
4555
46- % In-region intervals are consecutive crossings whose midpoint is inside.
47- totalLen = 0;
48- for k = 1:m - 1
49- tm = (ts(k) + ts(k + 1)) / 2;
50- if point_in_poly(px + tm * dx, py + tm * dy, vx, vy)
51- totalLen = totalLen + (ts(k + 1) - ts(k));
56+ if local
57+ % Local segment: the in-region interval straddling t = 0, i.e. bounded
58+ % by the nearest crossing on each side of the current (interior) point.
59+ tlo = -inf;
60+ thi = inf;
61+ for k = 1:m
62+ if ts(k) <= 0 && ts(k) > tlo
63+ tlo = ts(k);
64+ elseif ts(k) >= 0 && ts(k) < thi
65+ thi = ts(k);
66+ end
5267 end
53- end
54- if totalLen <= 0
55- continue
56- end
57-
58- % Pick a point uniformly across the union of in-region intervals.
59- u = totalLen * rand;
60- tpick = 0;
61- for k = 1:m - 1
62- tm = (ts(k) + ts(k + 1)) / 2;
63- if point_in_poly(px + tm * dx, py + tm * dy, vx, vy)
64- len = ts(k + 1) - ts(k);
65- if u <= len
66- tpick = ts(k) + u;
67- break
68+ if ~isfinite(tlo) || ~isfinite(thi) || thi <= tlo
69+ continue
70+ end
71+ tpick = tlo + (thi - tlo) * rand;
72+ else
73+ % In-region intervals are consecutive crossings whose midpoint is
74+ % inside; sample uniformly across their union.
75+ totalLen = 0;
76+ for k = 1:m - 1
77+ tm = (ts(k) + ts(k + 1)) / 2;
78+ if point_in_poly(px + tm * dx, py + tm * dy, vx, vy)
79+ totalLen = totalLen + (ts(k + 1) - ts(k));
80+ end
81+ end
82+ if totalLen <= 0
83+ continue
84+ end
85+ u = totalLen * rand;
86+ tpick = 0;
87+ for k = 1:m - 1
88+ tm = (ts(k) + ts(k + 1)) / 2;
89+ if point_in_poly(px + tm * dx, py + tm * dy, vx, vy)
90+ len = ts(k + 1) - ts(k);
91+ if u <= len
92+ tpick = ts(k) + u;
93+ break
94+ end
95+ u = u - len;
6896 end
69- u = u - len;
7097 end
7198 end
7299 px = px + tpick * dx;
hitandrun_sampler.mmodified+18−12View file
@@ -25,7 +25,7 @@ end
2525 % Build a region (convex by default) and draw N uniform samples from it.
2626 convex = true;
2727 [vx, vy] = make_region(convex);
28-[sx, sy] = sample_region(vx, vy, N, convex);
28+[sx, sy] = sample_region(vx, vy, N, convex, false);
2929
3030 % The figure app is the prebuilt single-file page, relative to the project root
3131 % (the current working directory when a top-level script is run).
@@ -39,14 +39,17 @@ uihtml(gl, 'HTMLSource', html, 'Data', pack_data(vx, vy, sx, sy, N, convex), ...
3939 end
4040
4141 function on_event(src, ev)
42-% Figure -> script. Two requests the controls send:
43-% 'resample' {n, x, y, convex} -> draw n fresh samples in the given region;
44-% reply with a 'samples' event.
45-% 'newRegion' {n, convex} -> build a new region (convex or non-convex),
46-% draw n samples; reply with a 'data' event.
42+% Figure -> script. Two requests the controls send (both carry the region type
43+% `convex` and, for non-convex regions, the `local` sampling mode):
44+% 'resample' {n, x, y, convex, local} -> draw n fresh samples in the given
45+% region; reply with a 'samples' event.
46+% 'newRegion' {n, convex, local} -> build a new region (convex or
47+% non-convex), draw n samples; reply with a
48+% 'data' event.
4749 d = ev.HTMLEventData;
4850 n = 10000;
4951 convex = true;
52+local = false;
5053 if isstruct(d)
5154 if isfield(d, 'n')
5255 n = max(1, round(d.n));
@@ -54,30 +57,33 @@ if isstruct(d)
5457 if isfield(d, 'convex')
5558 convex = logical(d.convex);
5659 end
60+ if isfield(d, 'local')
61+ local = logical(d.local);
62+ end
5763 end
5864 switch ev.HTMLEventName
5965 case 'resample'
6066 vx = d.x(:);
6167 vy = d.y(:);
62- [sx, sy] = sample_region(vx, vy, n, convex);
68+ [sx, sy] = sample_region(vx, vy, n, convex, local);
6369 sendEventToHTMLSource(src, 'samples', pack_samples(sx, sy, n));
6470 case 'newRegion'
6571 [vx, vy] = make_region(convex);
66- [sx, sy] = sample_region(vx, vy, n, convex);
72+ [sx, sy] = sample_region(vx, vy, n, convex, local);
6773 sendEventToHTMLSource(src, 'data', pack_data(vx, vy, sx, sy, n, convex));
6874 end
6975 end
7076
71-function [sx, sy] = sample_region(vx, vy, n, convex)
77+function [sx, sy] = sample_region(vx, vy, n, convex, local)
7278 % Convex regions use the fast JIT-compiled chord sampler; non-convex regions
73-% use the general multi-segment sampler.
79+% use the general sampler, in union (default) or local-segment mode.
7480 tic;
7581 if convex
7682 [sx, sy] = hit_and_run(vx, vy, n, 50);
7783 else
78- [sx, sy] = hit_and_run_general(vx, vy, n, 50);
84+ [sx, sy] = hit_and_run_general(vx, vy, n, 50, local);
7985 end
80-fprintf('sample_region: N=%d convex=%d in %.3f s\n', n, convex, toc);
86+fprintf('sample_region: N=%d convex=%d local=%d in %.3f s\n', n, convex, local, toc);
8187 end
8288
8389 function data = pack_data(vx, vy, sx, sy, N, convex)