"""Self-contained synthetic education data. No network, recipe or product claims.
Usage: python3 depth_analysis.py spc (or doe/pe/reliability/thermal/signal/ssd/pm).
Keep this script and the selected CSV(s) in one separate folder.
"""
import csv,json,math,statistics as st,sys
from collections import defaultdict
from pathlib import Path
ROOT=Path(__file__).resolve().parent

def read(name):
 with (ROOT/(name+'.csv')).open(newline='',encoding='utf-8') as f: rows=list(csv.DictReader(f))
 if not rows: raise ValueError('empty data')
 return rows

def num(row,key):
 v=float(row[key])
 if not math.isfinite(v): raise ValueError('nonfinite '+key)
 return v

def unique(rows,keys):
 seen=set()
 for r in rows:
  k=tuple(r[x] for x in keys)
  if k in seen: raise ValueError('duplicate key '+str(k))
  seen.add(k)

def spc():
 rows=read('spc');unique(rows,['order'])
 if any(int(r['order'])!=i+1 for i,r in enumerate(rows)): raise ValueError('time order must be contiguous')
 if [r['phase'] for r in rows]!=['baseline']*12+['monitor']*6: raise ValueError('baseline/monitor layout')
 b=[num(r,'value_nm') for r in rows[:12]];mr=[abs(y-x) for x,y in zip(b,b[1:])]
 center=st.mean(b);m=st.mean(mr);sigma=m/1.128
 if sigma<=0: raise ValueError('no baseline dispersion')
 lo,hi=center-3*sigma,center+3*sigma;mrhi=3.267*m
 previous=b[-1];result=[]
 for r in rows[12:]:
  x=num(r,'value_nm');move=abs(x-previous);previous=x
  result.append({'order':int(r['order']),'value_nm':x,'moving_range':move,'I_signal':x<lo or x>hi,'MR_signal':move>mrhi,'example_spec_pass':49<=x<=51})
 shifted=read('spc-mean-shift');unique(shifted,['order']);shift_result=[];previous=b[-1]
 if [int(r['order']) for r in shifted]!=list(range(13,19)):raise ValueError('shift time order')
 for r in shifted:
  x=num(r,'value_nm');move=abs(x-previous);previous=x
  shift_result.append({'order':int(r['order']),'value_nm':x,'I_signal':x<lo or x>hi,'MR_signal':move>mrhi})
 return {'baseline_mean':center,'baseline_MR_mean':m,'I_LCL':lo,'I_UCL':hi,'MR_UCL':mrhi,'monitor':result,'separate_mean_shift_monitor':shift_result}

def doe():
 rows=read('doe');unique(rows,['run']);unique(rows,['A','B','replicate'])
 groups=defaultdict(list)
 for r in rows:
  a,b=int(r['A']),int(r['B'])
  if a not in [-1,1] or b not in [-1,1]:raise ValueError('coded factors must be -1/+1')
  groups[a,b].append(num(r,'response_au'))
 if set(groups)!={(-1,-1),(-1,1),(1,-1),(1,1)} or {len(v) for v in groups.values()}!={3}:raise ValueError('balanced full factorial required')
 means={k:st.mean(v) for k,v in groups.items()};n=len(rows);grand=st.mean([num(r,'response_au') for r in rows]);effects={}
 for key in ['A','B','AB']:
  signs=lambda r:int(r['A'])*int(r['B']) if key=='AB' else int(r[key])
  effects[key]=st.mean([num(r,'response_au') for r in rows if signs(r)==1])-st.mean([num(r,'response_au') for r in rows if signs(r)==-1])
 sse=sum((x-means[k])**2 for k,v in groups.items() for x in v);df=n-4;mse=sse/df
 confirm=read('doe-confirmation');unique(confirm,['run']);c=defaultdict(list)
 for r in confirm:
  if int(r['B'])!=1 or int(r['A']) not in [-1,1]:raise ValueError('confirmation factors')
  c[int(r['A'])].append(num(r,'response_au'))
 if set(c)!={-1,1} or min(map(len,c.values()))<2:raise ValueError('confirmation groups')
 return {'cell_means':{str(k):v for k,v in means.items()},'grand_mean':grand,'factorial_effects':effects,'A_simple_B_minus':means[1,-1]-means[-1,-1],'A_simple_B_plus':means[1,1]-means[-1,1],'difference_of_differences':2*effects['AB'],'pure_error_SSE':sse,'pure_error_df':df,'MSE':mse,'effect_SE':math.sqrt(4*mse/n),'confirmation_A_difference_at_B_plus':st.mean(c[1])-st.mean(c[-1])}

