#!/usr/bin/env python3
# baseline_orbit_test.py
# QUESTION (baseline, before any perturbation): left to breathe with NOTHING done to it,
# does the SPATIAL field return to itself — is there a stable attractor it relaxes onto
# (=> old perturbations wash out, we effectively HAVE a baseline despite never resetting)
# or does it accumulate/drift forever (=> contaminated, no clean baseline)?
#
# The global telemetry shows a clean ~120-cycle coherence limit cycle. This asks whether
# the SPATIAL field (5561 coarse stream) is ALSO on a clean orbit, or churns underneath it.
#
# Read-only: subscribes to 5561, sends NOTHING to the daemon. No injection.
#
# Method — "does it return to itself":
#   Capture a long undisturbed stretch. For a reference frame at time t0, measure the
#   spatial distance to every later frame. If there's an attractor/orbit, distance should
#   DIP periodically (field comes back near where it was) — a recurrence signature.
#   If it's accumulating, distance grows monotonically and never returns.
#   Cross-check against the mass-leak: work on MEAN-SUBTRACTED, per-frame-normalized fields
#   so we test PATTERN return, not the slow density ramp.

import sys, io, struct, time
sys.stdout = io.TextIOWrapper(sys.stdout.buffer, encoding='utf-8', errors='replace')
import numpy as np
try:
    import zmq
except ImportError:
    print("pip install pyzmq --break-system-packages"); sys.exit(1)

TILES, CH = 32, 6
NVALS = TILES*TILES*CH
HDR = 16
FRAME_BYTES = HDR + NVALS*4
PORT = "tcp://localhost:5561"
CH_NAMES = ["rho","ux","uy","sxx","syy","sxy"]
N_CAPTURE = 1500          # ~150s at 100Hz => ~15000 cycles => >100 of the ~120-cycle global orbits
OUT_DATA = "baseline_orbit.npz"
OUT_REPORT = "baseline_orbit_report.txt"

def capture(n):
    ctx=zmq.Context(); s=ctx.socket(zmq.SUB)
    s.setsockopt(zmq.RCVTIMEO,5000); s.setsockopt_string(zmq.SUBSCRIBE,"")
    s.connect(PORT)
    frames=[]; cycles=[]
    print(f"Capturing {n} frames (read-only, NO injection)...", flush=True)
    t0=time.time()
    while len(frames)<n:
        try: buf=s.recv()
        except zmq.Again:
            print(f"  timeout after {len(frames)} — is 5561 streaming?"); break
        if len(buf)!=FRAME_BYTES or buf[:4]!=b"KGCF": continue
        cyc=struct.unpack_from("<I",buf,4)[0]
        vals=np.frombuffer(buf,dtype=np.float32,count=NVALS,offset=HDR).copy()
        frames.append(vals.reshape(TILES*TILES,CH)); cycles.append(cyc)
        if len(frames)%250==0: print(f"  {len(frames)}/{n} cycle={cyc}", flush=True)
    s.close(); ctx.term()
    print(f"Captured {len(frames)} in {time.time()-t0:.1f}s", flush=True)
    return np.array(frames), np.array(cycles)

def norm_frames(F):
    # per-frame, per-channel mean-subtract + scale -> tests PATTERN, removes mass-leak amplitude
    G=F.astype(np.float64).copy()  # (T, tiles, CH)
    mu=G.mean(axis=1,keepdims=True)              # per-frame per-channel mean
    G=G-mu
    sd=G.std(axis=1,keepdims=True); sd[sd<1e-12]=1.0
    G=G/sd
    return G.reshape(len(G),-1)                   # (T, tiles*CH)

