GPT Second Script

  
import numpy as np
import matplotlib.pyplot as plt
import wave
import scipy.signal as signal
from scipy.stats import entropy

filename = "clean.wav"   # change if needed

# -----------------------------
# Load Audio
# -----------------------------
with wave.open(filename, 'rb') as wf:
    framerate = wf.getframerate()
    nframes = wf.getnframes()
    audio = wf.readframes(nframes)

audio = np.frombuffer(audio, dtype=np.int16)

# Convert to mono if stereo
if len(audio.shape) > 1:
    audio = audio[:,0]

audio = audio.astype(float)

# -----------------------------
# STFT
# -----------------------------
f, t, Zxx = signal.stft(audio, fs=framerate, nperseg=4096)
magnitude = np.abs(Zxx)

# -----------------------------
# Helper: Band Energy
# -----------------------------
def band_energy(center_freq, width=10):
    idx = np.where((f >= center_freq - width) & (f <= center_freq + width))[0]
    return magnitude[idx].mean(axis=0)

# =====================================================
# OUTPUT 1: Original 3 Frequencies
# =====================================================
freqs_original = [1492, 1987, 2255]
energies_original = [band_energy(freq) for freq in freqs_original]

plt.figure()
for i, freq in enumerate(freqs_original):
    plt.plot(t, energies_original[i], label=f"{freq} Hz")
plt.title("Original 3-Band Energy")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
plt.savefig("output_1_original3.png")
plt.close()

# =====================================================
# OUTPUT 2: Auto-Detect 3 Most Dominant Narrow Bands
# =====================================================
mean_spectrum = magnitude.mean(axis=1)
top_indices = np.argsort(mean_spectrum)[-3:]
auto_freqs = f[top_indices]
auto_energies = [band_energy(freq) for freq in auto_freqs]

plt.figure()
for i, freq in enumerate(auto_freqs):
    plt.plot(t, auto_energies[i], label=f"{round(freq,2)} Hz")
plt.title("Auto-Detected Top 3 Bands")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.legend()
plt.savefig("output_2_auto_detected.png")
plt.close()

# =====================================================
# OUTPUT 3: Mutual Exclusivity
# =====================================================
stacked = np.vstack(auto_energies)
dominance = np.argmax(stacked, axis=0)

exclusivity_score = np.mean([
    np.sum(dominance == i) / len(dominance)
    for i in range(3)
])

with open("output_3_mutual_exclusivity.txt", "w") as f_out:
    f_out.write(f"Dominance Distribution: {np.bincount(dominance)}\n")
    f_out.write(f"Exclusivity Score (ideal ~0.33): {exclusivity_score}\n")

# =====================================================
# OUTPUT 4: Dominance Rotation
# =====================================================
transitions = np.sum(dominance[:-1] != dominance[1:])
transition_rate = transitions / len(dominance)

with open("output_4_dominance_rotation.txt", "w") as f_out:
    f_out.write(f"Total transitions: {transitions}\n")
    f_out.write(f"Transition rate: {transition_rate}\n")

# =====================================================
# OUTPUT 5: Correlation Matrix
# =====================================================
corr_matrix = np.corrcoef(auto_energies)

with open("output_5_correlation.txt", "w") as f_out:
    f_out.write(str(corr_matrix))

# =====================================================
# OUTPUT 6: Quantization / Block Structure
# =====================================================
block_lengths = []
current = dominance[0]
length = 1

for d in dominance[1:]:
    if d == current:
        length += 1
    else:
        block_lengths.append(length)
        current = d
        length = 1
block_lengths.append(length)

with open("output_6_blocks.txt", "w") as f_out:
    f_out.write(f"Mean block length: {np.mean(block_lengths)}\n")
    f_out.write(f"Std block length: {np.std(block_lengths)}\n")

# =====================================================
# OUTPUT 7: 725 CLUE CHECK
# =====================================================
energy_725 = band_energy(725)

plt.figure()
plt.plot(t, energy_725)
plt.title("725 Hz Energy Over Time")
plt.xlabel("Time (s)")
plt.ylabel("Amplitude")
plt.savefig("output_7_725.png")
plt.close()

with open("output_7_725_stats.txt", "w") as f_out:
    f_out.write(f"Mean energy at 725 Hz: {np.mean(energy_725)}\n")
    f_out.write(f"Std energy at 725 Hz: {np.std(energy_725)}\n")

# =====================================================
# OUTPUT 8: Ternary Decode With 725 Timing Offset
# =====================================================

# Convert dominance states into symbol stream
symbols_full = dominance.copy()

# Interpret 725 as timing offset in frames
offset = 725

if offset >= len(symbols_full):
    raise ValueError("Offset 725 exceeds symbol length.")

# Extract symbols starting at offset
symbols = symbols_full[offset:]

# Downsample to reduce frame jitter (symbol grouping)
symbol_stride = 10  # adjust if needed
symbols_grouped = symbols[::symbol_stride]

# Limit N to avoid extremely large integers
N = min(20, len(symbols_grouped))
symbols_used = symbols_grouped[:N]

# Explicit ternary decode using formula
ternary_value = 0
for i in range(N):
    ternary_value += symbols_used[i] * (3 ** i)

# Stability test across nearby offsets
offsets_to_test = [724, 725, 726]
stability_results = {}

for off in offsets_to_test:
    if off < len(symbols_full):
        temp = symbols_full[off::symbol_stride][:N]
        value = sum(temp[i] * (3 ** i) for i in range(len(temp)))
        stability_results[off] = value

with open("output_8_ternary_decode.txt", "w") as f_out:
    f_out.write(f"Offset used: 725 frames\n")
    f_out.write(f"Symbol stride: {symbol_stride}\n")
    f_out.write(f"Symbols used: {symbols_used.tolist()}\n")
    f_out.write(f"N used: {N}\n")
    f_out.write(f"Decoded ternary value: {ternary_value}\n")
    f_out.write("\nStability test (neighbor offsets):\n")
    for k, v in stability_results.items():
        f_out.write(f"Offset {k}: {v}\n")

print("Updated ternary decode complete.")