def pe():
 rows=read('pe');unique(rows,['unit','site']);unique(rows,['unit','order']);units=defaultdict(dict);sites=defaultdict(list);orders=defaultdict(list)
 if any(len({r[k] for r in rows})!=1 for k in ['product','program','tester','temperature_C']):raise ValueError('unmatched conditions')
 for r in rows:
  if r['site'] not in ['S1','S2'] or r['fail'] not in ['0','1'] or r['order'] not in ['1','2']:raise ValueError('invalid site/result/order')
  units[r['unit']][r['site']]=int(r['fail']);sites[r['site']].append(int(r['fail']));orders[r['order']].append(int(r['fail']))
 if any(set(v)!={'S1','S2'} for v in units.values()):raise ValueError('unpaired units')
 patterns=defaultdict(list)
 for u,v in units.items():patterns[f"S1={v['S1']},S2={v['S2']}"].append(u)
 return {'unique_units':len(units),'attempts':len(rows),'site_fail_fraction':{k:sum(v)/len(v) for k,v in sites.items()},'order_fail_fraction':{k:sum(v)/len(v) for k,v in orders.items()},'paired_patterns':dict(patterns)}

def reliability():
 rows=read('reliability');unique(rows,['unit']);obs=[]
 for r in rows:
  t=num(r,'time_h')
  if t<=0 or r['event'] not in ['0','1']:raise ValueError('time/event')
  obs.append((t,int(r['event'])))
 survival=1.;steps=[]
 for t in sorted({t for t,e in obs}):
  n=sum(x>=t for x,e in obs);d=sum(x==t and e==1 for x,e in obs);c=sum(x==t and e==0 for x,e in obs)
  if d:survival*=1-d/n
  steps.append({'time_h':t,'at_risk_before':n,'failures':d,'censored':c,'survival_estimate':survival})
 return {'units':len(obs),'failures':sum(e for t,e in obs),'censored':sum(not e for t,e in obs),'KM_steps':steps,'extrapolation':None}

def thermal():
 rows=read('thermal');unique(rows,['case']);out=[]
 for r in rows:
  p=num(r,'power_W');a=num(r,'junction_case_K_W');b=num(r,'case_ambient_K_W')
  if p<0 or a<=0 or b<=0:raise ValueError('power/resistance')
  predicted=num(r,'ambient_C')+p*(a+b)
  out.append({'case':r['case'],'equivalent_K_W':a+b,'single_path_junction_C':predicted,'observed_minus_model_C':num(r,'observed_junction_C')-predicted})
 return out

def signal():
 rows=read('signal');unique(rows,['case']);out=[]
 for r in rows:
  ui=num(r,'UI_ps');w=num(r,'eye_width_ps');h=num(r,'eye_height_mV');rw=num(r,'width_requirement_ps');rh=num(r,'height_requirement_mV')
  if ui<=0 or not 0<=w<=ui or min(h,rw,rh)<0:raise ValueError('eye quantities')
  out.append({'case':r['case'],'width_fraction_UI':w/ui,'time_margin_ps':w-rw,'voltage_margin_mV':h-rh})
 return out

def ssd():
 rows=read('ssd');unique(rows,['scenario','request']);groups=defaultdict(list)
 for r in rows:
  v=num(r,'latency_ms')
  if v<0:raise ValueError('negative latency')
  groups[r['scenario']].append(v)
 def percentile(xs,p):return xs[math.ceil(p*len(xs))-1]
 return {k:{'requests':len(v),'mean_ms':st.mean(v),'p50_ms':percentile(sorted(v),.5),'p95_ms':percentile(sorted(v),.95),'p99_ms':percentile(sorted(v),.99),'max_ms':max(v)} for k,v in groups.items()}

def pm():
 rows=read('pm');unique(rows,['wafer']);groups=defaultdict(list)
 for r in rows:
  if r['phase'] not in ['before','after']:raise ValueError('phase')
  if num(r,'particle_count')<0:raise ValueError('negative count')
  groups[r['phase']].append(r)
 if set(groups)!={'before','after'}:raise ValueError('phase missing')
 return {k:{'wafers':len(v),'mean_pressure_error_au':st.mean(num(r,'pressure_error_au') for r in v),'mean_particles':st.mean(num(r,'particle_count') for r in v)} for k,v in groups.items()}

MODES={f.__name__:f for f in [spc,doe,pe,reliability,thermal,signal,ssd,pm]}
if __name__=='__main__':
 if len(sys.argv)!=2 or sys.argv[1] not in MODES:raise SystemExit('Choose: '+', '.join(MODES))
 print(json.dumps(MODES[sys.argv[1]](),ensure_ascii=False,indent=2,allow_nan=False))