95 lines
3.0 KiB
Python
95 lines
3.0 KiB
Python
import struct, os
|
|
|
|
d = r"C:\Project\MAG160C\analysis\pairs_recheck_20260811"
|
|
npix = 19200
|
|
nsegs = 3
|
|
|
|
def load16(name, n):
|
|
b = open(os.path.join(d, name), "rb").read()
|
|
return struct.unpack_from("<%dH" % n, b, 0)
|
|
|
|
def loadi16(name, n):
|
|
b = open(os.path.join(d, name), "rb").read()
|
|
return struct.unpack_from("<%dh" % n, b, 0)
|
|
|
|
f20 = load16("pair_000.f20", npix)
|
|
ref = load16("pair_000.ref", npix)
|
|
raw = load16("pair_000.raw", npix)
|
|
|
|
# working tables (from build-artifacts, verified = interp of EP2/EP3)
|
|
gain = load16(r"C:\Project\MAG160C\build-artifacts\mag160c_official_nuc_gain.bin", npix*nsegs*2)
|
|
thr = loadi16(r"C:\Project\MAG160C\build-artifacts\mag160c_official_nuc_thr.bin", npix*(nsegs-1))
|
|
|
|
# blind records from DDT (all endpoints identical; use EP2)
|
|
ddt = open(r"C:\Project\MAG160C\analysis\ida\Core160043865.ddt", "rb").read()
|
|
thr_ep0 = 0x80
|
|
block = (nsegs-1)*npix*2 + 1024 + nsegs*npix*4
|
|
blind_off = thr_ep0 + 6*block
|
|
recs = []
|
|
for i in range(32):
|
|
r = blind_off + i*40
|
|
target, typ = struct.unpack_from("<II", ddt, r)
|
|
neigh = struct.unpack_from("<8I", ddt, r+8)
|
|
recs.append((target, typ, neigh))
|
|
print("blind recs:", len(recs), "first:", recs[0])
|
|
|
|
def blind_apply(buf, recs):
|
|
out = list(buf)
|
|
for (t, typ, ng) in recs:
|
|
if typ == 3:
|
|
v = out[ng[0]] + out[ng[1]] + out[ng[2]]
|
|
out[t] = v // 3
|
|
elif typ == 4:
|
|
v = out[ng[0]] + out[ng[1]] + out[ng[2]] + out[ng[3]]
|
|
out[t] = v >> 2
|
|
elif typ == 5:
|
|
v = out[ng[0]] + out[ng[1]] + out[ng[2]] + out[ng[3]] + out[ng[4]]
|
|
out[t] = v // 5
|
|
elif typ == 6:
|
|
v = sum(out[ng[j]] for j in range(6))
|
|
out[t] = v // 6
|
|
elif typ == 7:
|
|
v = sum(out[ng[j]] for j in range(7))
|
|
out[t] = v // 7
|
|
elif typ == 8:
|
|
v = sum(out[ng[j]] for j in range(8))
|
|
out[t] = v >> 3
|
|
return out
|
|
|
|
# NUC with exact official math
|
|
nuc = [0]*npix
|
|
for i in range(npix):
|
|
d2 = (int(f20[i]) - int(ref[i])) >> 1
|
|
seg = 0
|
|
if d2 > thr[i*2]:
|
|
seg = 1
|
|
if d2 > thr[i*2+1]:
|
|
seg = 2
|
|
p = (seg*npix + i)*2
|
|
v = (gain[p] * d2) >> 12 + gain[p+1] # CAREFUL: precedence! official: off + ((gain*d2)>>12)
|
|
v = gain[p+1] + ((gain[p] * d2) >> 12)
|
|
if v < 0: v = 0
|
|
if v > 65535: v = 65535
|
|
nuc[i] = v
|
|
|
|
# stats before blind
|
|
diff = sum(1 for i in range(npix) if nuc[i] != raw[i])
|
|
sad = sum(abs(nuc[i]-raw[i]) for i in range(npix))
|
|
print(f"NUC only: diff={diff}/{npix} MAE={sad/npix:.4f}")
|
|
|
|
# apply blind comp
|
|
out = blind_apply(nuc, recs)
|
|
diff2 = sum(1 for i in range(npix) if out[i] != raw[i])
|
|
sad2 = sum(abs(out[i]-raw[i]) for i in range(npix))
|
|
print(f"NUC+blind: diff={diff2}/{npix} MAE={sad2/npix:.4f}")
|
|
|
|
# list mismatch pixels
|
|
mism = [(i, nuc[i], out[i], raw[i]) for i in range(npix) if out[i] != raw[i]]
|
|
print("mismatches:", len(mism))
|
|
for m in mism[:15]:
|
|
print(" pix", m)
|
|
|
|
# blind target pixels check
|
|
for (t, typ, ng) in recs:
|
|
print(f" rec target={t} type={typ} neigh={list(ng[:typ-2])}")
|