nx_motion.nx source

1// nx_motion.nx -- frame-to-frame motion (Lucas-Kanade optical flow). 2// 3// Brightness-constancy + small-motion assumption: 4// I(x, y, t) = I(x + u, y + v, t + 1) 5// Taylor expand: 6// I_x * u + I_y * v + I_t = 0 7// Over an N-pixel window this is N equations in 2 unknowns; solve 8// via 2x2 least-squares (normal equations + Cramer's rule). 9// 10// What this unlocks: 11// + motion detection (frame diff) 12// + tracking (per-pixel flow vectors) 13// + camera-motion estimation (mean flow direction) 14// + scene-cut detection (flow magnitude spike) 15// + activity recognition substrate 16// 17// Returns flow in Q8 sub-pixel precision (1 unit = 1/256 pixel). 18// 5x5 integration window default. 19// 20// genealogy_id: lucas_kanade_1981 + horn_schunck_1981 21// lineage_id: gradient_based_optical_flow + 2x2_least_squares 22 23// nx_safety_envelope: 24// intended_use: AUTO_APPLIED -- primitive-specific tuning queued 25// sil_target: SIL1 26// evidence: [bulk_applied_2026-05-16, see-file-comment-for-detail] 27// verdict: NOT_YET_EVALUATED 28 29import "syscalls.nx" 30import "nx_image.nx" 31 32const NX_MOTION_Q: i64 = 256 // Q8 sub-pixel 33const NX_MOTION_WIN: i64 = 5 // 5x5 integration 34 35// ===== Frame difference ================================================= 36// 37// Absolute difference between two same-size grayscale images. 38 39func nx_motion_frame_diff(a: *Image, b: *Image) -> *Image { 40 let w: i64 = a.width 41 let h: i64 = a.height 42 let dst: *Image = nx_image_alloc(w, h, 1) 43 var y: i64 = 0 44 while y < h { 45 var x: i64 = 0 46 while x < w { 47 let av: i64 = nx_image_get(a, x, y, 0) 48 let bv: i64 = nx_image_get(b, x, y, 0) 49 var d: i64 = av - bv 50 if d < 0 { d = -d } 51 if d > 255 { d = 255 } 52 nx_image_set(dst, x, y, 0, d) 53 x = x + 1 54 } 55 y = y + 1 56 } 57 return dst 58} 59 60// ===== Temporal derivative ============================================ 61// 62// I_t = frame2 - frame1. Signed output in ImageS64. 63 64func nx_motion_temporal_grad(frame1: *Image, frame2: *Image) -> *ImageS64 { 65 let w: i64 = frame1.width 66 let h: i64 = frame1.height 67 let dst: *ImageS64 = nx_image_s64_alloc(w, h) 68 var y: i64 = 0 69 while y < h { 70 var x: i64 = 0 71 while x < w { 72 let v: i64 = nx_image_get(frame2, x, y, 0) - nx_image_get(frame1, x, y, 0) 73 nx_image_s64_set(dst, x, y, v) 74 x = x + 1 75 } 76 y = y + 1 77 } 78 return dst 79} 80 81// ===== Lucas-Kanade per-pixel flow ===================================== 82// 83// Computes flow (u, v) at each pixel by solving the 2x2 normal 84// equations over a 5x5 window of spatial+temporal gradients. 85// 86// Pixels with insufficient texture (det close to zero) get u=v=0. 87// 88// Returns two ImageS64 maps: flow_u (Q8) and flow_v (Q8). 89 90struct FlowField { 91 u: *ImageS64, // Q8 horizontal flow 92 v: *ImageS64, // Q8 vertical flow 93} 94 95func nx_motion_lucas_kanade(frame1: *Image, frame2: *Image) -> *FlowField { 96 let w: i64 = frame1.width 97 let h: i64 = frame1.height 98 99 // Compute spatial + temporal gradients on frame1. 100 let Ix: *ImageS64 = nx_image_sobel_x(frame1) 101 let Iy: *ImageS64 = nx_image_sobel_y(frame1) 102 let It: *ImageS64 = nx_motion_temporal_grad(frame1, frame2) 103 104 let flow_u: *ImageS64 = nx_image_s64_alloc(w, h) 105 let flow_v: *ImageS64 = nx_image_s64_alloc(w, h) 106 107 let half: i64 = NX_MOTION_WIN / 2 108 var y: i64 = half 109 while y < h - half { 110 var x: i64 = half 111 while x < w - half { 112 // Accumulate over 5x5 window. 113 var sxx: i64 = 0 114 var sxy: i64 = 0 115 var syy: i64 = 0 116 var sxt: i64 = 0 117 var syt: i64 = 0 118 var dy: i64 = -half 119 while dy <= half { 120 var dx: i64 = -half 121 while dx <= half { 122 let ix: i64 = nx_image_s64_get(Ix, x + dx, y + dy) 123 let iy: i64 = nx_image_s64_get(Iy, x + dx, y + dy) 124 let it: i64 = nx_image_s64_get(It, x + dx, y + dy) 125 sxx = sxx + ix * ix 126 sxy = sxy + ix * iy 127 syy = syy + iy * iy 128 sxt = sxt + ix * it 129 syt = syt + iy * it 130 dx = dx + 1 131 } 132 dy = dy + 1 133 } 134 // 2x2 system: [sxx sxy] [u] [-sxt] 135 // [sxy syy] [v] = [-syt] 136 // det = sxx * syy - sxy^2 137 let det: i64 = sxx * syy - sxy * sxy 138 if det != 0 { 139 // u = (syy * (-sxt) - sxy * (-syt)) / det 140 // v = (sxx * (-syt) - sxy * (-sxt)) / det 141 // Scale up by Q for sub-pixel precision. 142 let u_num: i64 = (sxy * syt - syy * sxt) * NX_MOTION_Q 143 let v_num: i64 = (sxy * sxt - sxx * syt) * NX_MOTION_Q 144 nx_image_s64_set(flow_u, x, y, u_num / det) 145 nx_image_s64_set(flow_v, x, y, v_num / det) 146 } 147 x = x + 1 148 } 149 y = y + 1 150 } 151 152 let field: *FlowField = (sys_mmap(16)) as *FlowField 153 field.u = flow_u 154 field.v = flow_v 155 return field 156} 157 158// ===== Mean flow ======================================================== 159// 160// Mean (u, v) over high-confidence pixels (those where |flow| >= threshold). 161// Useful for global camera-motion estimate. 162 163func nx_motion_mean_flow(flow: *FlowField, mag_threshold: i64, 164 out_u: *i64, out_v: *i64) -> i64 { 165 let w: i64 = flow.u.width 166 let h: i64 = flow.u.height 167 var sum_u: i64 = 0 168 var sum_v: i64 = 0 169 var n: i64 = 0 170 var y: i64 = 0 171 while y < h { 172 var x: i64 = 0 173 while x < w { 174 let uv: i64 = nx_image_s64_get(flow.u, x, y) 175 let vv: i64 = nx_image_s64_get(flow.v, x, y) 176 var au: i64 = uv 177 if au < 0 { au = -au } 178 var av: i64 = vv 179 if av < 0 { av = -av } 180 if au + av >= mag_threshold { 181 sum_u = sum_u + uv 182 sum_v = sum_v + vv 183 n = n + 1 184 } 185 x = x + 1 186 } 187 y = y + 1 188 } 189 if n == 0 { 190 out_u[0] = 0 191 out_v[0] = 0 192 return 0 193 } 194 out_u[0] = sum_u / n 195 out_v[0] = sum_v / n 196 return n 197}