import os
import sys
import time
import subprocess
import ctypes
import numpy as np
import cv2

# We'll try to use sounddevice for audio as it's more robust for pip installs
try:
    import sounddevice as sd
except ImportError:
    print("Warning: sounddevice not found. Audio features will be disabled.")
    print("Please run: pip install sounddevice numpy opencv-python")
    sd = None

# ---------------------------------------------------------------------
# Configuration
# ---------------------------------------------------------------------
GRID_W = 36
GRID_H = 28
TOTAL_PIXELS = GRID_W * GRID_H
SAMPLE_RATE = 48000
FFT_SIZE = 2048
CAVITY_LO_BIN = 5
CAVITY_HI_BIN = 25
EMA_ALPHA = 0.15
CONTROL_UPDATE_HZ = 30.0

# ---------------------------------------------------------------------
# Inline C Code for performance-critical grid processing
# ---------------------------------------------------------------------
C_CODE = r"""
#include <math.h>
#include <string.h>
#include <stdio.h>

typedef struct {
    float winding;
    float hyperfine;
    float entropy;
    float accuracy;
    int isEPR;
    int isViolation;
    float avgMotion;
} Metrics;

void process_grid(const float* gray, float* prev_luma, int w, int h, 
                  float cavityEnergy, Metrics* out, char* statusMsg) {
    float totalDX = 0.0f, totalDY = 0.0f;
    float totalMotion = 0.0f;
    float lumaSum = 0.0f;
    int totalPixels = w * h;

    float spinTheta[1008]; // Max size for 36*28
    float spinMag[1008];

    // Entropy and Spin Analysis
    float totalMag = 0.0f;
    for (int i = 0; i < totalPixels; i++) {
        spinTheta[i] = 0.0f;
        spinMag[i] = 0.0f;
    }

    for (int y = 0; y < h; ++y) {
        for (int x = 0; x < w; ++x) {
            int idx = y * w + x;
            float lum = gray[idx];
            lumaSum += lum;

            if (y > 0 && y < h - 1 && x > 0 && x < w - 1) {
                float left  = gray[y * w + (x - 1)];
                float right = gray[y * w + (x + 1)];
                float up    = gray[(y - 1) * w + x];
                float down  = gray[(y + 1) * w + x];

                float dx = right - left;
                float dy = down - up;
                float theta = atan2f(dy, dx);
                float mag = sqrtf(dx*dx + dy*dy);

                float motion = fabsf(lum - prev_luma[idx]);
                totalMotion += motion;

                totalDX += cosf(theta) * mag;
                totalDY += sinf(theta) * mag;

                spinTheta[idx] = theta;
                spinMag[idx] = mag;
            }
            prev_luma[idx] = lum;
        }
    }

    float avgMotion = totalMotion / totalPixels;
    float meanLuma = lumaSum / totalPixels;

    // Winding number
    float circulation = 0.0f;
    for (int y = 2; y < h - 2; ++y) {
        for (int x = 2; x < w - 2; ++x) {
            int idx = y * w + x;
            float t = spinTheta[idx];
            circulation += (sinf(t) * totalDX - cosf(t) * totalDY);
        }
    }
    out->winding = fabsf(circulation) / (totalPixels * 2.0f);
    if (out->winding > 1.0f) out->winding = 1.0f;

    // Hyperfine (Scaled to match C++ logic)
    out->hyperfine = avgMotion * cavityEnergy * 8.0f;
    out->avgMotion = avgMotion;

    // Entropy
    for (int i = 0; i < totalPixels; i++) totalMag += spinMag[i] + 0.001f;
    float entropyVal = 0.0f;
    for (int i = 0; i < totalPixels; i++) {
        float p = (spinMag[i] + 0.001f) / totalMag;
        if (p > 0.0f) entropyVal -= p * log2f(p);
    }
    float maxEnt = log2f((float)totalPixels);
    out->entropy = entropyVal;
    float normEntropy = entropyVal / maxEnt;
    if (normEntropy > 1.0f) normEntropy = 1.0f;

    // Accuracy
    out->accuracy = (1.0f - normEntropy) * (out->hyperfine + 0.01f) * 4.0f;
    if (out->accuracy > 1.0f) out->accuracy = 1.0f;
    if (out->accuracy < 0.0f) out->accuracy = 0.0f;

    // EPR detection (Simplified)
    float leftMotion = 0.0f, rightMotion = 0.0f;
    for (int y = 0; y < h; y++) {
        for (int x = 0; x < w; x++) {
            float delta = fabsf(prev_luma[y*w+x] - meanLuma);
            if (x < w/2) leftMotion += delta;
            else rightMotion += delta;
        }
    }
    float leftCorr = leftMotion / (totalPixels/2);
    float rightCorr = rightMotion / (totalPixels/2);
    
    out->isEPR = (fabsf(leftCorr - rightCorr) < 0.05f && cavityEnergy > 0.15f && avgMotion > 0.02f && out->hyperfine > 0.25f);
    
    out->isViolation = 0;
    if (avgMotion > 0.08f && cavityEnergy < 0.03f) out->isViolation = 1;
    else if (avgMotion < 0.02f && cavityEnergy > 0.25f) out->isViolation = 1;

    if (out->isEPR) strcpy(statusMsg, "EPR BRIDGE DETECTED");
    else if (out->isViolation) strcpy(statusMsg, "SKISS-VIOLATION");
    else if (normEntropy > 0.85f) strcpy(statusMsg, "SKISS-INCOMPLETE");
    else strcpy(statusMsg, "SKISS-COMPLETE (Bonded)");
}
"""

