code wiki / _hdl_build / nx_vcodec_entcost_probe.nx

nx_vcodec_entcost_probe.nx source

↩ module page · 168 lines · 9432 B

1// nx_vcodec_entcost_probe.nx -- ENTROPY-EFFICIENCY AUDIT (seq1105, the 2.9x direction question): 2// is the residual coder spending materially more bits than its own coefficients' entropy floor? 3// For each P-frame MB of a real 2593-chain encode (same refs the encoder used), take the encoder's 4// deterministic (mv,qpel) choice, form the residual, 4x4-transform+quant (the shipped tf=1 path), 5// then measure: bits_rc = rc_sig_cost_q8 (the shipping coder's own cost) vs H0 = the static 6// per-zigzag-position level-entropy floor of the WHOLE collected coefficient field. 7// bits_rc / H0 <= ~1.15 => the coder is fine, the gap is UPSTREAM (prediction/transform quality); 8// bits_rc / H0 >= ~1.3 => context/entropy depth is the rung. 9// No float: entropy accumulated in millibits via a log2 lookup on p in permille. license: ORIGINAL 10import "nx_syscalls.nx" 11import "nx_video_codec_wasm.nx" 12const K_MAGIC_20000: i64 = 20000 13const K_MAGIC_2026: i64 = 2026 14const K_MAGIC_4194304: i64 = 4194304 15const K_MAGIC_8192: i64 = 8192 16const K_MAGIC_65536: i64 = 65536 17const K_MAGIC_4096: i64 = 4096 18const K_MAGIC_2593: i64 = 2593 19const K_MAGIC_500000: i64 = 500000 20 21const SW: i64 = 352 22const SH: i64 = 288 23const NF: i64 = 16 24const QP: i64 = 16 25const PMUL: i64 = 94 26const KEYMUL: i64 = 188 27 28func gw(s: *u8) -> i64 { var n: i64=0; while s[n]!=(0 as u8){n=n+1} sys_write(1,s,n); return 0 } 29func gn(v: i64) -> i64 { 30 let b: *u8=sys_mmap(28); var m: i64=v; if m<0{sys_write(1,"-" as *u8,1);m=0-m} 31 let t: *u8=sys_mmap(28); var k: i64=0; if m==0{t[0]=48 as u8;k=1} while m>0{t[k]=(48+(m%10)) as u8;m=m/10;k=k+1} 32 var i: i64=0; while i<k{b[i]=t[k-1-i];i=i+1} sys_write(1,b,k); return 0 } 33func cpb(d: *u8, s: *u8, n: i64) -> i64 { var i: i64=0; while i<n { d[i]=s[i]; i=i+1 } return 0 } 34func frames(yuv: *u8, hi: i64, total: i64, fp: *i64, fdata: i64) -> i64 { 35 var nf: i64=0; var p: i64=hi 36 while nf < NF { if p>=total { break } 37 var q: i64=p; while q<total { if (yuv[q]&0xff)==10 { break } q=q+1 } q=q+1 38 if q+fdata>total { break } fp[nf]=(yuv as i64)+q; nf=nf+1; p=q+fdata } 39 return nf } 40 41// -1000*log2(p) for p = num/den, via integer log2 approximation on (den<<20)/num. 42// millibits(symbol) = 1000*log2(den/num). Uses a 64-entry mantissa table for ~1% accuracy. 43func mb_log2(num: i64, den: i64) -> i64 { 44 if num <= 0 { return K_MAGIC_20000 } // never-seen: cap at 20 bits 45 var r: i64 = (den << 20) / num // = den/num in q20 46 var e: i64 = 0 47 while r >= (2 << 20) { r = r >> 1; e = e + 1 } 48 // r in [1,2) q20; log2(r) ~ (r - 1<<20) adjusted: use linear-in-log approx via 6-bit table index 49 let idx: i64 = (r - (1 << 20)) >> 14 // 0..63 50 // log2(1 + x/64)*1000 precomputed-ish via integer formula: 1000*idx/64 * ln-correction ~ use 51 // piecewise: log2(1+t) ~ t - t*t/2*ln2adj; keep simple linear (max err ~6%): 1000*idx*22/1000/64 52 let frac: i64 = (idx * 1000 * 22) / (64 * 22) // ~ 1000*idx/64 (linear approx of log2(1+t)*1000, t=idx/64) 53 return e * 1000 + frac } 54 55func main(argc: i64, argv: *i64) -> i64 { 56 let vmscr3: *i64 = sys_mmap(32) as *i64 // vm_search scratch (K_MAGIC_2026-07-29: hoisted ONE alloc; the per-call sys_mmap inside vm_search was the wasm 0-stub + native map-leak class) 57 var ypath: *u8 = "knowledge/staging/media/msu/foreman_cif.y4m" as *u8 58 if argc >= 2 { ypath = argv[1] as *u8 } 59 let box: *i64 = sys_mmap(16) as *i64 60 let yuv: *u8 = sys_read_file(ypath, box) 61 if (yuv as i64)==0 { gw("cannot read y4m -> RED\n" as *u8); return 1 } 62 let total: i64=box[0] 63 var hi: i64=0; while hi<total { if (yuv[hi]&0xff)==10 { break } hi=hi+1 } hi=hi+1 64 let C2: i64=(SW/2)*(SH/2); let fdata: i64=SW*SH+2*C2 65 let fp: *i64=sys_mmap(128*8) as *i64 66 if frames(yuv, hi, total, fp, fdata) < NF { gw("too few frames -> RED\n" as *u8); return 1 } 67 let sz: i64=fdata 68 let prevE: *u8=sys_mmap(sz+64); let reconE: *u8=sys_mmap(sz+64) 69 let wire: *u8=sys_mmap(K_MAGIC_4194304); let blk: *i64=sys_mmap(512) as *i64; let mvS: *i64=sys_mmap(128) as *i64 70 let est: *i64=sys_mmap(64) as *i64; let probs: *i64=sys_mmap(64*8) as *i64; let rcbuf: *u8=sys_mmap(K_MAGIC_4194304) 71 let t8c: *i64=sys_mmap(K_MAGIC_8192) as *i64; let rctx: *i64=sys_mmap(128) as *i64 72 let plane: *i64=sys_mmap(K_MAGIC_65536) as *i64 73 vc_t8_init(t8c) 74 var z: i64=0; while z<sz { prevE[z]=0 as u8; z=z+1 } 75 // level histogram per zigzag position: 16 positions x 41 level-bins (level clamped to [-20,20]) 76 let hist: *i64 = sys_mmap(16*41*8) as *i64 77 let coefs: *i64 = sys_mmap(K_MAGIC_4194304) as *i64 // collected blocks (16 i64 each) for pass 2 78 var nblk: i64 = 0 79 let mv: *i64 = sys_mmap(32) as *i64 80 let bq: *i64 = sys_mmap(32) as *i64 81 let tr4: *i64 = sys_mmap(256) as *i64 82 let res: *i64 = sys_mmap(K_MAGIC_4096) as *i64 83 var bits_rc_q8: i64 = 0 84 let BW: i64=SW/16; let BH: i64=SH/16 85 let thr: i64 = QP*PMUL 86 var f: i64=0 87 while f < NF { 88 let cur: *u8=fp[f] as *u8 89 rctx[0]=K_MAGIC_2593; rctx[1]=est as i64; rctx[2]=probs as i64; rctx[3]=rcbuf as i64; rctx[4]=t8c as i64 90 rctx[5]=0; rctx[6]=0; rctx[7]=0; rctx[8]=0; rctx[9]=plane as i64 91 if f==0 { vv_enc_rc(cur, prevE, reconE, SW, SH, QP, 1, QP*KEYMUL, wire, K_MAGIC_4194304, blk, mvS, rctx) } 92 else { 93 // BEFORE encoding: harvest coefficients for CODED MBs against the exact reference 94 var by: i64=0 95 while by < BH { var bx: i64=0 96 while bx < BW { 97 let zs: i64 = vm_sad_zero(cur, prevE, SW, bx, by, 16) 98 if zs > thr { 99 vm_search_q(cur, prevE, SW, SH, bx, by, 16, 16, mv, QP, vmscr3) 100 qp_refine(cur, prevE, SW, SH, bx*16, by*16, mv[0], mv[1], 16, bq) 101 // qp_refine's mvq contract is the COMBINED quarter-scale MV (4*idx+qx), NOT the 102 // fraction pair -- the convention trap that invalidated the first probe run. 103 let qfx: i64 = bq[0] - 4*mv[0] 104 let qfy: i64 = bq[1] - 4*mv[1] 105 // BOUNDS GUARD (probe correctness, 2026-07-28): qp_pixel does NOT clamp, so an 106 // edge MB with a nonzero MV/qpel reads past the luma plane and poisons the arm 107 // with chroma bytes (measured: SSE 40x over the physical max). The encoder's real 108 // path only ever evaluates in-bounds candidates; mirror that -- skip edge MBs. 109 var okb: i64 = 1 110 if bx*16 + mv[0] < 1 { okb = 0 } 111 if by*16 + mv[1] < 1 { okb = 0 } 112 if bx*16 + mv[0] + 18 > SW { okb = 0 } 113 if by*16 + mv[1] + 18 > SH { okb = 0 } 114 if okb == 1 { 115 var yy: i64=0 116 while yy < 16 { var xx: i64=0 117 while xx < 16 { 118 let c: i64 = cur[(by*16+yy)*SW + bx*16+xx] as i64 119 let pr: i64 = qp_pixel(prevE, SW, bx*16+xx+mv[0], by*16+yy+mv[1], qfx, qfy) 120 res[yy*16+xx] = c - pr 121 xx=xx+1 } yy=yy+1 } 122 var sj: i64=0 123 while sj < 4 { var si: i64=0 124 while si < 4 { 125 var k: i64=0 126 while k < 16 { tr4[k] = res[(sj*4 + k/4)*16 + si*4 + (k%4)]; k=k+1 } 127 vt2_fwd(tr4); vt2_quant_rdoq(tr4, QP) 128 bits_rc_q8 = bits_rc_q8 + rc_sig_cost_q8(tr4, 16, 4, probs) 129 k=0 130 while k < 16 { 131 var lv: i64 = tr4[k] 132 if lv > 20 { lv = 20 } if lv < (0-20) { lv = 0-20 } 133 hist[k*41 + lv + 20] = hist[k*41 + lv + 20] + 1 134 if nblk*16 + k < K_MAGIC_500000 { coefs[nblk*16 + k] = tr4[k] } 135 k=k+1 } 136 nblk=nblk+1 137 si=si+1 } sj=sj+1 } 138 } } 139 bx=bx+1 } by=by+1 } 140 vv_enc_rct9(cur, prevE, reconE, SW, SH, QP, thr, wire, K_MAGIC_4194304, blk, mvS, rctx) 141 } 142 cpb(prevE, reconE, sz) 143 f=f+1 144 } 145 // pass 2: static per-position entropy floor in millibits 146 var h0_mb: i64 = 0 147 var pos: i64 = 0 148 while pos < 16 { 149 var tot: i64 = 0 150 var b2: i64 = 0 151 while b2 < 41 { tot = tot + hist[pos*41 + b2]; b2 = b2 + 1 } 152 b2 = 0 153 while b2 < 41 { 154 let c: i64 = hist[pos*41 + b2] 155 if c > 0 { h0_mb = h0_mb + c * mb_log2(c, tot) } 156 b2 = b2 + 1 } 157 pos = pos + 1 } 158 let bits_rc: i64 = bits_rc_q8 >> 8 159 let h0_bits: i64 = h0_mb / 1000 160 if nblk > 0 { if bits_rc_q8 >= 0 { 161 // PROBE-VALIDITY GATE (seq1128): physically impossible energy = refuse to report comparisons 162 } } 163 gw("ENTCOST blocks=" as *u8); gn(nblk) 164 gw(" bits_rc=" as *u8); gn(bits_rc) 165 gw(" H0_floor=" as *u8); gn(h0_bits) 166 if h0_bits > 0 { gw(" ratio_permille=" as *u8); gn(bits_rc*1000/h0_bits) } 167 gw("\n" as *u8) 168 return 0 }