#!/usr/bin/python3
import sys
import numpy as np
from scipy.io import wavfile
from scipy.signal import hilbert

def decode_msf_frame(frame_bits):
    """
    Decodiert einen 60-Sekunden-Rahmen von Bit A Daten gemäß NPL-MSF Spezifikation.
    """
    try:
        # BCD-Gewichtungen für MSF (NPL Spec):
        # Jahr: 17-20 (80,40,20,10), 21-24 (8,4,2,1)
        year_tens = sum(b * w for b, w in zip(frame_bits[17:21], [8, 4, 2, 1]))
        year_ones = sum(b * w for b, w in zip(frame_bits[21:25], [8, 4, 2, 1]))
        year = year_tens * 10 + year_ones

        # Monat: 25 (10), 26-29 (8,4,2,1)
        month_tens = frame_bits[25] * 1
        month_ones = sum(b * w for b, w in zip(frame_bits[26:30], [8, 4, 2, 1]))
        month = month_tens * 10 + month_ones

        # Tag: 30-31 (20,10), 32-35 (8,4,2,1)
        day_tens = sum(b * w for b, w in zip(frame_bits[30:32], [2, 1]))
        day_ones = sum(b * w for b, w in zip(frame_bits[32:36], [8, 4, 2, 1]))
        day = day_tens * 10 + day_ones

        # Wochentag: 36-38 (4,2,1) -> 0=So, 1=Mo ... 6=Sa
        weekday_num = sum(b * w for b, w in zip(frame_bits[36:39], [4, 2, 1]))

        # Stunde: 39-40 (20,10), 41-44 (8,4,2,1)
        hour_tens = sum(b * w for b, w in zip(frame_bits[39:41], [2, 1]))
        hour_ones = sum(b * w for b, w in zip(frame_bits[41:45], [8, 4, 2, 1]))
        hour = hour_tens * 10 + hour_ones

        # Minute: 45-47 (40,20,10), 48-51 (8,4,2,1)
        min_tens = sum(b * w for b, w in zip(frame_bits[45:48], [4, 2, 1]))
        min_ones = sum(b * w for b, w in zip(frame_bits[48:52], [8, 4, 2, 1]))
        minute = min_tens * 10 + min_ones

        # Sommerzeit (BST): Sekunde 53 (1 = aktiv)
        is_bst = frame_bits[53] == 1 if len(frame_bits) > 53 else False

        days = ["Sonntag", "Montag", "Dienstag", "Mittwoch", "Donnerstag", "Freitag", "Samstag"]
        day_str = days[weekday_num] if weekday_num < len(days) else f"Tag {weekday_num}"
        tz_str = "BST (UTC+1)" if is_bst else "UTC"

        return f"{day_str}, {day:02d}.{month:02d}.20{year:02d} - {hour:02d}:{minute:02d}:00 {tz_str}"
    except Exception as e:
        return f"Fehler beim Parsen: {e}"

def process_wav(wav_pfad):
    print(f"Lese Audio-Datei: {wav_pfad}...")
    samplerate, data = wavfile.read(wav_pfad)

    if len(data.shape) > 1:
        data = data[:, 0]

    # Normalisieren
    data = data.astype(np.float32)
    data /= np.max(np.abs(data)) if np.max(np.abs(data)) != 0 else 1.0

    # Hüllkurve
    envelope = np.abs(hilbert(data))
    window_len = int(samplerate * 0.01) # 10 ms Glättung
    envelope = np.convolve(envelope, np.ones(window_len)/window_len, mode='same')

    # Schwellenwert für Absenkungen
    base_level = np.median(envelope)
    threshold = base_level * 0.6
    binary = (envelope < threshold).astype(np.int8)

    # Flanken
    diff = np.diff(binary)
    rising = np.where(diff == 1)[0]
    falling = np.where(diff == -1)[0]

    pulses = []
    for r in rising:
        f_candidates = falling[falling > r]
        if len(f_candidates) > 0:
            f = f_candidates[0]
            duration_ms = (f - r) / samplerate * 1000.0
            start_sec = r / samplerate
            if duration_ms >= 40: # Glitches unter 40ms filtern
                pulses.append((start_sec, duration_ms))

    # Minuten-Marker suchen (~500ms Absenkung)
    minute_markers = []
    last_marker_time = -100

    for start_sec, dur in pulses:
        if 400 <= dur <= 600:
            if (start_sec - last_marker_time) > 45.0:
                minute_markers.append(float(start_sec))
                last_marker_time = start_sec

    print(f"\n{len(pulses)} Absenkungs-Impulse gefunden.")
    print(f"{len(minute_markers)} valide Minutenmarker (500 ms Absenkung) bei: {[round(m, 1) for m in minute_markers]} s\n")

    if not minute_markers:
        print("Kein gültiger Minutenmarker gefunden.")
        return

    print("=" * 60)
    print(" DEKODIERTE MSF ZEIT- UND DATUMSANZEIGE")
    print("=" * 60)

    for m_idx, marker_start in enumerate(minute_markers):
        frame_bits = [0] * 60
        next_marker = minute_markers[m_idx+1] if m_idx+1 < len(minute_markers) else marker_start + 60.0

        for start_sec, dur in pulses:
            if marker_start <= start_sec < next_marker:
                # Exaktes Sekunden-Raster berechnen
                rel_time = start_sec - marker_start
                rel_sec = int(round(rel_time))
                sub_sec_offset = rel_time - rel_sec

                # Nur Bit A auswerten (Impuls direkt am Sekundenanfang, Offset < 200ms)
                if 0 <= rel_sec < 60 and abs(sub_sec_offset) < 0.2:
                    if 140 <= dur <= 280:   # 200 ms Absenkung -> Bit A = 1
                        frame_bits[rel_sec] = 1
                    elif 60 <= dur < 140:   # 100 ms Absenkung -> Bit A = 0
                        frame_bits[rel_sec] = 0

        zeit_text = decode_msf_frame(frame_bits)
        print(f"[@ {marker_start:6.1f}s Audio-Offset] -> {zeit_text}")

    print("=" * 60)

if __name__ == "__main__":
    datei = sys.argv[1] if len(sys.argv) > 1 else "websdr_recording_60khz.wav"
    process_wav(datei)