| """Frozen v3 experiment; ground truth is restricted to the audit path.""" |
| from pathlib import Path |
| import argparse,csv,json,time,sys,platform,hashlib |
| ROOT=Path(__file__).resolve().parents[1];sys.path.insert(0,str(ROOT)) |
| import numpy as np |
| from aureole import WorldMemory,freeze,draw,correct,exact_mse |
| from aureole.core import proposal_from_bound |
| from aureole.renderer import Scene,VisibilityPrior,receiver_grid,light_grid,unoccluded,physical_table |
| from aureole.certificates import CertificateMemory,visibility_certificate,enclosure |
| from aureole.innovation import prepare,eliminate,exact_risk_two |
|
|
|
|
| def iid(h,q,oracle,n,rng): |
| rows=np.flatnonzero(q.sum(1)>0);y=h.sum(1) |
| if not len(rows):return y,0 |
| cdf=np.minimum(np.cumsum(q[rows],1),1.0) |
| last=q.shape[1]-1-np.argmax(q[rows,::-1]>0,axis=1) |
| cdf[np.arange(q.shape[1])[None,:]>=last[:,None]]=1.0 |
| for _ in range(n): |
| j=(rng.random(len(rows))[:,None]>=cdf).sum(1) |
| f=oracle(rows,j) |
| y[rows]+=(f-h[rows,j])/q[rows,j,None]/n |
| return y,len(rows)*n |
|
|
|
|
| def audit_iid(truth,h,q,n=2): |
| r=truth-h |
| second=np.divide(r*r,q[...,None],out=np.zeros_like(r),where=q[...,None]>0).sum(1) |
| return np.maximum(second-r.sum(1)**2,0).mean(-1)/n |
|
|
|
|
| def summarize(rows,config): |
| summary=[];comparisons=[] |
| for phase in config['phases']: |
| for method in config['methods']: |
| rs=[r for r in rows if r['phase']==phase and r['method']==method] |
| summary.append({'phase':phase,'method':method, |
| 'expected_mse':float(np.mean([r['expected_mse'] for r in rs])), |
| 'observed_mse':float(np.mean([r['observed_mse'] for r in rs])), |
| 'rays':int(sum(r['rays'] for r in rs)), |
| 'rays_per_receiver':float(np.mean([r['rays']/r['receivers'] for r in rs])), |
| 'certified_fraction_before':float(np.mean([r['certified_fraction_before'] for r in rs])), |
| 'false_certificates':int(sum(r['false_certificates'] for r in rs)), |
| 'max_enclosure_violation':float(max(r['enclosure_violation'] for r in rs)), |
| 'cpu_median_ms':float(np.median([r['runtime_seconds'] for r in rs])*1000), |
| 'cpu_p99_ms':float(np.quantile([r['runtime_seconds'] for r in rs],.99)*1000), |
| 'memory_bytes':max(r['memory_bytes'] for r in rs)}) |
| for baseline in ('epoch_eliminate','v2_guarded','certificate_iid','constant_certificate','raw_importance'): |
| a=[];b=[] |
| for sid in config['scene_ids']: |
| a.append(np.mean([r['expected_mse'] for r in rows if r['phase']==phase and r['scene']==sid and r['method']=='certificate_eliminate'])) |
| b.append(np.mean([r['expected_mse'] for r in rows if r['phase']==phase and r['scene']==sid and r['method']==baseline])) |
| a,b=np.array(a),np.array(b) |
| if b.mean()<1e-20: |
| value=None;interval=None |
| else: |
| value=float(1-a.mean()/b.mean());rng=np.random.default_rng(9301) |
| indices=rng.integers(0,len(a),(10000,len(a))) |
| x=1-a[indices].mean(1)/np.maximum(b[indices].mean(1),1e-30) |
| interval=np.percentile(x,[2.5,97.5]).tolist() |
| comparisons.append({'phase':phase,'baseline':baseline,'contender':'certificate_eliminate', |
| 'expected_mse_reduction':value,'scene_bootstrap_95':interval}) |
| return summary,comparisons |
|
|
|
|
| def run(output): |
| protocol=ROOT/'experiments_innovation.json';cfg=json.loads(protocol.read_text()) |
| out=ROOT/output;out.mkdir(parents=True,exist_ok=True) |
| points=receiver_grid(*cfg['receiver_grid']);lights=light_grid(cfg['emitter_grid_side']) |
| H,W=cfg['receiver_grid'];vh,vw=cfg['viewport'];prior=VisibilityPrior(ROOT/'models/visibility_prior.npz') |
| schedule=[(p,k) for p,n in cfg['phases'].items() for k in range(n)] |
| rows=[];preparation=[];images={};started=time.perf_counter() |
| for sid in cfg['scene_ids']: |
| original=Scene.create(sid) |
| tic=time.perf_counter();p_all=prior(original.features(points[:,None,:],lights[None,:,:])) |
| preparation.append({'scene':sid,'prior_seconds':time.perf_counter()-tic}) |
| b0=unoccluded(points,lights);b1=unoccluded(points,lights,True,.7) |
| for seed in cfg['replicate_seeds']: |
| memory={} |
| for method in cfg['methods']: |
| if method in ('raw_importance','v2_guarded'):memory[method]=WorldMemory(len(points),len(lights),f'scene-{sid}') |
| else:memory[method]=CertificateMemory(points,lights,original.spheres,f'scene-{sid}', |
| 'epoch' if method=='epoch_eliminate' else 'unsafe' if method=='unsafe_eliminate' else 'margin') |
| rngs={m:np.random.default_rng(seed+sid*1000) for m in cfg['methods']} |
| for frame,(phase,k) in enumerate(schedule): |
| g=original.spheres.copy() |
| if phase in ('smooth_motion','jump'): |
| u=k+1 if phase=='smooth_motion' else cfg['phases']['smooth_motion'] |
| g[:,0]+=np.array([.004,-.003,.002])*u |
| g[:,1]+=np.array([-.002,.003,.002])*u |
| g[:,3]+=np.array([.0002,-.0001,.00015])*u |
| scene=Scene(sid,g) |
| if phase=='jump':scene=scene.changed() |
| x0=2+(k%4);y0=8 |
| ids=(np.arange(y0,y0+vh)[:,None]*W+np.arange(x0,x0+vw)[None,:]).ravel() |
| b=(b0 if phase in ('cold','warm','revisit') else b1)[ids] |
| |
| truth=physical_table(scene,points[ids],lights,b);target=truth.sum(1) |
| vis_truth=scene.visibility(points[ids,None,:],lights[None,:,:]) |
| for method in cfg['methods']: |
| tic=time.perf_counter();mem=memory[method];rng=rngs[method] |
| false=0;known=np.zeros(b.shape[:2],bool);values=np.zeros_like(known,float) |
| if method in ('raw_importance','v2_guarded'): |
| v=np.zeros_like(p_all[ids]) if method=='raw_importance' else mem.predict(ids,p_all[ids]) |
| q=proposal_from_bound(b,v,mem.trusted(ids),active=method=='v2_guarded') |
| snap=freeze(b*v[...,None],q);j=draw(snap,2,rng) |
| observed=scene.visibility(points[ids,None,:],lights[j]) |
| y=correct(snap,j,b[np.arange(len(ids))[:,None],j]*observed[...,None]) |
| if method=='v2_guarded':mem.commit(ids,j,observed,revise_on_conflict=True) |
| count=len(ids)*2;h=snap.control |
| else: |
| mem.begin_geometry(scene.spheres) |
| values,known=mem.lookup(ids) |
| p=np.full_like(p_all[ids],.5) if method=='constant_certificate' else p_all[ids] |
| score=proposal_from_bound(b,p,active=True) |
| h,q=prepare(b,p,values,known,score) |
| count_box=[0] |
| def oracle(rr,j): |
| |
| vv,mm=visibility_certificate(scene,points[ids[rr]],lights[j]) |
| mem.commit(ids[rr],j,vv,mm);count_box[0]+=len(rr) |
| return b[rr,j]*vv[:,None] |
| if method=='certificate_iid':y,_=iid(h,q,oracle,2,rng) |
| else:y,_=eliminate(h,q,oracle,2,rng) |
| count=count_box[0] |
| runtime=time.perf_counter()-tic |
| |
| false=int(np.count_nonzero(known & (values!=vis_truth))) |
| if method in ('raw_importance','v2_guarded'):risk=exact_mse(truth,snap,2) |
| else: |
| simulated=np.where(known[...,None],h,truth) |
| bias=(simulated-truth).sum(1) |
| risk=(audit_iid(simulated,h,q) if method=='certificate_iid' else exact_risk_two(simulated,h,q))+(bias*bias).mean(-1) |
| lo,hi=enclosure(b,values,known) |
| violation=float(max(np.max(lo-target),np.max(target-hi),0)) |
| rows.append({'scene':sid,'seed':seed,'phase':phase,'phase_frame':k,'frame':frame,'method':method, |
| 'expected_mse':float(risk.mean()),'observed_mse':float(np.mean((y-target)**2)), |
| 'rays':count,'receivers':len(ids),'runtime_seconds':runtime, |
| 'memory_bytes':0 if method=='raw_importance' else mem.nbytes, |
| 'certified_fraction_before':float(known.mean()),'false_certificates':false, |
| 'enclosure_violation':violation,'enclosure_mean_width':float(np.mean(hi-lo)), |
| 'negative_channel_fraction':float(np.mean(y<0))}) |
| if sid==cfg['scene_ids'][0] and seed==cfg['replicate_seeds'][0] and k==0: |
| images[f'{phase}_{method}']=y.reshape(vh,vw,3);images[f'{phase}_reference']=target.reshape(vh,vw,3) |
| print(f'completed scene {sid}; {len(rows)} frame-method records',flush=True) |
| with (out/'innovation_raw.csv').open('w',newline='') as f: |
| writer=csv.DictWriter(f,fieldnames=list(rows[0]));writer.writeheader();writer.writerows(rows) |
| summary,comparisons=summarize(rows,cfg) |
| report={'protocol':cfg,'protocol_sha256':hashlib.sha256(protocol.read_bytes()).hexdigest(), |
| 'elapsed_seconds':time.perf_counter()-started,'environment':{'python':platform.python_version(),'numpy':np.__version__,'device':'CPU'}, |
| 'records':len(rows),'online_segment_queries':sum(r['rays'] for r in rows), |
| 'summary':summary,'comparisons':comparisons,'prior_preparation':preparation, |
| 'limitations':['Not matched time or memory; certificate queries perform additional clearance arithmetic.', |
| 'Finite direct-light renderer with 36 emitters and fixed receiver points.', |
| 'Current geometry is authoritative; this is not a solution to truly unobservable geometry changes.', |
| 'CPU batch timings exclude common prior preparation, analytic bound preparation and offline reference audits.', |
| 'Unsafe ablation deliberately violates validity; its expected MSE includes bias.', |
| 'No retraining, GPU benchmark, commercial upscaler or unified SR/RR/FG validation.']} |
| (out/'innovation_report.json').write_text(json.dumps(report,indent=2)+'\n') |
| np.savez_compressed(out/'innovation_frames.npz',**images) |
| print(json.dumps({'records':len(rows),'queries':report['online_segment_queries'],'primary':[x for x in comparisons if x['phase']=='smooth_motion' and x['baseline']=='epoch_eliminate']},indent=2)) |
|
|
|
|
| if __name__=='__main__': |
| p=argparse.ArgumentParser();p.add_argument('--output',default='innovation_reproduced');run(p.parse_args().output) |
|
|