def main():
    F,C=capture(N_CAPTURE)
    if len(F)<300:
        print("Too few frames; aborting."); sys.exit(1)
    dcyc=np.diff(C)
    print(f"cycle step median={np.median(dcyc):.0f} (expect ~10)", flush=True)

    G=norm_frames(F)         # (T, D)
    T=len(G)

    lines=[]; 
    def out(s): print(s,flush=True); lines.append(s)
    out("="*64)
    out("BASELINE ORBIT / RETURN TEST — undisturbed spatial field (5561)")
    out("="*64)
    out(f"frames={T}  cycles {C[0]}..{C[-1]} (~{C[-1]-C[0]} cycles)")

    # --- RECURRENCE: distance from several reference frames to all later frames ---
    # Look for periodic DIPS (field returns near a past state).
    refs=[0, T//4, T//2]
    out("")
    out("RETURN SIGNATURE (does distance-from-reference dip periodically?):")
    dip_periods=[]
    for r in refs:
        d=np.sqrt(np.mean((G-G[r])**2,axis=1))   # distance from ref to all frames
        # ignore the trivial zero at r; look forward
        fwd=d[r+5:]
        if len(fwd)<50: continue
        # find local minima (returns) and their spacing in frames
        mins=[]
        for i in range(2,len(fwd)-2):
            if fwd[i]<fwd[i-1] and fwd[i]<fwd[i+1] and fwd[i]<np.median(fwd)*0.85:
                mins.append(i)
        if len(mins)>=2:
            spac=np.diff(mins)
            # convert frame-spacing to cycle-spacing (~10 cyc/frame)
            cyc_per=np.median(spac)*np.median(dcyc)
            dip_periods.append(cyc_per)
            out(f"  ref@frame{r}: {len(mins)} returns, ~{cyc_per:.0f} cycles between returns, "
                f"min dist reached {fwd[mins].min():.3f} (start dist {d[r]:.3f}->{fwd.min():.3f})")
        else:
            out(f"  ref@frame{r}: NO periodic returns found (distance does not dip back)")

    # --- ACCUMULATION vs BOUNDED: does distance-from-start grow forever or stay bounded? ---
    d0=np.sqrt(np.mean((G-G[0])**2,axis=1))
    # linear trend of distance over time
    tt=np.arange(T)
    slope=np.polyfit(tt,d0,1)[0]
    early=d0[:T//5].mean(); late=d0[-T//5:].mean()
    out("")
    out("ACCUMULATION CHECK (pattern distance from start over time):")
    out(f"  early mean dist {early:.3f} | late mean dist {late:.3f} | trend slope {slope:+.2e}/frame")
    bounded = abs(late-early)/max(early,1e-9) < 0.15 and abs(slope) < (early/ (T*3))
    out(f"  {'BOUNDED — pattern distance stays flat: field orbits, does NOT accumulate' if bounded else 'GROWING — pattern distance trends up: field is drifting/accumulating'}")

    # --- global coherence-analog on coarse field for cross-check with telemetry ~120-cyc ---
    out("")
    if dip_periods:
        out(f"SPATIAL return period ~{np.median(dip_periods):.0f} cycles "
            f"(telemetry global coherence orbit was ~120 cyc — compare).")
    out("")
    out("=== BASELINE VERDICT ===")
    has_orbit = len(dip_periods)>0
    if has_orbit and bounded:
        out("  ATTRACTOR PRESENT: the undisturbed spatial field RETURNS to near past states")
        out("  on a stable period and does not accumulate. => It relaxes onto an attractor;")
        out("  old perturbations plausibly wash out; we effectively HAVE a baseline despite")
        out("  never resetting. The natural churn is a bounded orbit = a clean noise floor")
        out("  to run an on/off perturbation test against.")
    elif not has_orbit and not bounded:
        out("  ACCUMULATING / NON-RETURNING: the spatial field does not come back to itself")
        out("  and drifts. => Jason's contamination worry is supported: no clean baseline;")
        out("  history accumulates. A perturbation test needs each trial from a fresh daemon,")
        out("  or must measure deltas over windows short vs the drift.")
    else:
        out("  MIXED: orbit and accumulation signals disagree — report both, do not force a verdict.")
        out("  (e.g. bounded but no clean periodic return = wanders on a bounded manifold.)")

    np.savez(OUT_DATA, frames=F.astype(np.float32), cycles=C, d0=d0)
    with open(OUT_REPORT,"w") as f: f.write("\n".join(lines)+"\n")
    print(f"\nWrote {OUT_REPORT} and {OUT_DATA}. DONE.", flush=True)

if __name__=="__main__":
    main()