#!/usr/bin/env python3
import sys, os, struct
import numpy as np

SR = 44100
TONES = [110, 220, 440, 880, 1760, 3520, 7040, 14080]

def read_wav(path):
    d = open(path, 'rb').read()
    i = 12
    ch = bps = None
    data_off = data_sz = None
    while i + 8 <= len(d):
        cid = d[i:i+4]; sz = struct.unpack('<I', d[i+4:i+8])[0]
        if cid == b'fmt ':
            ch = struct.unpack('<H', d[i+8:i+10])[0]
            bps = struct.unpack('<H', d[i+16:i+18])[0]
        if cid == b'data':
            data_off = i + 8; data_sz = sz
        i += 8 + sz + (sz & 1)
    raw = d[data_off:data_off+data_sz]
    if bps == 32:
        a = np.frombuffer(raw, dtype='<f4')
        return a.reshape(-1, ch)[:, 0], ch
    if bps == 24:
        b = np.frombuffer(raw, dtype=np.uint8).reshape(-1, ch*3)
        v = b[:, 0] | (b[:,1].astype(np.int32) << 8) | (b[:,2].astype(np.int32) << 16)
        v = (v ^ (1 << 23)) - (1 << 23)
        return v / 8388608.0, ch
    a = np.frombuffer(raw, dtype='<i2')
    return a.reshape(-1, ch)[:, 0] / 32768.0, ch

def tone_db(x, f, sr=SR):
    n = len(x)
    w = np.hanning(n*2-1) if False else None
    # use Goertzel over entire file
    w0 = 2*np.pi*f/sr
    c = 2*np.cos(w0)
    s0=s1=s2=0.0
    for v in x:
        s0 = v + c*s1 - s2
        s2 = s1; s1 = s0
    p = s1*s1 + s2*s2 - c*s1*s2
    return 10*np.log10(p + 1e-30) - 20*np.log10(n)  # normalize by length

def main():
    outs = sys.argv[1:]
    colw = 12
    n0 = None
    datas = {}
    for p in outs:
        y, ch = read_wav(p)
        datas[p] = y
        if n0 is None or len(y) > n0: n0 = len(y)
    print(f"{'Hz':>5}" + ''.join(f" {os.path.basename(p)[:colw]:>{colw}}" for p in outs))
    for f in TONES:
        row = f"{f:5d}"
        base = None
        for p in outs:
            d = tone_db(datas[p], f)
            if base is None: base = d
            row += f" {d:>{colw}.1f}"
        print(row, f"  (rel ref {outs[0]})")

if __name__ == '__main__':
    main()