class Metrics(ctypes.Structure):
    _fields_ = [
        ("winding", ctypes.c_float),
        ("hyperfine", ctypes.c_float),
        ("entropy", ctypes.c_float),
        ("accuracy", ctypes.c_float),
        ("isEPR", ctypes.c_int),
        ("isViolation", ctypes.c_int),
        ("avgMotion", ctypes.c_float),
    ]

# ---------------------------------------------------------------------
# Compilation logic
# ---------------------------------------------------------------------
def get_c_lib():
    c_file = "mapsh_core.c"
    so_file = "./mapsh_core.so"
    
    with open(c_file, "w") as f:
        f.write(C_CODE)
    
    try:
        subprocess.check_call(["gcc", "-O3", "-shared", "-o", so_file, "-fPIC", c_file, "-lm"])
    except Exception as e:
        print(f"Error compiling C core: {e}")
        return None
    
    lib = ctypes.CDLL(so_file)
    lib.process_grid.argtypes = [
        ctypes.POINTER(ctypes.c_float),
        ctypes.POINTER(ctypes.c_float),
        ctypes.c_int,
        ctypes.c_int,
        ctypes.c_float,
        ctypes.POINTER(Metrics),
        ctypes.c_char_p
    ]
    return lib

# ---------------------------------------------------------------------
# Main Application Class
# ---------------------------------------------------------------------
class MAPHS:
    def __init__(self):
        self.lib = get_c_lib()
        if not self.lib:
            print("Falling back to (slower) Python implementation is not yet implemented.")
            sys.exit(1)
            
        self.prev_luma = np.zeros(TOTAL_PIXELS, dtype=np.float32)
        self.cavity_energy = 0.0
        self.gain_ema = 0.0
        self.control_output = 0.5
        self.last_gain = -float('inf')
        self.accuracy_ema = 0.0
        
        # Audio buffer
        self.audio_buffer = np.zeros(FFT_SIZE, dtype=np.float32)

    def audio_callback(self, indata, frames, time, status):
        if status:
            print(status)
        # Shift and append
        self.audio_buffer = np.roll(self.audio_buffer, -frames)
        self.audio_buffer[-frames:] = indata[:, 0]
        
        # FFT energy normalized by FFT_SIZE to match C++ logic
        fft_data = np.abs(np.fft.rfft(self.audio_buffer)) / FFT_SIZE
        bins = fft_data[CAVITY_LO_BIN:CAVITY_HI_BIN+1]
        self.cavity_energy = np.mean(bins)

    def compute_control(self, gain_db):
        if self.last_gain > -float('inf'):
            if gain_db > self.last_gain:
                self.control_output += 0.01
            else:
                self.control_output -= 0.01
        self.control_output = max(0.0, min(1.0, self.control_output))
        self.last_gain = gain_db
        return self.control_output

    def run(self):
        print("MAPHS Python - Black Hole Energy Extraction")
        cap = cv2.VideoCapture(0)
        if not cap.isOpened():
            print("Error: Could not open camera.")
            return

        if sd:
            try:
                stream = sd.InputStream(callback=self.audio_callback, channels=1, 
                                        samplerate=SAMPLE_RATE, blocksize=512)
                stream.start()
            except Exception as e:
                print(f"Audio init failed: {e}. Cavity energy will be simulated.")
                self.cavity_energy = 0.05
                sd_active = False
            else:
                sd_active = True
        else:
            print("Running without audio input. Simulating baseline cavity energy.")
            self.cavity_energy = 0.05
            sd_active = False

        try:
            while True:
                ret, frame = cap.read()
                if not ret:
                    break
                
                # If no mic, we still need some energy for hyperfine to be non-zero
                if not sd_active:
                    self.cavity_energy = 0.05 + 0.02 * np.random.random()

                # Video processing
                small = cv2.resize(frame, (GRID_W, GRID_H), interpolation=cv2.INTER_AREA)
                gray = cv2.cvtColor(small, cv2.COLOR_BGR2GRAY).astype(np.float32) / 255.0
                gray_flat = gray.flatten()
                
                metrics = Metrics()
                status_msg = ctypes.create_string_buffer(256)
                
                self.lib.process_grid(
                    gray_flat.ctypes.data_as(ctypes.POINTER(ctypes.c_float)),
                    self.prev_luma.ctypes.data_as(ctypes.POINTER(ctypes.c_float)),
                    GRID_W, GRID_H,
                    float(self.cavity_energy),
                    ctypes.byref(metrics),
                    status_msg
                )

                # EMA and Gain
                super_gain = metrics.hyperfine * (1.0 + metrics.winding)
                self.gain_ema = (1.0 - EMA_ALPHA) * self.gain_ema + EMA_ALPHA * super_gain
                gain_db = 10.0 * np.log10(self.gain_ema) if self.gain_ema > 1e-6 else -100.0
                
                # Accuracy EMA
                self.accuracy_ema = (1.0 - 0.1) * self.accuracy_ema + 0.1 * metrics.accuracy
                
                # Estimated Energy (Power) in arbitrary "Extraction Units"
                # Base it on the gain and the cavity energy
                energy_ext = self.gain_ema * 100.0 
                
                control = self.compute_control(gain_db)
                
                # Output
                msg = status_msg.value.decode()
                print(f"\r[CTRL] mirror={control:.3f} | Gain={gain_db:5.1f}dB | Acc={self.accuracy_ema*100:4.1f}% | Energy={energy_ext:6.2f} mE | {msg:<22}", end="")
                
                # Visualization (Optional)
                cv2.imshow("MAPHS Monitor", cv2.resize(frame, (640, 480)))
                if cv2.waitKey(1) & 0xFF == ord('q'):
                    break

                time.sleep(1.0 / CONTROL_UPDATE_HZ)

        except KeyboardInterrupt:
            print("\nInterrupted by user.")
        finally:
            cap.release()
            cv2.destroyAllWindows()
            if sd:
                stream.stop()
            print("\nShutdown.")

if __name__ == "__main__":
    app = MAPHS()
    app.run()
