From 4c24f28a8818a9a8d386319b946dbb21f9f77bca Mon Sep 17 00:00:00 2001 From: Aron Novak Date: Thu, 11 Jun 2026 10:34:00 +0200 Subject: [PATCH 1/5] Add universal viewer with TC002C Duo radiometric support New src/thermalcam.py auto-detects InfiRay/Topdon USB thermal cameras and decodes real temperatures: TC001 (256x384, /64) and TC002C Duo (512x484, /16, which also yields a clean 512x384 image), with a labelled relative fallback otherwise. Adds docs/TC002C-DUO.md (frame layout + reverse-engineered vendor protocol), requirements.txt, and an experimental command module. Co-Authored-By: Claude Opus 4.8 (1M context) --- README.md | 56 ++++ docs/TC002C-DUO.md | 169 ++++++++++ requirements.txt | 7 + src/thermal_protocol.py | 146 ++++++++ src/thermalcam.py | 712 ++++++++++++++++++++++++++++++++++++++++ 5 files changed, 1090 insertions(+) create mode 100644 docs/TC002C-DUO.md create mode 100644 requirements.txt create mode 100644 src/thermal_protocol.py create mode 100644 src/thermalcam.py diff --git a/README.md b/README.md index 409cad6..cb25716 100644 --- a/README.md +++ b/README.md @@ -6,6 +6,62 @@ https://www.eevblog.com/forum/thermal-imaging/infiray-and-their-p2-pro-discussio Check out Leo's Github here: https://github.com/LeoDJ/P2Pro-Viewer/tree/main +## Universal viewer — `src/thermalcam.py` + +The original program (`src/tc001v4.2.py`, below) is hard-coded for the **Topdon +TC001** frame format. `src/thermalcam.py` is a camera-agnostic rewrite that also +works with newer Topdon / InfiRay models such as the **TC002C Duo**, which expose +a different USB frame geometry. + +What it adds: + +- **Auto-detects** the thermal camera on the USB bus — no `--device` needed. +- **Auto-detects the frame geometry** instead of assuming `256x384`, so it + handles the TC001 (`256x384`) and the TC002C Duo family (`256x392`, …). +- **Real error handling** (the original has none). +- **Real temperatures.** It auto-detects the camera's radiometric mode and + decodes true °C — the TC001 at `256x384` and the **TC002C Duo at `512x484`** + (which also gives a crisp 512×384 image). Cameras with no 16-bit mode fall back + to a clearly-labelled *relative* scale. The absolute reading can be offset- + calibrated with `--temp-offset` / the `[` `]` keys. See `docs/TC002C-DUO.md`. +- Same niceties: colormaps, HUD, recording, snapshots, scaling, blur, contrast, + hot/cold spot tracking — plus fixed-pattern-noise removal and a contrast stretch + for cameras that don't pre-AGC their preview. + +### Install + +```bash +sudo apt-get install python3-opencv # Debian/Ubuntu/Raspberry Pi +# or: pip install -r requirements.txt +``` + +### Run + +```bash +python3 src/thermalcam.py # auto-detect everything +python3 src/thermalcam.py --device /dev/video2 +python3 src/thermalcam.py --resolution 256x384 # e.g. force TC001 radiometric mode +python3 src/thermalcam.py --selftest 20 # headless: save a snapshot + diagnostics, no GUI +``` + +Useful options: `--temp-offset N` to correct the absolute °C, `--temp-scale +{auto,64,16}` / `--temp-order {le,be}` to force a radiometric decode, +`--no-stretch` / `--no-destripe` / `--smooth` for the image cleanups. + +### Keys + +`a/z` blur · `s/x` min-max threshold · `d/c` scale · `f/v` contrast · `m` colormap +· `h` HUD · `n` swap bands · `g` temporal smoothing · `e/w` fullscreen on/off · +`r/t` record/stop · `p` snapshot · `[`/`]` temperature-offset calibration · `q`/ESC quit + +### Camera support + +| Camera | Image | Temperature | +|--------|-------|-------------| +| Topdon TC001 | ✅ | ✅ real °C (`256x384`, /64) | +| Topdon TC002C Duo | ✅ hi-res 512×384 | ✅ real °C (`512x484`, /16; offset-calibratable) | +| Other InfiRay UVC | ✅ likely | auto-detected if a radiometric mode is present | + ## Introduction diff --git a/docs/TC002C-DUO.md b/docs/TC002C-DUO.md new file mode 100644 index 0000000..a9825a0 --- /dev/null +++ b/docs/TC002C-DUO.md @@ -0,0 +1,169 @@ +# Topdon TC002C Duo (and other InfiRay UVC cameras) on Linux + +Notes from getting `src/thermalcam.py` working with a **Topdon TC002C Duo**, and a +reverse-engineering log of the radiometric (true-temperature) protocol pulled +out of the Android app. Written so the work can be continued. + +The TC002C Duo is a dual-sensor (thermal + visible) camera built on the InfiRay +`iruvc` engine — the same family as the TC001 this repo originally targeted, but +it behaves differently over USB. + +USB id: `2bdf:0102` (`Topdon`). It enumerates as a standard UVC camera +(`/dev/videoN`), so `thermalcam.py` auto-detects it. + +## Frame geometry + +`v4l2-ctl -d /dev/video2 --list-formats-ext` reports a partly-garbage table (the +firmware exposes bogus descriptors like `4x12305`, `60x3299`). The real, usable +YUYV modes are: + +| Resolution | What it contains | +|-----------|------------------| +| `256x192` | thermal preview only (8-bit) | +| `256x392` | stacked: two 8-bit thermal planes (**no** 16-bit data) | +| `512x384` | 2× thermal preview, 8-bit only | +| `512x484` | **16-bit temperature band (top) + 512x384 thermal image** ← use this | +| `644x384`, `520x192` | other banded/fused layouts | + +**The `512x484` mode is the one to use.** It is a stack of horizontal bands: + +``` +rows 0 - 95 : 16-bit temperature data (256x192, little-endian, /16 Kelvin), + packed two sensor rows per 512-wide line -> 96*512*2 = 256*192*2 bytes +rows 96 - 97 : telemetry +rows 98 -481 : 512x384 thermal AGC image (8-bit, what we display) +rows 482 -483 : telemetry +``` + +Read naively as one YUYV image the temperature band looks like a green/torn +stripe (its high byte is mis-read as chroma) — which is exactly why it was first +mistaken for garbage. Decoded as 16-bit it is a clean radiometric frame. + +Two gotchas cost real debugging time and are handled in `thermalcam.py`: + +* **OpenCV buffering** hands back stale, mid-frame-desynced buffers. Set + `CAP_PROP_BUFFERSIZE = 1`. (`v4l2-ctl --stream-mmap` never desyncs.) +* **Opening an unsupported resolution desyncs the next stream.** Enumerate the + advertised modes (`VIDIOC_ENUM_FRAMESIZES`) and never probe one the camera + doesn't list. The occasional desynced frame is also detected (its shifted + bytes decode to absurd temperatures) and skipped. + +## Temperature + +Decode the band as little-endian 16-bit and apply the standard InfiRay TPD +formula: + +``` +temp_C = raw / 16 - 273.15 +``` + +This gives a stable, spatially-correct field (the warm/cool structure matches the +scene exactly). The **absolute** value can be offset by several °C because we are +not running the camera's shutter-based NUC or ambient/emissivity correction — use +`--temp-offset` / the `[` `]` keys to dial it in against a known reference (e.g. +your own skin ~34 °C). The original TC001 uses the same scheme at `256x384` but +with a `/64` divisor; `thermalcam.py` knows both via `KNOWN_LAYOUTS`. + +> Earlier notes in this file claimed the Duo exposes *no* 16-bit data over UVC. +> That was wrong — it was based on only testing `256x392`. The `512x484` mode +> carries full radiometric data with no vendor command required. + +## Vendor command protocol (still useful for image quality) + +True temperatures need no command, but the **image** still shows the sensor's +uncorrected non-uniformity (blotchy on flat scenes) because we can't trigger the +camera's shutter NUC over plain UVC. That — and finer temperature calibration — +is what the vendor command protocol below is for. It was reverse-engineered from +the Android app's native libraries before the `512x484` mode was found. + +## The vendor command protocol (reverse-engineered) + +The Android APK's Kotlin/Java only holds the InfiRay SDK *wrapper* +(`com.energy.iruvc.*`); the real protocol is in the native `.so` libraries, +which live in the **`split_config.arm64_v8a.apk`** split (pull from the phone via +`adb`, see below). The relevant libs: + +- `libircmd.so` — the `ir_cmd` command set + temperature math +- `libUSBUVCCamera.so` — a libuvc + libusb fork; `UVCCamera::sendCommand` is the transport +- `libUSBDualCamera.so` — the Duo dual-sensor path +- `libadvirtemp.so`, `libdualcalibration.so` — Y16→°C + Duo calibration + +### Transport: `UVCCamera::sendCommand` + +Commands are **USB vendor control transfers** (not UVC-class), i.e. plain +`libusb_control_transfer`: + +``` +sendCommand(this, bmRequestType, bRequest, wValue, wIndex, data, wLength) + -> libusb_control_transfer(handle, bmRequestType, bRequest, + wValue, wIndex, data, wLength, timeout=1000) +``` + +| Operation | bmRequestType | bRequest | wValue | +|-----------|--------------|----------|--------| +| register WRITE | `0x41` (OUT, vendor, interface) | `0x45` (`'E'`) | `0x0078` | +| register READ | `0xc1` (IN, vendor, interface) | `0x44` (`'D'`) | `0x0078` | + +`wIndex` is a **device register**: + +| Register | Meaning | +|----------|---------| +| `0x9d00` | command descriptor (opcode + arg length) | +| `0x9d08` | command payload (bulk data) | +| `0x1d00`, `0x1d08`, `0x1d40` | command parameters | +| `0x1d04`, `0x1d08` | result read-back | +| `0x200` | status register (poll: bit0 = busy, retry; error if >3) | + +### Worked example: `preview_start` + +Decoded from `Java_com_energy_iruvc_ircmd_LibIRCMD_preview_1start` in `libircmd.so`: + +``` +1. WRITE 0x9d00 <- 0f c1 00 00 00 00 00 08 # opcode 0xc10f, 8-byte arg block +2. WRITE 0x1d08 <- [fps, + (source==1 ? 0x80 : 0x00), # source: 0=SENSOR, 1=FIX_PATTERN + width_hi, width_lo, # big-endian + height_hi, height_lo, # big-endian + mode, # StartPreviewMode (0=VOC_DVP_MODE) + path] # PreviewPathChannel +3. POLL 0x200 (read 1 byte) until (status & 1) == 0; abort if status > 3 +``` + +`thermal_protocol.py` in this directory implements this transport and command. +**It is experimental and unverified on hardware** — see the safety notes there. + +## Optional next steps (nice-to-have, not required) + +Basic temperatures and a clean hi-res image already work from the `512x484` +stream with no commands. The vendor protocol would let us go further: + +* **Shutter NUC** — periodically trigger the sensor's flat-field correction + (`shutter_manual_switch` / `shutter_sta_set`) to remove the slowly-drifting + blotchy non-uniformity that's visible on flat scenes. +* **Calibrated °C** — apply the camera's emissivity / ambient / NUC correction + (`libadvirtemp.so` + calibration assets `nuc_table_high.bin`, `tau*.bin`, + `kt_high.bin`, `bt_high.bin`, `dual_calibration_parameters2.bin`) instead of the + raw linear formula, removing the absolute offset. +* **Visible / fused image** — the `644x384` and other modes appear to carry the + Duo's visible camera for thermal+visible fusion (`libUSBDualCamera.so`). + +To pursue these, either continue the static RE the way `preview_start` was +decoded, or capture the real app's USB traffic (phone with `usbmon`+root, or +Waydroid on Linux with the camera passed through). None of it is needed for the +working radiometric viewer. + +Once the stream delivers genuine 16-bit data, `thermalcam.py` already +auto-detects it (`auto_detect_radiometric`) and shows real °C — or force it with +`--temp-scale`. + +## Pulling the native libs from the phone + +```bash +adb shell pm path com.topdon.topInfrared # find base + split APKs +adb pull <.../split_config.arm64_v8a.apk> /tmp/split.apk +unzip -o /tmp/split.apk 'lib/*' -d /tmp/topdon_libs +ls /tmp/topdon_libs/lib/arm64-v8a/ # *.so +``` + +The libs are ARM64 Android (won't run on an x86_64 laptop) but disassemble fine +(Ghidra/radare2), which is how the protocol above was recovered. diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..ffe746e --- /dev/null +++ b/requirements.txt @@ -0,0 +1,7 @@ +# Core viewer (src/thermalcam.py) +opencv-python>=4.5 +numpy>=1.20 + +# Only for the experimental radiometric-unlock work (src/thermal_protocol.py). +# Not needed for the viewer. +# pyusb>=1.2 diff --git a/src/thermal_protocol.py b/src/thermal_protocol.py new file mode 100644 index 0000000..f9f4f53 --- /dev/null +++ b/src/thermal_protocol.py @@ -0,0 +1,146 @@ +#!/usr/bin/env python3 +""" +thermal_protocol.py - EXPERIMENTAL InfiRay/Topdon vendor command transport. + +Reverse-engineered from the Topdon Android app's native libraries +(libircmd.so / libUSBUVCCamera.so). See docs/TC002C-DUO.md for the full writeup. + +NOTE: real temperatures already work in thermalcam.py with NO commands at all +(via the 512x484 radiometric mode). This module is NOT needed for that. It exists +for the optional next steps - shutter NUC, full °C calibration, the visible/fused +image - which do need the vendor command channel. + +>>> THIS IS UNVERIFIED ON HARDWARE AND NOT WIRED INTO thermalcam.py. <<< +By default it runs DRY (prints the control transfers it WOULD send). Pass --run +to actually talk to the camera, at your own risk. + +Safety: + * Only the safe register reads and the `preview_start` command are implemented. + * The flash/OEM/erase commands from the SDK are deliberately NOT implemented. + * A USB replug resets the camera if it ends up in a bad state. + +Linux caveat: these are USB *vendor* control transfers (bmRequestType 0x41/0xc1), +not UVC-class controls, so the kernel uvcvideo driver does not expose them. To +send them, libusb must talk to the device directly, which conflicts with +streaming the same interface via v4l2 at the same time. Expect to claim/detach +the interface here and do the video capture through libuvc rather than v4l2 once +this path works. That integration is the open task. +""" + +import argparse +import sys +import time + +VID, PID = 0x2BDF, 0x0102 + +# Transport constants (see docs/TC002C-DUO.md). +REQTYPE_WRITE = 0x41 # OUT | vendor | interface +REQTYPE_READ = 0xC1 # IN | vendor | interface +BREQUEST_WRITE = 0x45 # 'E' +BREQUEST_READ = 0x44 # 'D' +WVALUE = 0x0078 +TIMEOUT_MS = 1000 + +# Registers. +REG_CMD = 0x9D00 # command descriptor +REG_PAYLOAD = 0x9D08 # bulk payload +REG_PARAM = 0x1D08 # command parameters +REG_STATUS = 0x0200 # status poll + + +class InfiRayDevice: + def __init__(self, dry_run=True): + self.dry_run = dry_run + self.dev = None + self.usb = None + + def open(self): + import usb.core # pyusb; only needed for a real run + self.usb = usb.core + dev = usb.core.find(idVendor=VID, idProduct=PID) + if dev is None: + raise RuntimeError(f"No {VID:04x}:{PID:04x} device found") + self.dev = dev + if not self.dry_run: + # The kernel uvcvideo driver owns the interface; vendor control + # transfers need it detached. (This is what stops simultaneous v4l2.) + for cfg in dev: + for intf in cfg: + n = intf.bInterfaceNumber + if dev.is_kernel_driver_active(n): + dev.detach_kernel_driver(n) + return self + + def write_reg(self, reg, data): + """Vendor WRITE of `data` bytes to register `reg`.""" + data = bytes(data) + if self.dry_run: + print(f" WRITE reg=0x{reg:04x} <- {data.hex(' ')}") + return len(data) + return self.dev.ctrl_transfer(REQTYPE_WRITE, BREQUEST_WRITE, WVALUE, reg, data, TIMEOUT_MS) + + def read_reg(self, reg, length): + """Vendor READ of `length` bytes from register `reg`.""" + if self.dry_run: + print(f" READ reg=0x{reg:04x} ({length} bytes) [dry-run -> zeros]") + return bytes(length) + return bytes(self.dev.ctrl_transfer(REQTYPE_READ, BREQUEST_READ, WVALUE, reg, length, TIMEOUT_MS)) + + def poll_status(self, tries=1000): + """Poll REG_STATUS until the command completes (bit0 clear). Returns the + last status byte; values > 3 indicate an error per the SDK.""" + for _ in range(tries): + status = self.read_reg(REG_STATUS, 1)[0] + if self.dry_run: + return 0 + if status & 1: # busy + continue + if (status >> 1) & 1 and status > 3: + raise RuntimeError(f"command error, status=0x{status:02x}") + return status + raise TimeoutError("status poll timed out") + + def preview_start(self, width, height, fps=25, mode=0, source=0, path=0): + """The decoded `preview_start` command (see docs/TC002C-DUO.md). + + NOTE: starts the sensor preview; it does not by itself switch the UVC + stream to 16-bit. The data-flow-mode command is still to be recovered. + """ + # 1. command descriptor: opcode 0xc10f, 8-byte arg block. + self.write_reg(REG_CMD, bytes([0x0F, 0xC1, 0, 0, 0, 0, 0, 0x08])) + # 2. parameter block (width/height big-endian). + params = bytes([ + fps & 0xFF, + 0x80 if source == 1 else 0x00, + (width >> 8) & 0xFF, width & 0xFF, + (height >> 8) & 0xFF, height & 0xFF, + mode & 0xFF, + path & 0xFF, + ]) + self.write_reg(REG_PARAM, params) + # 3. wait for completion. + return self.poll_status() + + +def main(argv=None): + p = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) + p.add_argument("--run", action="store_true", + help="Actually send to the camera (default: dry-run, prints only).") + p.add_argument("--width", type=int, default=256) + p.add_argument("--height", type=int, default=384) + args = p.parse_args(argv) + + dev = InfiRayDevice(dry_run=not args.run) + if args.run: + try: + dev.open() + except Exception as exc: # pyusb missing, no perms, device busy, ... + sys.exit(f"open failed: {exc}\n(Try: pip install pyusb; run as root; unplug other users.)") + print(f"preview_start({args.width}x{args.height}){' [DRY-RUN]' if not args.run else ''}:") + dev.preview_start(args.width, args.height) + print("done.") + return 0 + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/src/thermalcam.py b/src/thermalcam.py new file mode 100644 index 0000000..2dc43a7 --- /dev/null +++ b/src/thermalcam.py @@ -0,0 +1,712 @@ +#!/usr/bin/env python3 +""" +thermalcam.py - Universal Linux viewer for InfiRay-based USB thermal cameras. + +Originally written by Les Wright for the Topdon TC001 (see tc001v4.2.py), this +version is a camera-agnostic rewrite that also works with newer Topdon/InfiRay +models such as the TC002C Duo. + +These cameras pack a frame as a stack of horizontal bands: a high-res thermal +*image* band (8-bit AGC, neutral chroma) plus a 16-bit *temperature* band (raw +radiometric data, which shows up with extreme chroma when mis-read as YUYV). The +band geometry differs per model: + + * TC001 256x384 : image 256x192 + temp 256x192 (16-bit, /64 Kelvin) + * TC002C Duo 512x484 : temp 256x192 (16-bit, /16 Kelvin) + image 512x384 + +This viewer auto-detects the camera, picks a mode that carries real temperature +data, decodes it, and shows real degrees C. If no radiometric mode is available +it falls back to a clearly-labelled relative scale. See docs/TC002C-DUO.md. +""" + +import argparse +import fcntl +import glob +import os +import re +import struct +import sys +import time + +try: + import cv2 + import numpy as np +except ImportError as exc: # pragma: no cover - environment guard + sys.exit( + f"Missing dependency: {exc.name}. Install with:\n" + " sudo apt-get install python3-opencv (Debian/Ubuntu/Raspberry Pi)\n" + " pip install opencv-python numpy (everything else)" + ) + +# USB vendor IDs known to ship InfiRay-based thermal cameras. 0x2BDF is Topdon. +KNOWN_THERMAL_VIDS = {0x2BDF, 0x0BDA, 0x1514, 0x3474} + +KELVIN = 273.15 # temp_C = raw/scale - 273.15 + +# Exact band geometry per capture resolution, keyed by (width, height): +# temp_rows, image_rows, sensor (w,h), scale (raw units per Kelvin), byte order. +# Fixed geometry is far more robust than per-frame chroma detection, which breaks +# on the occasional desynced/metadata-header frame these cameras emit. +KNOWN_LAYOUTS = { + (256, 384): ((192, 384), (0, 192), (256, 192), 64.0, "le"), # TC001 / P2 family + (512, 484): ((0, 96), (98, 482), (256, 192), 16.0, "le"), # TC002C Duo +} +# Resolutions to try for genuine temperature data, in order of preference. +RADIOMETRIC_MODES = [(256, 384), (512, 484)] +# Plain image modes (no temperature), used only if no radiometric mode works. +IMAGE_MODES = [(512, 384), (256, 392), (256, 192)] + +COLORMAPS = [ + (cv2.COLORMAP_JET, "Jet"), + (cv2.COLORMAP_INFERNO, "Inferno"), + (cv2.COLORMAP_HOT, "Hot"), + (cv2.COLORMAP_MAGMA, "Magma"), + (cv2.COLORMAP_PLASMA, "Plasma"), + (cv2.COLORMAP_BONE, "Bone"), + (cv2.COLORMAP_VIRIDIS, "Viridis"), + (cv2.COLORMAP_RAINBOW, "Rainbow"), +] + + +def is_raspberrypi(): + try: + with open("/sys/firmware/devicetree/base/model", "r") as m: + return "raspberry pi" in m.read().lower() + except OSError: + return False + + +def _sysfs_usb_ids(video_node): + """Return (vid, pid, name) for a /dev/videoN node by walking sysfs, or None.""" + name = os.path.basename(video_node) + base = f"/sys/class/video4linux/{name}" + try: + modalias = open(os.path.join(base, "device", "modalias")).read().strip() + match = re.search(r"v([0-9A-Fa-f]{4})p([0-9A-Fa-f]{4})", modalias) + if not match: + return None + vid, pid = int(match.group(1), 16), int(match.group(2), 16) + except OSError: + return None + try: + product = open(os.path.join(base, "name")).read().strip() + except OSError: + product = "" + return vid, pid, product + + +def find_thermal_device(): + """Locate the thermal camera's /dev/videoN, preferring known thermal VIDs.""" + nodes = sorted(glob.glob("/dev/video*"), + key=lambda p: int(p[10:]) if p[10:].isdigit() else 999) + fallback = None + for node in nodes: + ids = _sysfs_usb_ids(node) + if ids is None: + continue + vid, _pid, product = ids + if vid in KNOWN_THERMAL_VIDS and _can_capture(node): + return node + if fallback is None and "USB Camera" in product and _can_capture(node): + fallback = node + return fallback + + +def supported_resolutions(device): + """Set of (w, h) YUYV modes the device advertises, via VIDIOC_ENUM_FRAMESIZES. + + Used to avoid even *opening* an unsupported resolution: doing so desyncs these + cameras' subsequent streams. Returns None if enumeration fails. + """ + yuyv = ord("Y") | (ord("U") << 8) | (ord("Y") << 16) | (ord("V") << 24) + # _IOWR('V', 74, sizeof(struct v4l2_frmsizeenum)==44) + vidioc_enum_framesizes = (3 << 30) | (44 << 16) | (ord("V") << 8) | 74 + try: + fd = os.open(device, os.O_RDWR) + except OSError: + return None + sizes = set() + try: + for idx in range(128): + buf = struct.pack("III", idx, yuyv, 0) + b"\x00" * 32 + try: + res = fcntl.ioctl(fd, vidioc_enum_framesizes, buf) + except OSError: + break + _i, _fmt, ftype, w, h = struct.unpack("IIIII", res[:20]) + if ftype != 1: # 1 == V4L2_FRMSIZE_TYPE_DISCRETE + break + sizes.add((w, h)) + finally: + os.close(fd) + return sizes or None + + +def _can_capture(node): + cap = cv2.VideoCapture(node, cv2.CAP_V4L2) + try: + if not cap.isOpened(): + return False + cap.set(cv2.CAP_PROP_CONVERT_RGB, 0) + ok, _ = cap.read() + return ok + finally: + cap.release() + + +# -- frame parsing -------------------------------------------------------- +def _largest_run(mask): + """Return (start, end) of the longest contiguous True run in a 1-D mask.""" + best_len = best = 0 + start = None + for i, v in enumerate(mask): + if v and start is None: + start = i + elif not v and start is not None: + if i - start > best_len: + best_len, best = i - start, (start, i) + start = None + if start is not None and len(mask) - start > best_len: + best = (start, len(mask)) + return best or None + + +def classify_bands(raw, swap=False): + """Split a raw YUYV frame into (image_band, temp_band) row ranges. + + The image band reads as a real YUYV picture (chroma near the neutral 128); + the temperature band is 16-bit data that reads as extreme chroma. Telemetry + rows (luma pinned to 0/255) are excluded. Either range may be None. + """ + chroma = raw[..., 1].mean(axis=1) + luma = raw[..., 0].mean(axis=1) + telemetry = (luma < 8) | (luma > 248) + is_image = (np.abs(chroma - 128.0) < 40) & ~telemetry + is_temp = ~is_image & ~telemetry + image_band = _largest_run(is_image) + temp_band = _largest_run(is_temp) + if swap: + image_band, temp_band = temp_band, image_band + return image_band, temp_band + + +class Frame: + """One parsed frame: a display image and (optionally) a temperature map, + both as float32 arrays at the image band's native resolution. + + `valid` is False when a known radiometric layout produced an implausible + temperature field - i.e. the camera emitted a desynced/garbage frame that the + caller should skip rather than display. + """ + + def __init__(self, image, temp, is_real, valid=True): + self.image = image # (H, W) luma + self.temp = temp # (H, W) degrees C, or None + self.is_real = is_real # True when temp is genuine radiometric + self.valid = valid + self.h, self.w = image.shape + + +def _plausible(temp): + """True for a physically sensible temperature field. Rejects desynced frames, + whose shifted image/telemetry bytes decode to a big chunk of absurd values, + while still allowing a small genuinely-hot region (a soldering iron, a flame).""" + median = float(np.median(temp)) + if not -25.0 <= median <= 160.0: + return False + garbage = float(np.mean((temp < -40.0) | (temp > 300.0))) + return garbage < 0.02 + + +def _decode_temp(raw, temp_rows, sensor, scale, order): + """Decode a temperature band into a sensor-resolution °C array, or None.""" + t0, t1 = temp_rows + sw, sh = sensor + band = np.ascontiguousarray(raw[t0:t1]).reshape(-1) + lo, hi = band[0::2].astype(np.uint16), band[1::2].astype(np.uint16) + u16 = (lo + (hi << 8)) if order == "le" else (hi + (lo << 8)) + if u16.size < sw * sh: + return None + temp = u16[:sw * sh].astype(np.float32) / float(scale) - KELVIN + temp = temp.reshape(sh, sw) + return temp if _plausible(temp) else None + + +def parse_frame(raw, scale=None, order="le", swap=False): + """Decode a raw YUYV frame into a Frame (image + optional temperature).""" + h, w, _ = raw.shape + layout = KNOWN_LAYOUTS.get((w, h)) + if layout is not None: + temp_rows, image_rows, sensor, def_scale, def_order = layout + if swap: + temp_rows, image_rows = image_rows, temp_rows + i0, i1 = image_rows + image = raw[i0:i1, :, 0].astype(np.float32) + temp = _decode_temp(raw, temp_rows, sensor, scale or def_scale, order or def_order) + if temp is None: + # Known radiometric layout but bad data this frame -> skip it. + return Frame(image, None, False, valid=False) + if temp.shape != image.shape: + temp = cv2.resize(temp, (image.shape[1], image.shape[0]), interpolation=cv2.INTER_NEAREST) + return Frame(image, temp, True) + + # Unknown resolution: fall back to chroma-based band detection (best effort). + image_band, temp_band = classify_bands(raw, swap) + if image_band is None: + return Frame(raw[..., 0].astype(np.float32), None, False) + i0, i1 = image_band + image = raw[i0:i1, :, 0].astype(np.float32) + if temp_band is None or not scale: + return Frame(image, None, False) + ih = i1 - i0 + npix = (temp_band[1] - temp_band[0]) * w + sensor = (w, ih) if ih * w <= npix else (w // 2, ih // 2) + temp = _decode_temp(raw, temp_band, sensor, scale, order) + if temp is None: + return Frame(image, None, False) + if temp.shape != image.shape: + temp = cv2.resize(temp, (image.shape[1], image.shape[0]), interpolation=cv2.INTER_NEAREST) + return Frame(image, temp, True) + + +class ThermalApp: + def __init__(self, args): + self.args = args + self.is_pi = is_raspberrypi() + + self.alpha = args.contrast + self.colormap = 0 + self.blur = 0 + self.threshold = 2.0 + self.hud = True + self.recording = False + self.video_out = None + self.rec_start = 0.0 + self.elapsed = "00:00:00" + self.snaptime = "None" + self.fullscreen = False + self.swap = args.swap_halves + self.smooth = args.smooth + self._ema = None + + self.temp_scale = None # raw-units-per-Kelvin, or None for relative + self.temp_order = args.temp_order + self.temp_offset = args.temp_offset + self.radiometric = False + + self.cap = None + self.native_w = self.native_h = 0 + self.scale = args.scale # finalised once native size is known + + # -- camera lifecycle ------------------------------------------------- + @staticmethod + def _open_at(device, w, h): + cap = cv2.VideoCapture(device, cv2.CAP_V4L2) + if not cap.isOpened(): + return None, None + cap.set(cv2.CAP_PROP_CONVERT_RGB, 0) + cap.set(cv2.CAP_PROP_FOURCC, cv2.VideoWriter_fourcc(*"YUYV")) + cap.set(cv2.CAP_PROP_FRAME_WIDTH, w) + cap.set(cv2.CAP_PROP_FRAME_HEIGHT, h) + # Without this OpenCV hands back stale buffered frames that are desynced + # mid-frame (the band layout shifts and the parse fails). + cap.set(cv2.CAP_PROP_BUFFERSIZE, 1) + # Read a few frames: the first ones after a mode change can be warm-up + # garbage that wouldn't classify correctly. + raw = None + for _ in range(5): + ok, frame = cap.read() + raw = ThermalApp._as_raw(frame) if ok else None + if raw is None or raw.shape[1] != w or raw.shape[0] != h: + cap.release() + return None, None + return cap, raw + + def _probe_radiometric(self, cap, raw, scale): + """True if this mode yields genuine temperature data. Checks the probe + frame plus a few live ones, since the occasional frame desyncs.""" + if parse_frame(raw, scale, self.temp_order, self.swap).is_real: + return True + for _ in range(10): + ok, frame = cap.read() + r = self._as_raw(frame) if ok else None + if r is not None and parse_frame(r, scale, self.temp_order, self.swap).is_real: + return True + return False + + def open_camera(self): + device = self.args.device or find_thermal_device() + if device is None: + sys.exit( + "No thermal camera found. Plug it in and check it appears in\n" + " v4l2-ctl --list-devices\n" + "then pass it explicitly, e.g. --device /dev/video2" + ) + if device.isdigit(): + device = "/dev/video" + device + + forced = self.args.temp_scale != "auto" + forced_scale = float(self.args.temp_scale) if forced else None + + # Only probe modes the camera actually advertises - opening an + # unsupported resolution desyncs the stream that follows it. + supported = supported_resolutions(device) + ok = (lambda wh: supported is None or wh in supported) + + pick = None # (cap, frame, scale, order) + if self.args.resolution: + w, h = (int(x) for x in self.args.resolution.lower().split("x")) + cap, raw = self._open_at(device, w, h) + if cap is not None: + pick = (cap, raw, forced_scale, self.temp_order) + else: + # Prefer a mode that yields genuine temperature data. Probe a few + # frames since the occasional one desyncs. + for w, h in RADIOMETRIC_MODES: + if not ok((w, h)): + continue + cap, raw = self._open_at(device, w, h) + if cap is None: + continue + scale = forced_scale if forced else KNOWN_LAYOUTS[(w, h)][3] + if self._probe_radiometric(cap, raw, scale): + pick = (cap, raw, scale, self.temp_order) + break + cap.release() + if pick is None: + for w, h in IMAGE_MODES: + if not ok((w, h)): + continue + cap, raw = self._open_at(device, w, h) + if cap is not None: + pick = (cap, raw, forced_scale, self.temp_order) + break + + if pick is None: + sys.exit(f"{device}: no usable YUYV mode. Try --resolution 256x192.") + + self.cap, raw, self.temp_scale, self.temp_order = pick + self.radiometric = self.temp_scale is not None + if not self.radiometric: + self.temp_scale = None + frame = parse_frame(raw, self.temp_scale, self.temp_order, self.swap) + self.native_h, self.native_w = frame.h, frame.w + if not self.args.scale: # auto: aim for ~768px wide + self.scale = max(1, min(5, round(768 / self.native_w))) + print(f"Opened {device}: frame {raw.shape[1]}x{raw.shape[0]}, " + f"image {self.native_w}x{self.native_h}, " + f"temperature={'REAL °C' if self.radiometric else 'relative (uncalibrated)'}") + if self.radiometric: + print(f" radiometric: scale=1/{self.temp_scale:g} K, order={self.temp_order}, " + f"offset={self.temp_offset:+.1f} °C (adjust with [ and ])") + else: + print(" no 16-bit data in this stream - see docs/TC002C-DUO.md") + + @staticmethod + def _as_raw(frame): + if frame is None: + return None + if frame.ndim == 3 and frame.shape[2] == 2: + return frame + if frame.ndim == 2 and frame.shape[1] % 2 == 0: + return frame.reshape(frame.shape[0], frame.shape[1] // 2, 2) + return None + + # -- per-frame processing -------------------------------------------- + @staticmethod + def _destripe(luma): + """Remove per-row and per-column fixed-pattern offsets while keeping + smooth gradients (uncorrected previews carry 2D FPN that a stretch + turns into stripes; we can't run the camera's shutter NUC over UVC).""" + row = np.median(luma, axis=1) + luma = luma - (row - cv2.blur(row.reshape(-1, 1), (1, 9)).ravel())[:, None] + col = np.median(luma, axis=0) + return luma - (col - cv2.blur(col.reshape(1, -1), (9, 1)).ravel())[None, :] + + def _display_luma(self, image): + """Turn the raw image band into a stretched, denoised 8-bit luma.""" + luma = image.copy() + if self.smooth > 0: + if self._ema is None or self._ema.shape != luma.shape: + self._ema = luma + else: + self._ema = self.smooth * self._ema + (1.0 - self.smooth) * luma + luma = self._ema + if not self.args.no_destripe: + luma = self._destripe(luma) + if not self.args.no_stretch: + lo, hi = np.percentile(luma, [1, 99]) + span = max(hi - lo, 45.0) + luma = np.clip((luma - lo) * (255.0 / span), 0, 255) + return luma.astype(np.uint8) + + def render(self, raw): + frame = parse_frame(raw, self.temp_scale, self.temp_order, self.swap) + if not frame.valid: + return None # desynced/garbage frame - caller reuses the last one + is_real = frame.is_real + if is_real: + temp_map = frame.temp + self.temp_offset + unit = "C" + else: + temp_map = self.args.rel_gain * frame.image + self.args.rel_offset + unit = "lvl" + + h, w = frame.h, frame.w + center = temp_map[h // 2, w // 2] + inner = temp_map[2:-2, 2:-2] + max_pos = tuple(p + 2 for p in np.unravel_index(np.argmax(inner), inner.shape)) + min_pos = tuple(p + 2 for p in np.unravel_index(np.argmin(inner), inner.shape)) + tmax, tmin = float(temp_map[max_pos]), float(temp_map[min_pos]) + tavg = float(temp_map.mean()) + + luma = self._display_luma(frame.image) + luma = cv2.convertScaleAbs(luma, alpha=self.alpha) + new_w, new_h = w * self.scale, h * self.scale + luma = cv2.resize(luma, (new_w, new_h), interpolation=cv2.INTER_CUBIC) + if self.blur > 0: + luma = cv2.blur(luma, (self.blur, self.blur)) + + cmap, cmap_name = COLORMAPS[self.colormap] + heatmap = cv2.applyColorMap(luma, cmap) + + self._draw_crosshair(heatmap, new_w, new_h, center, unit) + self._draw_markers(heatmap, max_pos, min_pos, tmax, tmin, tavg, unit, w, h, new_w, new_h) + if self.hud: + self._draw_hud(heatmap, cmap_name, tavg, unit, is_real) + return heatmap + + def _draw_crosshair(self, img, w, h, center_temp, unit): + cx, cy = w // 2, h // 2 + for color, thick in (((255, 255, 255), 2), ((0, 0, 0), 1)): + cv2.line(img, (cx, cy + 20), (cx, cy - 20), color, thick) + cv2.line(img, (cx + 20, cy), (cx - 20, cy), color, thick) + label = f"{center_temp:.1f} {unit}" + cv2.putText(img, label, (cx + 10, cy - 10), cv2.FONT_HERSHEY_SIMPLEX, 0.45, (0, 0, 0), 2, cv2.LINE_AA) + cv2.putText(img, label, (cx + 10, cy - 10), cv2.FONT_HERSHEY_SIMPLEX, 0.45, (0, 255, 255), 1, cv2.LINE_AA) + + def _draw_markers(self, img, max_pos, min_pos, tmax, tmin, tavg, unit, w, h, nw, nh): + sx, sy = nw / w, nh / h + for (row, col), col_bgr, value in ((max_pos, (0, 0, 255), tmax), (min_pos, (255, 0, 0), tmin)): + if abs(value - tavg) < self.threshold: + continue + x, y = int(col * sx), int(row * sy) + cv2.circle(img, (x, y), 5, (0, 0, 0), 2) + cv2.circle(img, (x, y), 5, col_bgr, -1) + text = f"{value:.1f} {unit}" + cv2.putText(img, text, (x + 10, y + 5), cv2.FONT_HERSHEY_SIMPLEX, 0.45, (0, 0, 0), 2, cv2.LINE_AA) + cv2.putText(img, text, (x + 10, y + 5), cv2.FONT_HERSHEY_SIMPLEX, 0.45, (0, 255, 255), 1, cv2.LINE_AA) + + def _draw_hud(self, img, cmap_name, tavg, unit, is_real): + cv2.rectangle(img, (0, 0), (180, 130), (0, 0, 0), -1) + rows = [ + f"Avg: {tavg:.1f} {unit}", + f"Threshold: {self.threshold:.0f}", + f"Colormap: {cmap_name}", + f"Blur: {self.blur} Scale: {self.scale}", + f"Contrast: {self.alpha:.1f}", + f"Snapshot: {self.snaptime}", + ] + for i, text in enumerate(rows): + cv2.putText(img, text, (8, 14 + i * 14), cv2.FONT_HERSHEY_SIMPLEX, 0.4, (0, 255, 255), 1, cv2.LINE_AA) + mode = "RADIOMETRIC" if is_real else "RELATIVE (uncal)" + mode_color = (0, 255, 0) if is_real else (0, 165, 255) + cv2.putText(img, mode, (8, 14 + len(rows) * 14), cv2.FONT_HERSHEY_SIMPLEX, 0.4, mode_color, 1, cv2.LINE_AA) + rec_color = (40, 40, 255) if self.recording else (200, 200, 200) + cv2.putText(img, f"REC {self.elapsed}", (8, 14 + (len(rows) + 1) * 14), cv2.FONT_HERSHEY_SIMPLEX, 0.4, rec_color, 1, cv2.LINE_AA) + + # -- recording / snapshots ------------------------------------------- + def start_recording(self, frame_size): + now = time.strftime("%Y%m%d--%H%M%S") + self.video_out = cv2.VideoWriter(now + "-output.avi", cv2.VideoWriter_fourcc(*"XVID"), 25, frame_size) + self.recording = True + self.rec_start = time.time() + + def stop_recording(self): + self.recording = False + self.elapsed = "00:00:00" + if self.video_out is not None: + self.video_out.release() + self.video_out = None + + def snapshot(self, heatmap): + now = time.strftime("%Y%m%d-%H%M%S") + cv2.imwrite(f"thermal-{now}.png", heatmap) + self.snaptime = time.strftime("%H:%M:%S") + print(f"Saved thermal-{now}.png") + + # -- main loop -------------------------------------------------------- + def run(self): + self.open_camera() + win = "Thermal" + cv2.namedWindow(win, cv2.WINDOW_GUI_NORMAL) + cv2.resizeWindow(win, self.native_w * self.scale, self.native_h * self.scale) + self._print_keys() + while True: + ok, frame = self.cap.read() + if not ok or frame is None: + time.sleep(0.05) + continue + raw = self._as_raw(frame) + if raw is None: + continue + heatmap = self.render(raw) + if heatmap is None: # desynced frame; show the previous one + continue + cv2.imshow(win, heatmap) + if self.recording and self.video_out is not None: + self.elapsed = time.strftime("%H:%M:%S", time.gmtime(time.time() - self.rec_start)) + self.video_out.write(heatmap) + if self.handle_key(cv2.waitKey(1) & 0xFF, heatmap): + break + self.cleanup() + + def handle_key(self, key, heatmap): + if key in (ord("q"), 27): + return True + elif key == ord("a"): + self.blur += 1 + elif key == ord("z"): + self.blur = max(0, self.blur - 1) + elif key == ord("s"): + self.threshold += 1 + elif key == ord("x"): + self.threshold = max(0, self.threshold - 1) + elif key == ord("d"): + self.scale = min(5, self.scale + 1) + self._resize_window() + elif key == ord("c"): + self.scale = max(1, self.scale - 1) + self._resize_window() + elif key == ord("f"): + self.alpha = min(3.0, round(self.alpha + 0.1, 1)) + elif key == ord("v"): + self.alpha = max(0.0, round(self.alpha - 0.1, 1)) + elif key == ord("m"): + self.colormap = (self.colormap + 1) % len(COLORMAPS) + elif key == ord("h"): + self.hud = not self.hud + elif key == ord("n"): + self.swap = not self.swap + elif key == ord("g"): + self.smooth = 0.0 if self.smooth > 0 else (self.args.smooth or 0.5) + elif key == ord("w"): + self._set_fullscreen(False) + elif key == ord("e"): + self._set_fullscreen(True) + elif key == ord("p"): + self.snapshot(heatmap) + elif key == ord("r") and not self.recording: + self.start_recording((heatmap.shape[1], heatmap.shape[0])) + elif key == ord("t"): + self.stop_recording() + elif key == ord("]"): # calibration: nudge temperature offset / relative gain + if self.radiometric: + self.temp_offset = round(self.temp_offset + 0.5, 1) + else: + self.args.rel_gain = round(self.args.rel_gain + 0.05, 3) + elif key == ord("["): + if self.radiometric: + self.temp_offset = round(self.temp_offset - 0.5, 1) + else: + self.args.rel_gain = round(self.args.rel_gain - 0.05, 3) + return False + + def _resize_window(self): + if not self.fullscreen and not self.is_pi: + cv2.resizeWindow("Thermal", self.native_w * self.scale, self.native_h * self.scale) + + def _set_fullscreen(self, on): + self.fullscreen = on + cv2.setWindowProperty("Thermal", cv2.WND_PROP_FULLSCREEN, + cv2.WINDOW_FULLSCREEN if on else cv2.WINDOW_NORMAL) + if not on: + self._resize_window() + + def cleanup(self): + if self.recording: + self.stop_recording() + if self.cap is not None: + self.cap.release() + cv2.destroyAllWindows() + + def _print_keys(self): + print( + "\nKey bindings:\n" + " a/z blur +/- s/x min-max label threshold +/-\n" + " d/c scale +/- f/v contrast +/-\n" + " m cycle colormap h toggle HUD\n" + " n swap bands g toggle temporal smoothing\n" + " e/w fullscreen on/off\n" + " r/t record / stop p snapshot\n" + " [/] temperature offset (or relative gain) calibration\n" + " q/ESC quit\n" + ) + + # -- headless self-test ---------------------------------------------- + def selftest(self, frames): + self.open_camera() + last = None + for i in range(frames): + ok, frame = self.cap.read() + if not ok or frame is None: + continue + raw = self._as_raw(frame) + if raw is None: + continue + f = parse_frame(raw, self.temp_scale, self.temp_order, self.swap) + rendered = self.render(raw) + if rendered is not None: + last = rendered + if i == frames - 1 and f.valid: + if f.is_real: + t = f.temp + self.temp_offset + print(f"frame {i}: REAL temps center={t[f.h//2, f.w//2]:.1f}C " + f"min={t.min():.1f} max={t.max():.1f}") + else: + print(f"frame {i}: relative (no radiometric data)") + self.cap.release() + if last is None: + print("SELFTEST FAILED: no frames decoded") + return 1 + cv2.imwrite("selftest-thermal.png", last) + print(f"SELFTEST OK: wrote selftest-thermal.png ({last.shape[1]}x{last.shape[0]})") + return 0 + + +def parse_args(argv=None): + p = argparse.ArgumentParser(description="Universal Linux thermal camera viewer (TC001, TC002C Duo, InfiRay).") + p.add_argument("--device", help="Video device, e.g. /dev/video2 or 2. Default: auto-detect.") + p.add_argument("--resolution", help="Force capture WxH, e.g. 512x484. Default: best radiometric mode.") + p.add_argument("--scale", type=int, default=0, help="Display upscaling 1-5 (0 = auto from image size).") + p.add_argument("--contrast", type=float, default=1.0, help="Initial contrast (0.0-3.0).") + p.add_argument("--temp-scale", default="auto", + help="Radiometric raw-units-per-Kelvin divisor: 'auto', 64 (TC001), or 16 (Duo).") + p.add_argument("--temp-order", default="le", choices=["le", "be"], help="16-bit byte order.") + p.add_argument("--temp-offset", type=float, default=0.0, + help="°C added to decoded temperatures (calibrate against a known reference; live keys [ ]).") + p.add_argument("--rel-gain", type=float, default=0.2, help="Relative-mode level->pseudo-temp gain.") + p.add_argument("--rel-offset", type=float, default=0.0, help="Relative-mode offset.") + p.add_argument("--no-stretch", action="store_true", help="Disable the display contrast stretch.") + p.add_argument("--no-destripe", action="store_true", help="Disable fixed-pattern-noise removal.") + p.add_argument("--smooth", type=float, default=0.5, help="Temporal smoothing 0.0-0.9 (0 disables). Live: 'g'.") + p.add_argument("--swap-halves", action="store_true", help="Swap which band is image vs data. Live: 'n'.") + p.add_argument("--selftest", type=int, metavar="N", help="Headless: grab N frames, save a snapshot, exit.") + return p.parse_args(argv) + + +def main(argv=None): + args = parse_args(argv) + app = ThermalApp(args) + if args.selftest: + return app.selftest(args.selftest) + try: + app.run() + except KeyboardInterrupt: + app.cleanup() + print("\nBye.") + return 0 + + +if __name__ == "__main__": + sys.exit(main()) From 33fb0ff1aba133c363339220c511212522125530 Mon Sep 17 00:00:00 2001 From: Aron Novak Date: Thu, 11 Jun 2026 10:55:44 +0200 Subject: [PATCH 2/5] Add MJPEG bridge, orientation controls, reusable frame API - src/thermal_bridge.py: serve the parsed/colormapped thermal image as an MJPEG-over-HTTP stream so any app (OpenCV, ffmpeg, browser) can consume it like a webcam, no dependency. - thermalcam: --rotate / --flip (live 'o' key) and a clean colormap_frame() method (overlay-free) used by the bridge. Co-Authored-By: Claude Opus 4.8 (1M context) --- README.md | 26 ++++++ src/thermal_bridge.py | 178 ++++++++++++++++++++++++++++++++++++++++++ src/thermalcam.py | 50 ++++++++++-- 3 files changed, 249 insertions(+), 5 deletions(-) create mode 100644 src/thermal_bridge.py diff --git a/README.md b/README.md index cb25716..692723b 100644 --- a/README.md +++ b/README.md @@ -54,6 +54,32 @@ Useful options: `--temp-offset N` to correct the absolute °C, `--temp-scale · `h` HUD · `n` swap bands · `g` temporal smoothing · `e/w` fullscreen on/off · `r/t` record/stop · `p` snapshot · `[`/`]` temperature-offset calibration · `q`/ESC quit +### Feed it to other apps (MJPEG bridge) + +`src/thermal_bridge.py` re-serves the parsed, colormapped thermal image as an +MJPEG-over-HTTP stream, so any app can consume it like a webcam — no kernel +module, no root, no dependency: + +```bash +python3 src/thermal_bridge.py --colormap inferno --rotate 90 # serves http://127.0.0.1:8090 +``` + +```python +import cv2 +cap = cv2.VideoCapture("http://127.0.0.1:8090/stream.mjpg") # OpenCV, ffmpeg, browsers, ... +``` + +Options mirror the viewer (`--colormap`, `--rotate`, `--flip`, `--width/--height`, +`--temp-scale`, `--device`, …). Open `http://127.0.0.1:8090/` in a browser for a +live preview. + +### Orientation + +The sensor is landscape (4:3). Rotate/mirror in software with `--rotate +{0,90,180,270}` and `--flip {none,h,v}`, or live with the `o` key. (The InfiRay +sensor also has a hardware mirror/flip, reachable only via the vendor command +protocol — software rotation is simpler.) + ### Camera support | Camera | Image | Temperature | diff --git a/src/thermal_bridge.py b/src/thermal_bridge.py new file mode 100644 index 0000000..67c8ceb --- /dev/null +++ b/src/thermal_bridge.py @@ -0,0 +1,178 @@ +#!/usr/bin/env python3 +""" +thermal_bridge.py - Serve the thermal camera as an MJPEG-over-HTTP stream. + +Re-exposes an InfiRay/Topdon thermal camera (parsed + colormapped by thermalcam) +as a plain MJPEG stream that ANY app can open like a webcam - no kernel module, +no root, no reboot, no code dependency: + + python3 src/thermal_bridge.py --colormap inferno --rotate 90 + + # then, in any consumer: + OpenCV: cv2.VideoCapture("http://127.0.0.1:8090/stream.mjpg") + ffmpeg: ffmpeg -i http://127.0.0.1:8090/stream.mjpg ... + browser: http://127.0.0.1:8090/ (preview page) + +This exists because the camera's default UVC output is a torn/garbage frame; the +bridge does the 512x484 radiometric parse and hands out a clean thermal image. +""" + +import argparse +import sys +import threading +from http.server import BaseHTTPRequestHandler, ThreadingHTTPServer + +import cv2 + +import thermalcam + +BOUNDARY = "frame" + +# Latest encoded JPEG, published by the capture thread and consumed by clients. +_latest = {"jpeg": None} +_new_frame = threading.Condition() +_stop = threading.Event() + + +def capture_loop(app, out_size, quality): + """Grab frames, colormap them, publish the newest JPEG. Skips desynced frames.""" + encode = [cv2.IMWRITE_JPEG_QUALITY, quality] + while not _stop.is_set(): + ok, frame = app.cap.read() + raw = app._as_raw(frame) if ok else None + if raw is None: + continue + bgr = app.colormap_frame(raw) + if bgr is None: # desynced frame; keep the previous one + continue + if out_size is not None: + bgr = cv2.resize(bgr, out_size, interpolation=cv2.INTER_AREA) + ok, buf = cv2.imencode(".jpg", bgr, encode) + if not ok: + continue + with _new_frame: + _latest["jpeg"] = buf.tobytes() + _new_frame.notify_all() + + +class Handler(BaseHTTPRequestHandler): + protocol_version = "HTTP/1.0" + + def log_message(self, *a): # quiet + pass + + def do_GET(self): + if self.path.startswith("/stream") or self.path.rstrip("/") == "": + if self.path.rstrip("/") == "": + return self._page() + return self._stream() + if self.path.startswith("/snapshot"): + return self._snapshot() + self.send_error(404) + + def _page(self): + html = (b"" + b"") + self.send_response(200) + self.send_header("Content-Type", "text/html") + self.send_header("Content-Length", str(len(html))) + self.end_headers() + self.wfile.write(html) + + def _snapshot(self): + with _new_frame: + _new_frame.wait(timeout=2.0) + jpeg = _latest["jpeg"] + if jpeg is None: + return self.send_error(503) + self.send_response(200) + self.send_header("Content-Type", "image/jpeg") + self.send_header("Content-Length", str(len(jpeg))) + self.end_headers() + self.wfile.write(jpeg) + + def _stream(self): + self.send_response(200) + self.send_header("Content-Type", f"multipart/x-mixed-replace; boundary={BOUNDARY}") + self.send_header("Cache-Control", "no-cache") + self.end_headers() + try: + while not _stop.is_set(): + with _new_frame: + _new_frame.wait(timeout=2.0) + jpeg = _latest["jpeg"] + if jpeg is None: + continue + self.wfile.write(b"--" + BOUNDARY.encode() + b"\r\n") + self.wfile.write(b"Content-Type: image/jpeg\r\n") + self.wfile.write(f"Content-Length: {len(jpeg)}\r\n\r\n".encode()) + self.wfile.write(jpeg) + self.wfile.write(b"\r\n") + except (BrokenPipeError, ConnectionResetError): + pass + + +def build_app(args): + """Open the camera through thermalcam with the requested settings.""" + names = [n.lower() for _c, n in thermalcam.COLORMAPS] + if args.colormap.lower() not in names: + sys.exit(f"--colormap must be one of: {', '.join(names)}") + + tc_argv = ["--scale", "1", "--rotate", str(args.rotate), "--flip", args.flip, + "--temp-scale", args.temp_scale, "--temp-order", args.temp_order] + if args.device: + tc_argv += ["--device", args.device] + if args.resolution: + tc_argv += ["--resolution", args.resolution] + if args.no_destripe: + tc_argv += ["--no-destripe"] + tc_argv += ["--smooth", str(args.smooth)] + + app = thermalcam.ThermalApp(thermalcam.parse_args(tc_argv)) + app.colormap = names.index(args.colormap.lower()) + app.open_camera() + return app + + +def main(argv=None): + p = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) + p.add_argument("--host", default="127.0.0.1", help="Bind address (0.0.0.0 to expose on the LAN).") + p.add_argument("--port", type=int, default=8090) + p.add_argument("--colormap", default="inferno", help="Palette (jet, inferno, hot, ...). Default inferno.") + p.add_argument("--width", type=int, help="Output width (default: camera native).") + p.add_argument("--height", type=int, help="Output height (default: camera native).") + p.add_argument("--quality", type=int, default=85, help="JPEG quality 1-100.") + p.add_argument("--device", help="Video device. Default: auto-detect.") + p.add_argument("--resolution", help="Force capture WxH.") + p.add_argument("--rotate", type=int, default=0, choices=[0, 90, 180, 270]) + p.add_argument("--flip", default="none", choices=["none", "h", "v"]) + p.add_argument("--temp-scale", default="auto") + p.add_argument("--temp-order", default="le", choices=["le", "be"]) + p.add_argument("--no-destripe", action="store_true") + p.add_argument("--smooth", type=float, default=0.5) + args = p.parse_args(argv) + + app = build_app(args) + out_size = (args.width, args.height) if args.width and args.height else None + + worker = threading.Thread(target=capture_loop, args=(app, out_size, args.quality), daemon=True) + worker.start() + + server = ThreadingHTTPServer((args.host, args.port), Handler) + url = f"http://{args.host}:{args.port}" + print(f"Thermal MJPEG stream at {url}/stream.mjpg (preview: {url}/ snapshot: {url}/snapshot.jpg)") + print("Open it from any app, e.g.:") + print(f" python3 -c \"import cv2; c=cv2.VideoCapture('{url}/stream.mjpg'); print(c.read()[0])\"") + try: + server.serve_forever() + except KeyboardInterrupt: + pass + finally: + _stop.set() + server.shutdown() + app.cap.release() + return 0 + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/src/thermalcam.py b/src/thermalcam.py index 2dc43a7..ba43f62 100644 --- a/src/thermalcam.py +++ b/src/thermalcam.py @@ -288,6 +288,8 @@ def __init__(self, args): self.swap = args.swap_halves self.smooth = args.smooth self._ema = None + self.rotate = args.rotate # 0/90/180/270 + self.flip = args.flip # none/h/v self.temp_scale = None # raw-units-per-Kelvin, or None for relative self.temp_order = args.temp_order @@ -440,19 +442,36 @@ def _display_luma(self, image): luma = np.clip((luma - lo) * (255.0 / span), 0, 255) return luma.astype(np.uint8) + def _orient(self, arr): + """Apply the current rotation/flip to an image or temperature array.""" + if arr is None: + return None + if self.rotate == 90: + arr = cv2.rotate(arr, cv2.ROTATE_90_CLOCKWISE) + elif self.rotate == 180: + arr = cv2.rotate(arr, cv2.ROTATE_180) + elif self.rotate == 270: + arr = cv2.rotate(arr, cv2.ROTATE_90_COUNTERCLOCKWISE) + if self.flip == "h": + arr = cv2.flip(arr, 1) + elif self.flip == "v": + arr = cv2.flip(arr, 0) + return arr + def render(self, raw): frame = parse_frame(raw, self.temp_scale, self.temp_order, self.swap) if not frame.valid: return None # desynced/garbage frame - caller reuses the last one is_real = frame.is_real + image = self._orient(frame.image) if is_real: - temp_map = frame.temp + self.temp_offset + temp_map = self._orient(frame.temp) + self.temp_offset unit = "C" else: - temp_map = self.args.rel_gain * frame.image + self.args.rel_offset + temp_map = self.args.rel_gain * image + self.args.rel_offset unit = "lvl" - h, w = frame.h, frame.w + h, w = image.shape center = temp_map[h // 2, w // 2] inner = temp_map[2:-2, 2:-2] max_pos = tuple(p + 2 for p in np.unravel_index(np.argmax(inner), inner.shape)) @@ -460,7 +479,7 @@ def render(self, raw): tmax, tmin = float(temp_map[max_pos]), float(temp_map[min_pos]) tavg = float(temp_map.mean()) - luma = self._display_luma(frame.image) + luma = self._display_luma(image) luma = cv2.convertScaleAbs(luma, alpha=self.alpha) new_w, new_h = w * self.scale, h * self.scale luma = cv2.resize(luma, (new_w, new_h), interpolation=cv2.INTER_CUBIC) @@ -476,6 +495,22 @@ def render(self, raw): self._draw_hud(heatmap, cmap_name, tavg, unit, is_real) return heatmap + def colormap_frame(self, raw): + """A clean colormapped BGR frame (oriented, no crosshair/HUD/markers), + for feeding other apps. Returns None for a desynced frame.""" + frame = parse_frame(raw, self.temp_scale, self.temp_order, self.swap) + if not frame.valid: + return None + image = self._orient(frame.image) + luma = self._display_luma(image) + luma = cv2.convertScaleAbs(luma, alpha=self.alpha) + nw, nh = image.shape[1] * self.scale, image.shape[0] * self.scale + luma = cv2.resize(luma, (nw, nh), interpolation=cv2.INTER_CUBIC) + if self.blur > 0: + luma = cv2.blur(luma, (self.blur, self.blur)) + cmap, _ = COLORMAPS[self.colormap] + return cv2.applyColorMap(luma, cmap) + def _draw_crosshair(self, img, w, h, center_temp, unit): cx, cy = w // 2, h // 2 for color, thick in (((255, 255, 255), 2), ((0, 0, 0), 1)): @@ -590,6 +625,8 @@ def handle_key(self, key, heatmap): self.swap = not self.swap elif key == ord("g"): self.smooth = 0.0 if self.smooth > 0 else (self.args.smooth or 0.5) + elif key == ord("o"): # cycle rotation 0->90->180->270 + self.rotate = (self.rotate + 90) % 360 elif key == ord("w"): self._set_fullscreen(False) elif key == ord("e"): @@ -637,7 +674,7 @@ def _print_keys(self): " d/c scale +/- f/v contrast +/-\n" " m cycle colormap h toggle HUD\n" " n swap bands g toggle temporal smoothing\n" - " e/w fullscreen on/off\n" + " o rotate 90 deg e/w fullscreen on/off\n" " r/t record / stop p snapshot\n" " [/] temperature offset (or relative gain) calibration\n" " q/ESC quit\n" @@ -691,6 +728,9 @@ def parse_args(argv=None): p.add_argument("--no-destripe", action="store_true", help="Disable fixed-pattern-noise removal.") p.add_argument("--smooth", type=float, default=0.5, help="Temporal smoothing 0.0-0.9 (0 disables). Live: 'g'.") p.add_argument("--swap-halves", action="store_true", help="Swap which band is image vs data. Live: 'n'.") + p.add_argument("--rotate", type=int, default=0, choices=[0, 90, 180, 270], + help="Rotate the image clockwise. Live: cycle with 'o'.") + p.add_argument("--flip", default="none", choices=["none", "h", "v"], help="Mirror the image horizontally/vertically.") p.add_argument("--selftest", type=int, metavar="N", help="Headless: grab N frames, save a snapshot, exit.") return p.parse_args(argv) From 2f49a730402198ae9307cfba6ecfc5abdfef3294 Mon Sep 17 00:00:00 2001 From: Aron Novak Date: Thu, 11 Jun 2026 11:12:08 +0200 Subject: [PATCH 3/5] Reject desynced/tiled frames; auto-recover the stream These cameras occasionally emit a horizontally-tiled/desynced frame. The frame validity check only looked at the temperature band, so such frames slipped through as garbage. Add an image-band alignment guard (mean-abs-difference vs a shifted copy; gradients and flat scenes are exempt) so the viewer and the MJPEG bridge never emit them. The bridge also reopens the camera after a run of unusable frames to recover from a stuck desync. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/thermal_bridge.py | 21 +++++++++++++-------- src/thermalcam.py | 38 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 51 insertions(+), 8 deletions(-) diff --git a/src/thermal_bridge.py b/src/thermal_bridge.py index 67c8ceb..a51f759 100644 --- a/src/thermal_bridge.py +++ b/src/thermal_bridge.py @@ -35,16 +35,21 @@ def capture_loop(app, out_size, quality): - """Grab frames, colormap them, publish the newest JPEG. Skips desynced frames.""" + """Grab frames, colormap them, publish the newest JPEG. Skips desynced/tiled + frames, and reopens the camera if the stream wedges into a stuck desync.""" encode = [cv2.IMWRITE_JPEG_QUALITY, quality] + skips = 0 while not _stop.is_set(): ok, frame = app.cap.read() raw = app._as_raw(frame) if ok else None - if raw is None: - continue - bgr = app.colormap_frame(raw) - if bgr is None: # desynced frame; keep the previous one + bgr = app.colormap_frame(raw) if raw is not None else None + if bgr is None: # bad frame; keep serving the previous one + skips += 1 + if skips % 60 == 0: # ~2-4 s of nothing usable -> recover + print(f"[bridge] {skips} unusable frames, reopening camera...", flush=True) + app.reopen() continue + skips = 0 if out_size is not None: bgr = cv2.resize(bgr, out_size, interpolation=cv2.INTER_AREA) ok, buf = cv2.imencode(".jpg", bgr, encode) @@ -160,9 +165,9 @@ def main(argv=None): server = ThreadingHTTPServer((args.host, args.port), Handler) url = f"http://{args.host}:{args.port}" - print(f"Thermal MJPEG stream at {url}/stream.mjpg (preview: {url}/ snapshot: {url}/snapshot.jpg)") - print("Open it from any app, e.g.:") - print(f" python3 -c \"import cv2; c=cv2.VideoCapture('{url}/stream.mjpg'); print(c.read()[0])\"") + print(f"Thermal MJPEG stream at {url}/stream.mjpg (preview: {url}/ snapshot: {url}/snapshot.jpg)", flush=True) + print("Open it from any app, e.g.:", flush=True) + print(f" python3 -c \"import cv2; c=cv2.VideoCapture('{url}/stream.mjpg'); print(c.read()[0])\"", flush=True) try: server.serve_forever() except KeyboardInterrupt: diff --git a/src/thermalcam.py b/src/thermalcam.py index ba43f62..06cb524 100644 --- a/src/thermalcam.py +++ b/src/thermalcam.py @@ -218,6 +218,27 @@ def _plausible(temp): return garbage < 0.02 +def _image_aligned(luma): + """False if the image band looks horizontally tiled - the signature of a + desynced frame whose line stride is wrong (content repeats at width/4 or /2). + + Tested by mean-absolute-difference against a shifted copy: tiled panels are + near-identical (diff ~0), whereas a normal scene - even a smooth left-to-right + temperature gradient - differs substantially when shifted. Flat/near-uniform + scenes are exempt.""" + s = float(luma.std()) + if s < 8.0: # too uniform to judge; treat as fine + return True + w = luma.shape[1] + for period in (w // 4, w // 2): + if period < 2: + continue + diff = float(np.abs(luma[:, period:] - luma[:, :-period]).mean()) + if diff < max(3.0, 0.12 * s): # panels ~identical at this period => tiled + return False + return True + + def _decode_temp(raw, temp_rows, sensor, scale, order): """Decode a temperature band into a sensor-resolution °C array, or None.""" t0, t1 = temp_rows @@ -242,6 +263,8 @@ def parse_frame(raw, scale=None, order="le", swap=False): temp_rows, image_rows = image_rows, temp_rows i0, i1 = image_rows image = raw[i0:i1, :, 0].astype(np.float32) + if not _image_aligned(image): + return Frame(image, None, False, valid=False) temp = _decode_temp(raw, temp_rows, sensor, scale or def_scale, order or def_order) if temp is None: # Known radiometric layout but bad data this frame -> skip it. @@ -388,6 +411,8 @@ def open_camera(self): sys.exit(f"{device}: no usable YUYV mode. Try --resolution 256x192.") self.cap, raw, self.temp_scale, self.temp_order = pick + self._device = device + self._wh = (raw.shape[1], raw.shape[0]) # for reopen() after a desync self.radiometric = self.temp_scale is not None if not self.radiometric: self.temp_scale = None @@ -404,6 +429,19 @@ def open_camera(self): else: print(" no 16-bit data in this stream - see docs/TC002C-DUO.md") + def reopen(self): + """Reopen the camera at the same resolution to recover from a stuck + desync (the stream can wedge into a misaligned state).""" + try: + self.cap.release() + except Exception: + pass + cap, _raw = self._open_at(self._device, *self._wh) + if cap is not None: + self.cap = cap + return True + return False + @staticmethod def _as_raw(frame): if frame is None: From 55f5a3416dcfe48c9fc64edf4db790032938dac8 Mon Sep 17 00:00:00 2001 From: Aron Novak Date: Tue, 16 Jun 2026 07:02:16 +0200 Subject: [PATCH 4/5] radiometry: emissivity correction + per-unit calibration profile MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Reproduce the InfiRay temperature-correction math on Linux from the Android libadvirtemp tables, so the camera can report true object °C off-device: - tables.py extracts the embedded float64 ems/distance/target-temp LUTs into adv_tables.npz via a self-contained ELF reader. - corrector.py: graybody emissivity / reflected-temp inversion (temp_correct), round-trip validated <1e-6 °C against the decoded inverse; temperature-dependent effective-emissivity LUT; Magnus vapour pressure. - profile.py: per-unit JSON profile (keyed by USB serial) with scene defaults and NUC/kt/bt slots filled at provisioning. Vendor-protocol read of per-unit NUC/kt/bt is still to come. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/radiometry/__init__.py | 12 +++ src/radiometry/adv_tables.npz | Bin 0 -> 26899 bytes src/radiometry/corrector.py | 186 ++++++++++++++++++++++++++++++++++ src/radiometry/profile.py | 100 ++++++++++++++++++ src/radiometry/tables.py | 146 ++++++++++++++++++++++++++ 5 files changed, 444 insertions(+) create mode 100644 src/radiometry/__init__.py create mode 100644 src/radiometry/adv_tables.npz create mode 100644 src/radiometry/corrector.py create mode 100644 src/radiometry/profile.py create mode 100644 src/radiometry/tables.py diff --git a/src/radiometry/__init__.py b/src/radiometry/__init__.py new file mode 100644 index 0000000..4f303cc --- /dev/null +++ b/src/radiometry/__init__.py @@ -0,0 +1,12 @@ +"""Radiometric calibration for InfiRay/Topdon UVC cameras (e.g. TC002C Duo). + +Reproduces the camera's factory true-temperature pipeline on Linux/Pi by +reimplementing the math from the Android native libraries in pure numpy. The +native libs are arm64/Bionic/JNI and cannot be loaded off-device, so the model's +correction lookup tables are extracted once (`tables.py` -> `adv_tables.npz`) and +the per-pixel correction is reimplemented (`corrector.py`). Per-unit factory +parameters (NUC / kt / bt) are read off each camera over the USB vendor protocol +and stored in a per-unit profile (`profile.py`). + +See docs/TC002C-DUO.md for the reverse-engineering background. +""" diff --git a/src/radiometry/adv_tables.npz b/src/radiometry/adv_tables.npz new file mode 100644 index 0000000000000000000000000000000000000000..111d05f70169a88035546627e86d5f17405bc1b5 GIT binary patch literal 26899 zcmeFWd0bL$+cw&)Woc$*FTzt6Mxx4-?z_x<olyzK|Mv)3{> zG;035mGd^v*IcT3)&FZVFLTOp(r9JVkF$+ui*lY{$?|#~T<|J5`9^MW{;Qlc$c6JM zFAQ;zl;_Jp(1V-fdhAWGRW zn8ukKSBYgKn1~KL~d~7p){fILch-h8+wAp&LYop)A zkSSWu{21$5;nnz%X`Q+y#4_;mecU!r#HW4F;Ijc&`Bc;E9z#j+*{gws!M@lhiBTbX zm%egq9Oj+rax38Kx>Xzz-@&54f-8dz{FZ2D>#7oOd-S*z$WOe+=S6=TM zq4lGS?tZNk;;f*a8x)254Ep%F)u<2~%$qyk&L2aHvV>CD#jM^njJ^2pfGOF?=`)HW zHDn^yR3x7&f!7l2OlLO9w3GwPJ|v0KUAO^q+?3>F$+8~3_By0S>8 zBxd1hpGTR$oJJpk?yvLZz62%hfjR;bZVq&i=Wp*5mtY%BP|No=_B)_7*TjHy+oNA& z7oG;SlR+Cy=GOh$ALes!qxO04F#SO3BY-!{@ZaG-K71(VANc<#DY@^TlzjC!^!FD0 zYYMWTE6oS4t!oSlDC&H^)e0e}-cB&iEm+anb>nMVOX>FHgHy&Q?`Qko7);x``{eVU z?I)d`cb{~=^1P=i=tjuN=UXS_)!T~-3c8_ky;+467BL`URdT7$?pBtNzF~K!I|qh9 zTGvpwQlcu^DV=^ciOn^sgZDlgF@4I(3LJ39Oq5qzr@PnL;jrzOp9#Xkig%=63{gs; z_EEQUouy#|@76#|m1$a`!D6{n}I82AZ+khtJ+5 z_@40oA?SW;OpbN;9!Jme^(~>ox(r3v<@?N+L%BMX<)wqReWq2(u~YXBomfi^yJt%N z^ILIa@r9(cr?L#N`=ORQgy5w3*89%|r;={CY<#e7Q^L8T<)Qn}xJU}0gnxT1p7iG( zS|0D4lwU%b2R)j+iHyJHXc0=fUj5nX!p2>+r`n1iJC?VcT>Z9i#p0~}9z&URyY=55 zSa)Cnxh1-20i^YpX$)fZtJ-r{&JI!hMT@vD;`r~D%Z%60u7#h7aa9K78+k5(@3{G} zVbNwmNwfAgYgSM=~pD}S; zRsQzS!REWWq%>e)hj_EEQ_|kw%gp{?L_Wt zi1mJlj7yi-IPWgFXumfSmsI4r)UmTD)nj|ER#})jZrG~y`K;PV1~Vp46lY+ft-nQ{ zLY~ZQs(bhNP~}eIwWANg9%bK}K1$qa@$xF=cYRh*?Z5wCf|C9#+r5! zZ4Mlg#H{3@*w8AjZd)b-7h*!_JfdtIer0^q0?2Ro9tJKAW( znh1XQ-Cy({>@qUprTS+kbT#CD{kPElrEvdn$=$y>MEf7iDGYPRYwKzSUvA2}p1&?A zux_%827bH3DXqqLG%qHEm6&HwWqJL8CG}e{F~8PeIxn5tS*EjDKfrPGb7~X`6T=$s zliTCU2_ZD(mx$bUy_91j?N}rG z^Kk3VFbfsf(PbYidS`}?$yVq)n#(AmWBQW4NzlLnxshVG#7;FZY8r~IFzdJCE~fi2 zzDp`dLf#X(om@|a8X<);DggaJ>fY{x#}TN;O0=(cugAf=M&DM99uHmLzr!Sd6n5{% z#MajKox9lkvva@fcV)fH(}6u(d-GY^TW1QaXZCu$WEXE4?|?|_-&uNx-GNj8;J8YI zW3V3_^8dv0kH7QR?!U0~HeYm3gJsW&!RK2~cAN}P8$WXU<(pRHg58D(mz7!X&&YUi zhi-A|L5fa+SH(!DbY9B&Qv{D^4@!MjJu5}1aJ5hoFN(@D-fDLIxSPK3caP*8T;|e( z-Pi3uz215K>C*Hg)(+WA9Y~R`6*@MI)%&bDH&6kIZ~Zd7-~3vs#@))_@9zQ$pk}|? zZ%A#w7G{~9RB{*b&0#nvY%zDSLWhS}AV!thaxFDT0I3M+ z*W()c@0i|vs%l?>H=(2WNWi6sAL_e?ZiPiuMECw-6krhaW_M}o*K>#a0<@-IpX)ty zrq>Y>@OaiOv3I|V7OZ?pPr@&QsfPVcsK=-*)N_Ypg%wp(rifY)8AQwJzId6NSV_(e zMz1e_Z{97)tO)AY4K^%K=yS*XkKxzmo@x1 z65;>tjf4Lec;x-<|9IrRcVm9cU$T1R$h_Q(uS!i*4#sq(o$b{k!y>C*My+3hK(wKHIquQ$CO-37Va>05lNZ*)s`PKEl`M=vQz zy{wJe_udBDUR&{nf4uG>07R|JT`VPPckTP=6JKC)+sCS};k?bSlmwhO*u*lZmV3tI z%VNme!Ta61k9qYWbDy+C@=v^Z-G)O(HK4d;r%z?tn8fjd?kiY8gV_eVn8kC)5}f%t z2h*-4q1y(RRTo4QwiV*{M=inM6dT{T`G=?F)p%9gxggd2e_wi5UqeR!fph)e%=bTp zwC!)=WBo74gS{R9BzZfx@=asD!G?+c4_ zGj44ZzrLbfbb)tu`GRLFR$eL5KSGM@Fv{Er^RtS(b33Ljt7zri=9HYCCAkKDn@VqJ z#i_Il;CGI&8!fe`dsgHwMr;(H*-`V!@QDBH;hPO;ZRHonweFS-EC21M@3uL;Ds4Qf zJ+QC<{Kv%8lPjaOQqnS(=Ni2dCpYQ-M*^N73SRmTV(b3j67auq7W;3&X8r~&`-hG5 z!o&0a7~T<3TKaX}v2_bQ0~Rgad1+Er#Ed#Mfr){9cEf1(LZ0hHuCl#b|`i@6=uR{=FyG4;MBRqJrNVSH@{vP|b z@Nhrhzn#l}tIhv7_%Ga!{VlxHFTpq6&aG}7eM$I2oGyEMc=hYm4Lv`PmE4P6xub6T zp$AVM9x6R_hdcE812x-@evw5z=9#%Y;EgS%Cfz!+Vc1SWlN&_Amqtm(g`-3fiV=lH zB83V;r=j%*o)Kc=PKY@CMex=Mr7F7jRELVnANVA{5$Inr;6?|)AlT@lttK^-y3rHI zKUbphsr`hS?o*(xIr%+pqxXoqy{*YWLopi$XYz+%fH3EMPj|NFXUdcJnixaycXTaJ zy5f8~ydAceg4=a&7Uf?$9pHL~=Q2TtA9w?E2EJ(bEK^=R&iHwdY}G-zRD41iDPlkR zdArKStmcfrOY~f*3+AUVFvmHhZ=4|hYW9LQKFNuZtu(YPe@vU{_mUwANZyAOpt%1Q z@g;uzgUm}wm_&9!KiUG|Nf%j(U4UOo=Rw=0Ov(q+XCCx6NN(&qu}+XnBwix+V7lqn>5XX~${DnMDvj@C)e3#@X`wVSK)RUE_ z(A>?U1R1H%*1<9kBu2Kc`~jutV?zq7dCBh z2xNuX-u)2oYWO_?X3?te4tkyZUX-hz9|~I6?BAVT^805lBflmUfPOHw>Ew+dPqXv5 zor$#Wu(z(xvMI^!8z$q)a~ns#K!2rb4tB|2gC2xirf7C=gNkWF&C_gbes&g=UkE7> z_0YnMaX+=rCA|V|UkhNx9uu8Jtz>F9Uvme2K|jNa>V6^b-~i2J+5YBX4`i0!rEIJy zAO+gNTKSDN8B3`rGqVB#t%8i1L}Zq0zH!*Jo7x&gEj&k+&9R{QW8eSIkb=1h25NL8HLp5fGU^kgYOgSWhhq-#Tv%AQ~h%$uF#4l$6 zbjpFJex(Iw2TQ6bN&pNy%EDl3ss>nbqd-(kiCI4V^KYB zlfp+@UdrW=s3>n2RQwko!7bhg3WPY;iI-2MbHbZD!87;Z0FLmp$ZM8l8U6B&7yseF z?!kzdH`NjXO)=x-Q?j!h(aJ;5y|)f(4(XqU5F_7>W(qT)j1e$R&VFOh)L}@w*Sp?s z9mIe7673_ai6 zPe6E;nX^a@4#jh=nGwtzA7$GT3~KEn;o*lxI7jv>{bEyj#J z^T&(=mNRiGVLJ}dH=S@EGc)rWg@o`va%7Yy7r9P(3%3+doSmH#1Rq`P#^h)I&eE&_HRU5IEhFA^cT> z6l}*J)B3!lEDQ<5`8lCIMp)+MhI|T7)$ZXk`;n}g+NHOS(wYvuZown6>-`rz)u)H4+)LqY~^WU@Em?QpdunWQdjV)R|`iJAz+;qBg1z-pEo z2CcwO?o*dSz`lbbSrD<0nJJ8;3`#HhSNM$5$^@+56mzxQiHC?M-Bds`_hABEqa)lJ zWusEo8)H$toWU%`1vCp^Fm5$dfcw-7&d3HAv(&k}={~y0Jx%9oML9Ch^Z*KJ5RR_v z3r|8Y#t^Qgi=JG9qB9My65ud1lh~V-02ZAS5ES)er%!05i1>OpUe@61quLS@K(QGo z%hSL1L!Y6U{seH4I%`jGQs(BsXN`=XBt(G_6pk@tjUov9z^XP0sZ=7GK@+nQcs@|oF?l$>m1U_Y!)b01NjsGC0SN3p3leC) zCbmyEY360KV&V4VFC8;N(i0iCc93^#VRnc@Q0WP}OnO=@IK@ygprD;~u?9#m#?6kb zQ2{YGYvjBui0lrwlFE&mtxcsr3-Kho6a-09XO_dc;47Q-HjVzGljW=SvcR*RGzup^ z42ngx%1)1QPhP&ngo2z?P1DHX94(qzlm!69`W_N7nT}>Z%&k4Evrlo5*j26(0j2;xcf8!M{+ngy>$gZfakCMs-jJ}lQ-Re>%8VHsl zTv1=YT9D3$n1^Tc#`sjbOULiS1-p#0+)jqdB0e zU7=+M_ID@BimY;^y=1r_AqO`gLr)+vcl>0P;GD`f#)!9JFcQd>HbtQ z?;c8qBP=(eS5ap`Jd`YdBe6ddFvPkHm#}wc5@s|_H5|YdDu~s-DjbBKUS$CY zg@BarWH_SzZS)G&2Vf>)`Ds*CqvC)%56D6_0lb9_>b?}2o3uY(3n$NQQFo;l+oc%o zUL3pwV$+AryB}V|99h>~4pHK*=V}U*;nE<=|z)uOz;lOp^d!t$opC#44vb2Njs{^ zV-Qv{W8fp4D9KVAWC|aJyb11~{0UU|eaJ$FXXBc^&ZvhsCv!%|SlP`T$NTT>QY-vL zPY%eU%)q?g$G)=wqRMyG5$G(53v+|>Cfq3$PnJiMz0$QH%SxW`C~qw2b(ELk-!avV zfI-3TGqQAm*N6*E8tNAOwL0{StXE-MoGCQZtasEqz#G#xEpMUbEML#EJq=!c%kozd zFZruQ&`*<2z1=~LzJ2@Dba-@QJFS6MypM4;NcKU-@hoEsDUh}_35`DTt}47icEHYw z72TIgs!(Z}LRIDw%+qGXSwVTlVOC}O+$#{?+*vNR_O|Kj*aKr`0Of+PKI=PvL}sx% zrdYwDv+lLaR8!5uE68$ljmX+zl+lG#l<*bGb5Rr2TZz)L6 zmGN5g8c56RYsXAc4x~4Rs=eM~rw85o@Z4TJZkX#us+7A%Nvw}+2J>e)Z^9rLVGP(;X5HF)IAX>J4^v`sU& zFeI}PAr=y!WnCUr5CF@=av5eA+QbL3gDbN<4R~89u}`hCG*n%n?!ok9v-%E>8 zIDOwcYAn3AD!kp(b_0AZ0~gkJ!T@hi1&w^k_>@_SAdQ4xOV zL=fgP(SG{=^t@el`e!1cN=6-!e^%i1xRjoB19}yYOXmL`+@dz1SjzC56b?096}8B zA8dnS(cE@gHt`hAL32{FI2nf!LS307&et-U!ygA9lF4%vC|4EH48uq!@!lbC5swc( z73J~}LP@qrP|UyqKejYcAR=y~=+r1WFAOt~!t#;hd2Y&&0*8({(*P8_S(z1?=L;&o<{#L*(?+&-6MOV#Z?}IJ%a%fwB*cL|TlRjO$$&D-J~yH&D*! zQ1_%VBW!3F(&$6yTeLDL5x~<#!|e>6WL|ueq<%~uDnvz!HT5h}gQCU!nO2$pngCD zI@~{~+^GCnGuWNOPhw0p_oN-qv-0UqMOZ#@b*Fyz+L|qIEAIETY0GA{Qzb`Oh+I)k zo6P4&R$YIMq5>D;N0_Cy$$i);;G%?8KAq`)_SH}x=rhpSQ}LpOQExc+S;4WM^$RFx zFwzau4(Gzv1L-mnsI^=!tfi}G=-i0U!=Z!o((~|t>e7+6DRRrZKhL7!rXth*k;Dy; z^OrN>X@~#ZGq@+kjOI7A#nC31I7r%GfiZ11VBKMAITEJ!ow#9&sqeS+Mp~p_o|o1x z8oSwl&BD1)WafLh_j$q-CM|u54QsT|HOuIJiG}aCD)Q&YG~5r_O5zKGS4I1Uv@k(+ zGOPJrOV`vR_snJ?eee*=5!#@EU$$4ya^%TR;`F*Q8)x1Kklm3!0rvRoT+kIrL@M+9 zP0HMq$ht4B5O-rqQg3*{;yuBX!&cVFps3#8amm5OyWPF+KIul?X#q;a?%Eg+lKL4^i{_V{H*o%c7+H`lS}zoZ_5jE9Eq zoMCuh_W6Q?`9+~z7)EdoC3KwN%cd;&;8P&U>OxQ_DINFw3xON1&wVa^=*pq&TU%2B zqxSa{1?49B$%CAeT0Zsi@wOVW+yG@=fo<127I6Pfd~iDgu;yEa|VS5TljS zrZDJ2`*k74&1c5xGng(c2?xO^K4BiG=-8$KkU_yIZs0pK@&aFbnL8LzH4`feE1j#8 zN!3@twFN)@_j9c_NY+1*KW;=F9|h%y;o2){H!17W58I0><{o#Wj&nd;q^l1Vg)0Br z&p(jVp4JDv3)dV*#t$YJXCXH&uM}cMtC2la-bt3^wFuE6b6Dz9@Oy7n|0T&y7NS(- z&?_^YMq+1)P6$SR4X=Gra+BXTK=6{>gjf=cKccE9#m+Ho%8`W%e~2kz9=|4Rl349W zJWk!8Vn?Mn$P&UQ3!X0by8$Fd&(ul#&!!2*P8w;XDjPzyG3y_bmixsl9kfr0FOUtL zO#`cBxfbfV6MU~>{d~cLp-c5e=kr_~CbFWckecpvNL54o9PX zS-CSv_G7u8{T};4LfU)|8`TZ=f^(1d3HzK6-=<(%MiqK}La9t`bH?4WW**S&iRV%p z`>{U^&85-?WCS}pjXLuMs+{rySB1y8)9LDI2#h&O!_c5O83z^BrkL~kB!NoM+(Ftv z6=`Y!GRp-bB=J_D0*fZ@U=p4{*V?JE(H%qI8OLi zjGT)N%c?^)3=uaWwJM53v+8ILQi7JGD_=g@j+$W;qhXNWasigyw7;?6)M^8EJw zW=nMf|3!Ob%U=Sc5t4nT20mvaTQ20NA<)siGJ+Xcf(PLm{5dr11knSk`?&t_h9_Z9 z-dId?;q@Urp>d`A6^3~cFST28ibE)x>F-WsMuPDX13}*KzVN99PfLleo{Yj_B@B0m z=vs9?8>8A>CH%PFvgWvZapg+tu-$(A9j-N6W-BCi^&h&^_@yem+iws1iGP9~_??Z3 z>Dd{P2JQCom4bO1g&hCF6`9tX-Tz!_UT*r)b6xY}&ko8>v2seCl6e&b--A1l!t_;q zJDV0R>Z%jbhh%{}7My6@Jn;RdLxV?BMnF#BV`j1_5YstS#b>Gu1>O_JDm57%L4BE@ zlRQAkhAseDf?n+ut>!;aDLx3X9}P-6*F((5nN0wwHQ0TeV6XBWC5h zXP8B-ZvND0{De>*;)>PqB9WnEc)qt|EhV^L>>cN{M zK%8K>S)wa~QCP$49I9qJ{R|Ov$y$=pSh=81G|3ux4N)zik%Gi6I{n5`4t%X){ zmW6t@tGp??xTI5JS<9cUz(df9enS$KtGQ&81L5Q*2l@Bpd7$vRE$&8AytM%ziF_tE+KpKvN4nmTSJMBQZ z_}Z7+jGY|P16yMocs*o7ZOn`fSK{5VX^Nga0NtFkk#$t*6lh$!n<9`y6c4x=-ICM% z-F!xR%V-}jqPW4VwJ_h#>uQm8grCiJz-05jU}UZpF`0lWXWTQ?_tE!qj($CeJPTPH zsV+fX9`mP7QPeG=^$d-JboTv%ulmz82Dm!pE5zw_VN;0WS5wE3vsdPX$~dZYOMTNt#?DHy2e*2jmTy2Y z*mW*nhHH`zTvsGXB0mRT6(Tg1cuwJWo7Mi>6(xE6cI!+y&7Ksrb-vrID>O&G~eP^I&p(s4^qU8|p8G@joth zvC)QT1B5wI7Pd}v76V$OYS1&AMW{Uf!3*i@C8x3o@&M?Fs{kiFzF&??}r@)c#^sH8An?q^Kr>Uc;gfiULE6 zya$$wv@3uS*VC<6QvnWZMU~u|I$~^{#D{Mo-Y&jeg{+*5ZMVJGF8!u5{#?5Zu%Q*> z=VBPV-dHEg(Ddj`_&!ss_1Fv~$(6`3<50Jcu1e-fnr!=&qD{tc?+o1XD$fgahtQ>} zfc|oS7k$g8Q9jBOqv);C$A`hguHU{K!`&mAp!K~x&<`~3{g>lC_(SjyW4(}iWOylb zf$3u(9bAHxPUX3V7JGZ+CvS|yQN`I{`cN&q^o4U9us*pOLA-+pAazX7&1H%`D-C8pqn>!INBzJC8fDegPt? zk-hWGK2#c!+|KY7Z=W-juvWHu`(9=zhxoa;Lsl9)OT1woeD<8F83vn7BC91TJwa#L zj8PfQe^k^XMILN)HWCFt;kb{^b*l3mPYYPtmOia@G z^?b8Qrz@f}MY3Pr_MRMC`3y@=!5IvL%wFBAI6)#Ffg!9CGIiR`_So7i`Mkv z!$V8cOW5e-u&O3y?R7PM%+`kB;!>P!19+KmOKmDxI4QZQKp5Akudvm22UJ-`G}IF? zt7RwKw2S4JW2pP{1um!va#LvFfa7*_JwpIuVk6W&94G zQIKuNAVN+>wh)2W%adyHv2s=8=8q=LN8mp*=v~&>@rO>5Hb`a_N&l#g(gfSiWW2C9Jl%~u!%a8iJPHTx2fK3Rx8xSL zL|tVw-^*Bez4l(CmDwt15yZD&L(^`r5f5?K;V$RYrl>ul>d4+k)Jgs@1VWT{OU`lR zZx5(fgM1`G9QoI>;4Q30UXhZlSIOy?rB*)r*^egu{Oz~eI`OZH1;V9pu2a4FK9aM< zn79a#8v1Ccc9Py;HtDd8c#ITp>Qh-jOz}JYaL~{+>{$i`#&@ujqK%LwlbWQPNhbDW z-QwhSQa<>ak2_#brJi=k4zJCqwUq|2L|b+K%*zbNA)WZXY1%}`URUmxSW$g!twEku zyZq34c!mPxf(Dq8ur`2b+%1x88Kbb5=QfQfpOvr1*(eRA^rbpO3=i;zERV__{?3-2 z(&9L+?{KO;-ln`HDSt*`Xv&PY`WU8^7F{UM$9!|&toRfGS(Ote$zl^KKbHzvk9aFh zeP~@WLY~%r>pDmQYZ*1(#WZ`1sA{}-umT>E4dJmK5zo=O-WHbDx`{A2vnKs$|kC*$m=(kk6Mp*&n=F zw-~Nvyux0ITM+RENqrvI7wK}VpKQ&{($TL57uMQ6r`6i4JJ!GHJy^`Bt9cE7gkDV4 zH`z)&IFTvmCsS> zb%k3MDh;NtLlihA2A`{#y*cphCgAKiOzbh_QArzdOwg+HPoy4Yp@jt8<0_ZNEoP&z zWVx)jGNBh566BS8FzHGaqiBHF-Hq+)2570D6Z6Bszoep&DbB6OxPvN%5k@8O&&9g~ z7%viZxZBX1qOvk28j>DRQK1mn;}(we=xdJYlk(TJdWbgO5)}D!M*f`rxovkIg48BF?C-*8-CL)( z*%WB6ZJ#@df>pnlXxh2gA&%69xxwCOY_Dx)Z*}mrjnQky^2dENm$XOqht z$+`=#!wpgDdBstC&;~W`*^Cc-el-4`7BGk`x(A)6^6Fz|S*xV0hF=S|HHUO(I%yN* zbv+1cYn1ad`XYG8%eQ$Sjy$T^@>73Ye+KJ#`Sz?U7JsV3RWpTGS~zn)&qrVJ{KU#N zxJ`4%#5HV2hU1R|hs)oJa%{Yjv@zng=wqQikt-vQ%t+dIj%-;CD#w6QoZaMEbP%g7 zEU%YkuS}23?p21M9$-$IUGNOK;%V*XsXpro?|DHB@X~51+(Y|NY+cNG5_VtuwXNi? z9!IW?-JJ8j%=D%@s|~&vw}%|7(_rOa0)Nw*W`1~gZ}W(cH>_SwG7|*f^m`7 zK;n)6s7oJmU8J>LGCz9@@u>7Dtd@QLMc+k3gs81X`F>a*d&z02r^#WHS)!K~vv+yp zrX`$AE%~Qpqe8?X8#C!q{=0?`=II-u+Q==LOOwy^*Ja`g0G0-8~IcZ5C&(%ki1{6z0%k z?@HZVoexVNW;Gct0IbC&L+?@1lB~0??5Y>8FGlOC#Kyzx+u_5aE@){nrzotG4fxKu zZ@nqx9Lqto`Yd6)TDI{Yy-#e4IjdfWTBb}BkJgBwd5rpWZPz2zJn?D4mJ$C1r`4`6 zqZWDT~E3xcR=qQJp&Y2Q(poIjeuM7Cac1`2F`{T zJ)PY8NE(kP2^~dL>?IlXJAxm+8g&81YksNO1c1&M9`G5)!<4zJWK3vo@$CjO>L z*T&C3P7ybNva8R9S`5?OY$P2!60|na;;qp} z4}Aal!=-LHo$5%*%Co&(<@~I?0m+OE&A~YTN?h5>0 z)Uqn+rHyH|j@LTXqMT=0mSS&qe3R^AcE}T_Z)EL|BgH!p8Af?YhuXC`rd)*_u2MFn zvn6`Er%A_xm&po*?cd1eSsbA!7|rVY!P{niBI;y&;Io80O^1QDU*#w!@;Us? z;>$0gqlNV$4*>3xshc26s5U5KW=lrraA#LJHi5VkSW6k7rp=tGY_Ye~fI$%9ei%olf1;*vA zSry7qw3JXQHL8mrcID?jD0flqj5^AmY`VND;Y9UC=`sG*;ld1t061cOSn_~m$k;5S zkNu2UNP)K#a^S_RY~_t%g9E%-znC?_>%EuRUR6|$y?&L1g zH44q3x62>->i~~vrfzY_32xqCrnIZc^QSDyGBNAeI7|u?t`VGhMsW^aX#ap;!9Lx# zb7Wy%f~8kjxY9THouWDwSH2H)ZoOe~i~3S`&?Tzldj~E59Cyg&s*ra@PAj+t_PAcc zSQ)Fh7fL^gOU|G1X1o99lOK>6O+Ln9$@I5!i|A_u{zYM_-O{taLkle)keT3sb630M z{72HYf{a%>Zy2ivUl?1$M05QJ47WYxnfXym6A)bEU4olV<|e&z+5}jPvPXNvNYS)< zz&?C(VZ2MC@#^sLHjV1UTOZZZ!)ufW`e3e~{aH-h>aVXZlJ<|0V?=qZRGy!(W;M7n zvcCy>70Iq*IVvwRc^~$XE=Hl>*^@4CvF4+Y3>MQGU*!DZ$%pQQ>ilthuM#z_X@~aP zn6W)vCN%=E)V0+Z3wVKCDt;<85>Oh*`JU~_t~pRrv5zDSNg8q1*i(4D!qThi0PQJZ zQT8&K*O*;>`YT)^qn>V#3-li>6vJcl26}^ zaUsFQ8oTnIF--qU$s-*V7p5HN@wSzxYs!hZY76O+ahRfftZFoK0qYpe*)naefqgoS z&Y2Lz3^^JyGDMp&XZ=j`J{Q0J?34EyObEwaBDj7)f-evOyBK<9tRKA+$7#97e<9#; zZgEaFR~*`vPi1j7IZ8WrwWeZn6)T!*8h}9KyuCm}HGZ5db^UQN-g*<*g9IOkhmX%7 z{4Ea`<}vIDyZGCxlxuA*hYPcnLG7){R&lP+2qq*u!p^b7uR{8Zwp594VowK^vjADo^otc1fH0DF}!F z)WarG7iH(c{t)DHlrNfBVHsz!yX3aJFd?a=@@9r7Gx`F-^kT77D$SEAy0jkuR#$Ls zcz2odOPMSoD>=W*VBDaqy+ttO-O=+MyLoVveSUadp65BGtIM;{M^o=M!{hBw{8nzU zxP>Pj=>ljXqiwUVDg*67E>#T#Kzh3e*QzYaYmmYRonWZ?mC4b&MqP@YLUaT zinR)DdpnRhj5O6ZU*hLu#xPWGs@6uOg)fcuAEs9K-3>C6O)0X zq{ZYT*u_-Jps`}(RPvGES`N?#!CFx^fJw+MdOdhv&v#2zDJ+fm8;}`gi1+zb-)FO)-f@hTaB!dIu9AlDE)%Xz{eO6{~TK>FlyyU z)GSB4Vl8OZFs~^;nOh1}AB0H2JUyA?cm%StdFN>&FzT@Q^3dfe+$y|Pik(GlINxVB zVNXJB!Yj%>;B}V-S0{-JH4n{bKNJzUQn8WGp47>A`=+&7?L$o%9BqnzR}}T5nXocd z(9|pnWzJn#&-n0>{LTyKGM6>gZ7^p!o;>gu*?Y*@lMSd)>DQsr0>)!K4@vka!Upl_ z7xbQb(cmI)Lpk~dL*Ty?=Ky^;H^#p31CuiNri@z}fqU_<_`)~yuRXJNO?pI)P`IXix z*hdPV3C;hclY^(JODo|s@l7Q+NUJqv!Y2F5+4yaGIsk2(6(a#q&`)wut5fg*`6738 zmn1AzUXf)m08Ql)@^Sy#wwP72@*wjy8hr!!mP_D2uC%^2$|w`$&Fx88Y_gT$Gv%D3 z+8*UIR5#?sj_XN`2PbCPjjsy7G{k5sNC$3_pN?~S!#F}-OkB)M($fJwLxEbRUt;o| zHyJJVSYnq1VXJ|ZQJN_;_Mn?y?9Fxb}TX}7k^lUe z7n+Ot$unXsRas1QjwJeH-@?u|6(#^T;S+ihq-U|w(b8Fx7HI`Wb9)DPE9 zbbU;o@ZfnTK5rlPQ~DzGqMU4`8}iav`#ao|oO6z`zvlUe2f8{Ycxlhr{!ma#;6XyjpsmxS zJJ?K^c|&6mok)2$MZ&C0GSO1Z0@_2FQXe_XumGYxx!H_ zo*LK|5>JTBvV-j!(;pUXv%Y7XH-Ju5Y~t^ZgexBLAC4`tDj#(_ntp`3()Q{A+zG!$ zx++4w)!md~ITtDXI@XtnSp?WDs)XHb4|tZa+u;W6@z7=Iu@q+`aBAe+v4!tIdI|EQ zKV7e5wUL`tJEgE$m!-HQQ=OX4(LO?_9?VAYF_;HCQgHHF=O0OSEJrY0?u;{-)2M2k zjZehaeL3gW)82jzdG`_cz6E}<8y}A~#CaR>5{Wt{uer~9B`07rA(_Ju@bT93Sx3~G zelE|-W|eeZhV$!k9kniA8*;>Fm)wP2&X;^POEkU~2}gQNtv|;4(O)b#t6Uwg3WLrfwoG3tREXQ}#vyqPC%TuoNcz(M2 zSOd@^r5$3`i_RYGUD8Tp9f&q@r}zFaOS+sNLP@l7mf75sPa zW88p=x}L=NUR4_(kH4e>nR9A)F=)f1?L};KlAdP$pmz_@%4XOC^*ya*w>aE0nRmr4 zU(>1*3>dl;@v$8`tg*n_kb2dz(5qt6*y}vk`HVmHR;pHsYyaZ*oq5iOfmq!a_CZMZ zs))Fn(BB?vioT7NqU1ofT9wzZqag8<(|qt6n=6$7tArE& z4(r?>{44XKE{`?E_?Ey=Nm|+R?e^O8_Tl#Sj30U?CP2{6@Y6rujU(R<278B;bA7)S z`+7?DdWsr6B{a{`(Yp8!9HM*fJ>`AZ0XBq*lUV11U zW>Cf&)@YUs-+!9AO*+!F{=WNG>x-2?$Nspturs6tcfrPFc%VuDzUx-9ZuSn{ui9bNzuq$KxxwF3s-fd&L zJ}OQ4DpoRLa7{NxZFXbR zdmI45rNP*{1r!Z%szwGZ+`)f%Jji2q5$!6$;xWJ;MinX2c-vgct9~cZ+upd5-lw z!!E;}NXnjsazLaA9Il2tsS~UI%=0f#x~X?;TwlCnY9)=gODDbf1B$LJEz(J@a9J~b@1x#BPr zQDFP?@G|SSkBu%guGlrYWcWS!^ zcNI+dvE4k6b&Tf%0-}#7*zDEqgwZWqS^g5_+u|Sx!;m7#z+8WD5v8vZJTNEeorCsm zzzlA{_H58hAkUpEF_*|w!)lukU5zvl@6t#Z;36iFgel{M1ZCNSjw^_4!VqJ4O$6!Z z_~k%Og%M|EN61qDn+oL;#yz+h&J?8Yu*cSE{0@@WED4!TY#M|sxNM^yfRo9yq^}|y z-Qw?qtB_kee2)Ww%F{5?q*lJ;H6MRcYGf9lGllA#)>z6oOq}(F>_pt7ngrp*b=_-N zBc%2meAbC#wqCs&X&V(jo`-~9x ze>HX`P)%Ojx)m+cm53;@0!aapT9i0cArO+PU{sXUp`tR0inK)tCLlv5nF<;pAjOD6 zqEIeQ5CtiNAz@GkK?0&-3;`1elO$otKmxo!-uK?R_TBs5TkQ3(khT7uv;TAcefD&| z?Q1hvC%aLfR!>h+J}M42Hc1^k=8T@(A_hHECJvxyPBu$yq^69OOq=R;b2DStTJ{#E zAH;G3J~;tB2}vcUj6Am+;@n@`+c*y40mDW)-7O*RMTJs$RnOd+x7k6g73}-rcT~8} zx&-BGSPmQrev9X2gn30yUi9T(zeA7Rp9c_fCausLmia^b?QnmpCz z7Ih%+OUZs^Yvb%F|1X4>^O-_^^bnsk!N(tuHq*U@^M=`zz9V?w-ry6xXfV_ke_b^- zmLR{Ifk#73$r5Xvw*jPw-&h}54NyX zk@kQGZQEt1+$b>AO{D-f3@q8*_CY)oXLHUPAUe$+i)o=H*l*%)%rkV^D&PT5P}zF| zeb$=s8|k=(2yH>(`b^GNw9ES!EV0EF{Tpe6ce)49GN01)jcFxBjFcy>E)wzTn<`!D z_dt!k%N@lseP#8|2#Ps|xQd{nTl?Z%lZ=Xyp83x%Y5a5XOd`67hb#tb=JFrM1(2 z4`@${1IoHrGY#+tgu<3B9`r)w=4ub^39<77&Z?%h=^2*r@pRqGXo4{Se4xa|5MZYV zM6=VUx_y%i5_DMzyDWH~z}08cu!4Umw0ODaEX5LWS_#bA@J=75q*ycldNE>z97%D^ zoOX6uKmJF9!(r};DoA%plV;I87A=qHEV_ox?6Px<3SB;<+XZ_AUNhh`zsa)i(Y@1u z(!G6hzmI$%6g&5=*yv$XfYj}@XVaLI_wT;h!0IR02hjy>Vyg*tbK9jYq2ig~= zr$Oj8$PAyfbo!@NT6Z3l5`H*3Gfi5jqxmpkDt$!g{{4v4h#o#$!H@YGh?F~AP6LQ= zxwy5}H(v^OS+4@pnZuA_&$*aW$fog%4%d8<4QvI`c+d~|T=lg1gf-I>e?q5X#*V)O z?P+w2U;ZF+J8RXO2RPu(>*JxOl)Qm$+7$7Z1)86heo*b#dE#RIOzH}*7h?dcRy z99OHu;0+}682nJp7&z%}ut^=-!)r_RR$Q6+oH`+4{zB*cqh;gtBv>NbKipp-tWWqE zxv^MzkmE{@z&(x(6J?zTT!zy?2x30^f*Y?9G1a)hPE$fdL023{Ze1n`FOs{f2Jdit z7nH!Y#n(p5-bwa9L-+0#p3F&p@I>=qpyJfD4^T)F;d_Ug^>~HUi?N;Tp>h|vwmh?? zLl}vlWn@Ax(%>iVT`P9)ejQ_?XRLbI|(>v=mdBVvfDzr>%S>SIAWoD2REYn zV#%+4^U0S84+^|dX8y+{uEgup4X+U|$=oC^-PcvO!%qbN$~!u2)~RP4aCGpK9W%dc zet7V!k4`Bj8=$a-1H(r0tSP?fcy`u_qFD=aBKw5myGN{chBIJ@vVKEdgu-;}?c2_!hwPL1jc z#r`Rtoo>I2lqPz=BUDDFf8fV5UNb|&&n;}z*3^Hkvp&~?PW94H_0s%aP>SCW_@j7* zqiy6cTH&NLUSwg7}VH6%VTp?R<5(+;CM|uR(P8xsALTS6wT?6W(b*{ zf9M+k{=oy}l@PZdRBPv1J`|x3s)@|Cb323kskUNB@ZRaX-~ejSS5X5b+kke$dDMQabu;;Z;qp~;C_oo6zBqi zev=L0vB1Oqi4}lT;ZSPtzVuf7-9*iuRE-@bll?1h;~y2RkS*=W;q_@(_-UG#X^W8> z{Pl<+i^=cm9VRq`E-+k()i8fPCWPbYEDPfnYU~%*jeu7Q-%2bfhC~~+@0exMtd#dr z+->Onpx7r32nu96|?;0ovEwhUcdB$ByPV_woYzg1pgwc36*xR(WLmV-3 zJImM!<9GwfYF8b>CG;#GVA>H;(dJe7F=9*5OFGpA!lsiDV!a)JopxSm4OlRi_28XF zWff_b`Xq|ZBE8qNwj`4sHOC$kL~`wv$)t`c{Xk6hOxR`sKE>3CcQF0Kw7P>(DTfOy zVf`vzsFs9%f*GvD_WVxgU41J}hIMMI=CxIk+zpn|6@XYT=4YowKGk3mY_AgwxyJ=P zjLlYmX!}^r0;(fmEW)=$SNH?*ubw~wl}kM$9^z5mC;p+9){N=skur1F%|M$A%^aEFjtQxY=WG(JQh8E}{7qpHZTV7Kz*M&#sql3~(gGB$jz_Q!p2!0FNSY*(An9 ze5nPL1NTE5)h3A*@fR}DH?T&zbj#z>Yt%UnW~Jz4x67)|mR)%aVgm zMYeQr2CvY9b(Zw-?at6+=&KT)%K!_6nrl-QdZ(yemfb3H>a;i08+R>RndgT5FX2{k zy~5U1^}HaNN1MaBZ*!vfQSQt_LP>>xbp>*}=I>FGH419RiSM-E5KM-Rfx3aBhN7(H z^yOP5W_47V#4Hyq+@>=^7l;QaNO}LSoTCgw!rni+8BZvO(hcDSxuEfu&K4|@v!`)n zPkU3ZYF1E56Ws0|RQV%fxrGSs*ooSTm))VQIIaKJl=ByMEwwbWRw=t_$C$(ALs2Lz~#Fz_|Xx#XHlB1vaLAs@ao;(2bf7rBP&ESOap3hcm+hpjF5@G-^_sL zHSQJ7xsHD5f6U)6jyE9G8tbX+{F+-@=W$ReD}z)4q!a@0pSZYO7M-ZNoWwf9$6Nwg zOZM~1c1WKCB(JI5(O%u#02q@19%5?O<8U>{MaWxEl=3p z8t;1V1kKH|S8UGnq>50tYOLbW#HUw$!3(+LVcJ0+^K z&llR`Q)Kc%Qb5F+g~57osW4-}15lwtxv&(Kn7fsDVu&ryD5bBic@Z%p_xuGm`ZI#~ zTx#=4lZsvid$xn%>?!Nt85d{ZK0b8L6~BSsSMgPpv5N1dq1qa$$mc^^1%p*@f0+9; zGVBAN+)^C*JrHaT7HX@;_DJRU5|Uq=-ZhF7^8L% zj_r~I&RX^&YqFUG#hth<>`hK`^k82V3GLMps{YQ<%&hvQRdS~{_aPog`VjY>nbx-L z<7<&2Dj&YZWgXFW&S%2YRWn*DSuM0C6!n@NhS+i;>})ZSRf0f4G`ZN99k@X4y} zE`fQXVoCk$dDc$OIxDx>^^3?BNqllJLm#}EPMwYSN+NAuAWqTe3g%A`|K;Ii5YPD^T?V z9D)~M+B(=lQx+#(ST&K59P}LdlD|F306_&*5w{Fu#mR7qA6)5nE8_Rg7b1E*5X(U) zzw-AX{uv7^wMQ%YVrgI~k`hY(`D4{(9mNC!pF%vG%HLyR zMauIc_!azE_f5yB5J8A4hRBtr51sp_)=<^lt%Y-X`Of}gFH~RabP(TLHTwYWL ziX*wG!hkSGRIr!6Bngc!O2OR7g(}2jv--TVYIcd1^yCS}UsH1$%Kw@K3F$5`m!_4g zQ_Cfd<-Cu5@Eqy=PnI-C}T2ZIS1~ z(3E(jXH2DhxD@5iq|C(Dl{GuX!Q8R=8p519mANJt4EV-TxFbFIZ~WH7sSe#0+0srO zY38y~^m3lGPubK5uODyyzPqHj;_M5k00sN7yWZ1C*;_eS#5u@%@r=9)i5_?DqSclC z%n-!-_f(_Wv^I>W5I?UVNC7UdmCag0Ah#)2bSCepXajhEdJ7u>^MO%gS1s7?X&AD{Sty+5~5kT&9kqvGtbu?Ao?***Ab6zdx84 z$B5cV=FSW}-P!&OS0hBhjJ2|?Y&E}JNGn$_m#Y)Y=XZqgj)o`;Lo^CL!?ufNSiu{B zF<&(3a#y z+(`g{Zk5_D3+nD=TfI(EHd(ccp3!nIC+-^GW2*=zrL$tgT=Zp)O%JN4hBCXgq-1>*M1)@Fw3 zycyuvkXrQDvesMK64w46vMbu(t7JFZfbFsr@&3(FjcmNWYqMs5nYR(d3b6Xv^)>bs zFa5@15V6F&(e>1H*Qt{I4$E#16pRC}+d<^$AS?=zJS-8|mk5@h{DX=ye*M@i`7nDxkJpSmz~#&ZLNczX?lDH|4sJZY|Jc+M7&zj<+H1nNLK-DlYyaJc*mS zjcK-80lk!NZ@(e$DLD*(dl$&@h~vtI7=`L7g=a&h3B3n~z zmPbyeY2N9HE;vjgPVr$+QbB(&8FTZiGGXQ3eO8Lag<6t6k|RG7 zE`5}4o4v6xH`G+^?KF+r!ygVME;&u>eUR2y#3PTyGKwi=oa`h$YuEjl*vTcB3S50Y z!tqY9lGhUEX05+=fQqmYM72z=oxZ6(Hwx75W+Wdj*~2kYoM~|S;FE~5R<7gP;0>d< z6CRCMG<()D4L~%nRYYUd_2~q0$TNCg%t|eI*#Hh&h>@o~^-gqqD|Q)j+9bT=dwlnA zoyv(pD?ztMw{bu-!$e!31ie}d!Z2oisKl*3YS_QX`TcFRE^!F}<2lGhP6VjW? zYe;|fh8k3N+oWs?$n2SY$ve-y-!3}0ybqY&FLo4VVfS%ts#Z^+H;amTHK`_4OON(s zT)hb(<^~LMJBB#Klh!Uy#0k6fgu5Jgkq(pjr=eJDI8Xo2tAitb zk+_f5f9^q9LYQ!m*$5dd0;aMCM@z2DS-I~33GeOtN|`T>XW65()$Y~TUsqyqm7V1o zFXCaVmQUOCPJg)f#Fx%6|EU<`wQZx>hjRtMuY}q~91#P(6!)0+pIhtCp;w zu*gJ91KjcE#rImvJljfhvfuDvOJ0}q>jAKUHP|m^UJK6Tew!COx@P8DN3fw&DtC>a zIUzz;SVl>!(HFlOmu#~4iq*7KkLpp#`uszXXa?eG`f#oD%Q&@0Z*tzG&EE1h=5$q9 zF#^LH(0_v9m1tuwd1XS+2a%&c4kVYWXUZo>`WE8KDdhI+GI{|sOAwon+(V)-aEgW9 z#aAYZCx5{u1JX^}S+yuSit-oZdMEE`9+c9)uW}F{ddgZSyftUki`@rtBb2mw$2nz! z;ET7kot^tjqt?|JId5055r+?ZzvOvZ9Yd{U0Gla$nu6{gp|D3jc93!gf) ziEiQx9NAGF3AQ>A3~`!uA^QRGUGCVAbG~v(ewPh21-=rsz6!Sh&*|!=v?5n3o&j}Z zO*UX^jeac3UIc2XyBQXU6XIXaicKl$=(b?XN*Y2&!k=60ikr98UqbBS8KYa_$Cy!b z_);HCM@0(rsL;-;sc-&rXztEgL$)p1d^9?GoqCH?yMB)|I#3jA#WD50tqqLWDkh9G z1bxbkMZBu2!wg7Y(~SK#Y05f^Yhj0l5D3FEFt1gN&Y&8xfP!hLKZJ22#zvpmyDgLD-VwOv|4@p|cM zRl_)lEyDpvy9=wUW)Kz(|LPCflf8XPLEk+7<3LWiq^3Ml2GkI9e#5@gls$vMT+HQ| zVNCOm<6O=>Pf(3IXu^p2(nmDuBkIiYY!0+WtrsF<2@BhXJF-bZ*}R~S)mEI9HQTFb z;@;uli=c-KFQ(M@(+zkuo!K{VVmvD$WPGr7a8<;t_+nWwJ;rLlTl&+&)_Fa_FvRao zh5XDjC$Jx)N;b9LdfAI{KSvd&^n1cYr&%;k0XSQv6X3&-?wOPqU>E@0?P2w-|t#1ZdXU$N&dj}UVP1g z&m^@x1WI0D>xssBMw+8QT8U-6Y^V`Fw4I>=PJzjH>U0&jvv7!F4~^HZlZO;b1*)e# z6+fi2=R8{uS%ZD>m}pbk)`+to0%Jwc&K-{eJ}x49Jqz)-ef3rQ<{u7sCGARU_qz@O z>M0P`KoySws8N-8wS$0cU%zQL-IgIuKaGlc>2kBpo|>lKr_axWRoVPgBk+Q7uLLm#H*5lB$i?39CL-h6dNLWFxA_CzA_i zUokyR3Z|8^q*&Ak@9ssn4NlNS|5S3icsYaTISm4E>vNuBZdjVBUv1_+L$}#NZrf>P`e{kTxjcEChNd^`r$0?F8Bkhopcuknf;X>u+)^@r1aT)@l8_y9CmWTV z9a-kUMTPK4A=2Qp{WSKT_CK!*+I?4W5BgK5x#F+K{b?r07p|vClg6qOX}`ft1D?-w zI8%8IC=3G(_h|d5ikqLUx6dyqO3$$U8z68zjGF7HxShhiOoWsPlNYC_MA=>iJKBGb z$;)+7?SKDx?#u#Q?(-1Tj0Da)IJD;XieV#z5yuP)mtI? ziR*k8z0M6rK=i^qN?cw16^Ak}YJ+>TBy7Tfxie|)C5o9UPje4o9ud}A) z3ny!(Z{p+$veY|DlPRUE^E;16=|B$jNh|t zcP~|SKJLKYXLS>vUt!XGz!sDj)yHn!1 zMr=CVop@O~3qiUgqNFbCIX}A8*C|hjZ^eP-!}(2b7qcyoEB0pbr~!wGKg=AF0V3C{ zRa^iB_sbRA^gi#kqU)KW>!0$?Kky(wU9WwDy#2pAul$vpu1kro(}=F?$fnP@tpmRO zXF%eA*UROv9CeL1bPY6g%``qAUl#D~|1!Q$aF>7OrfYeiYg?deVet9*f>-^Yhd`}EqyG9c7Q0seL;C+08y$|b4wqSn(Y)a^HbH<;|AWQ9c@{c2V;wxP4yyQb zE`NZ&_~PGz8~;@#9SW-sRaJ+u`Z2j?O( AqW}N^ literal 0 HcmV?d00001 diff --git a/src/radiometry/corrector.py b/src/radiometry/corrector.py new file mode 100644 index 0000000..472dc87 --- /dev/null +++ b/src/radiometry/corrector.py @@ -0,0 +1,186 @@ +"""Radiometric correction pipeline: raw Y16 -> true object temperature (°C). + +Reimplements the InfiRay libadvirtemp math in vectorised numpy. The full firmware +pipeline is: + + raw_u16 --[NUC + kt/bt]--> apparent temp --[ems LUT]--> effective emissivity + --[temp_correct: emissivity / reflected-temp / transmission]--> object °C + --[distance/atmosphere]--> reported °C + +The emissivity/reflected-temp/transmission inversion (`emissivity_correct`, the +firmware's `temp_correct`) and its exact inverse (`reverse_temp_correct`) are fully +recovered and validated here by round-trip. The NUC + kt/bt stage that turns raw +counts into the apparent temperature needs the camera's per-unit factory +parameters (read over the vendor protocol, see `vendor.py`/`profile.py`); until a +profile supplies them, `raw_to_apparent` falls back to the legacy linear decode and +the corrector still applies the emissivity correction the legacy path never did. +""" + +import numpy as np + +from . import tables as _tables + +KELVIN = 273.15 +EMS_Q14 = 16384 # emissivity fixed-point scale (2^14), per the firmware +_TAU_FLOOR = 1e-4 # firmware guards tau < 1e-4 / emissivity == 0 + + +def emissivity_correct(t_meas, emissivity, t_refl=20.0, tau=1.0): + """True object temperature from an apparent/measured temperature. + + Inverts the graybody radiance-mixing model (firmware `temp_correct`): the + measured temperature carries the object's own emission plus reflected ambient, + attenuated by atmospheric transmission ``tau``:: + + Tmeas^4 = (1 - tau*e) * Trefl^4 + tau*e * Tobj^4 (T in K, W proportional T^4) + => Tobj = ([Tmeas^4 - (1 - tau*e) * Trefl^4] / (tau*e)) ** 1/4 + + ``emissivity`` is quantised to the camera's Q14 fixed point (e*16384) to match + the firmware. ``t_meas`` may be a scalar or a numpy array; ``t_refl`` is the + reflected apparent temperature (°C, ~ambient). Returns °C. + """ + e_int = np.rint(np.clip(emissivity, _TAU_FLOOR, 1.0) * EMS_Q14) + f = tau * e_int / EMS_Q14 # tau * emissivity + meas4 = (np.asarray(t_meas, dtype=np.float64) + KELVIN) ** 4 + refl4 = (t_refl + KELVIN) ** 4 + obj4 = np.maximum(meas4 - (1.0 - f) * refl4, 0.0) / f + return np.sqrt(np.sqrt(obj4)) - KELVIN + + +def reverse_temp_correct(t_obj, emissivity, t_refl=20.0, tau=1.0): + """Apparent temperature a body of true temperature ``t_obj`` would read. + + The exact inverse of `emissivity_correct` (firmware `reverse_temp_correct`); + kept so the forward correction can be validated by round-trip without the lib. + """ + e_int = np.rint(np.clip(emissivity, _TAU_FLOOR, 1.0) * EMS_Q14) + f = tau * e_int / EMS_Q14 + obj4 = (np.asarray(t_obj, dtype=np.float64) + KELVIN) ** 4 + refl4 = (t_refl + KELVIN) ** 4 + meas4 = (1.0 - f) * refl4 + f * obj4 + return np.sqrt(np.sqrt(meas4)) - KELVIN + + +def compatible_emissivity(target_temp_C, org_ems, tables, version=2): + """Effective in-band emissivity for an object at ``target_temp_C``. + + The user-set broadband emissivity differs from what the sensor sees in its + spectral band, and the gap varies with object temperature; the firmware looks + the corrected value up in a 2-D table over (target temperature, set emissivity) + -- `read_compatible_ems_version1`. We reproduce it as a clamped bilinear + interpolation of the extracted `ems_correct_table_v{version}` (axes in Kelvin / + emissivity). ``target_temp_C`` may be a numpy array (per-pixel); ``org_ems`` is + the scalar scene emissivity. + """ + temp_axis = tables[f"target_temp_list_of_ems_table_v{version}"] # Kelvin + ems_axis = tables[f"org_ems_list_of_ems_table_v{version}"] + table = tables[f"ems_correct_table_v{version}"] # (n_temp, n_ems) + # Interpolate the columns at the scalar set-emissivity -> a curve over temp, + # then interpolate that curve per pixel. np.interp clamps outside the axes, + # matching the firmware's range guard (it errors out of range; we hold the edge). + j = np.interp(np.clip(org_ems, ems_axis[0], ems_axis[-1]), + ems_axis, np.arange(ems_axis.size)) + j0 = int(np.floor(j)); j1 = min(j0 + 1, ems_axis.size - 1); fj = j - j0 + col = table[:, j0] * (1.0 - fj) + table[:, j1] * fj + tk = np.asarray(target_temp_C, dtype=np.float64) + KELVIN + return np.interp(tk, temp_axis, col) + + +def vapor_pressure(rh, temp_C): + """Water-vapour partial pressure (Pa) via the firmware's Magnus form. + + `calculate_vapor_pressure`: 611.2 * exp(17.67 * T / (T + 243.5)) * RH, with the + firmware's [800, 3000] Pa clamp. Feeds the atmospheric-transmission chain; at + indoor range that chain leaves tau ~ 1, so this is provided for completeness and + the pipeline defaults to tau = 1. + """ + es = 611.2 * np.exp(17.67 * temp_C / (temp_C + 243.5)) * rh + return np.clip(es, 800.0, 3000.0) + + +class RadiometricCorrector: + """Apply the radiometric correction to a frame given scene + unit parameters. + + Scene parameters (``emissivity``, ``t_refl``, ``tau``) default to sensible + indoor values; per-unit NUC/kt/bt come from a `profile` once provisioned. + """ + + def __init__(self, emissivity=0.95, t_refl=20.0, tau=1.0, profile=None, + temp_scale=16.0, tables=None, use_ems_lut=False, ems_version=2): + self.emissivity = emissivity + self.t_refl = t_refl + self.tau = tau + self.profile = profile + self.temp_scale = temp_scale + self.use_ems_lut = use_ems_lut + self.ems_version = ems_version + self._tables = tables + if use_ems_lut and self._tables is None: + self._tables = _tables.load() + + def raw_to_apparent(self, raw_u16): + """Raw Y16 counts -> apparent temperature (°C). + + With a provisioned profile this will apply the per-unit NUC + kt/bt kernel + (firmware `temp_measure_with_NUC_value`); without one it is the legacy + linear decode, which is spatially correct but carries the uncalibrated + absolute offset. + """ + if self.profile is not None and self.profile.has_nuc(): + raise NotImplementedError("NUC/kt/bt kernel lands with vendor reads (Phase 3)") + return np.asarray(raw_u16, dtype=np.float64) / self.temp_scale - KELVIN + + def apparent_to_object(self, t_apparent): + """Apply the emissivity / reflected-temp / transmission correction. + + With ``use_ems_lut`` the per-pixel effective emissivity is looked up from + the temperature-dependent table; otherwise the scalar scene emissivity is + used directly (the rigorously round-trip-validated path). + """ + emis = self.emissivity + if self.use_ems_lut: + emis = compatible_emissivity(t_apparent, self.emissivity, + self._tables, self.ems_version) + return emissivity_correct(t_apparent, emis, self.t_refl, self.tau) + + def apply(self, raw_u16): + """Full raw -> corrected object °C for a frame (or scalar).""" + return self.apparent_to_object(self.raw_to_apparent(raw_u16)) + + +def _self_test(): + """Round-trip: emissivity_correct must invert reverse_temp_correct exactly.""" + rng = np.random.default_rng(0) + t_obj = rng.uniform(-20, 120, 100000) + for emis in (0.5, 0.8, 0.95, 0.98, 1.0): + for t_refl in (-10.0, 20.0, 35.0): + for tau in (0.7, 0.9, 1.0): + meas = reverse_temp_correct(t_obj, emis, t_refl, tau) + back = emissivity_correct(meas, emis, t_refl, tau) + err = np.nanmax(np.abs(back - t_obj)) + assert err < 1e-6, f"round-trip {emis=} {t_refl=} {tau=}: {err}" + # Physical sanity: a perfect emitter (e=1) needs no correction. With e<1 the + # camera blends in reflected ambient, so a body HOTTER than its surroundings is + # under-read (correct up) and a body COOLER than its surroundings is over-read + # (correct down). + assert abs(emissivity_correct(30.0, 1.0, 20.0, 1.0) - 30.0) < 1e-9 + assert emissivity_correct(30.0, 0.9, 25.0, 1.0) > 30.0 # warm body, cool room + assert emissivity_correct(20.0, 0.9, 25.0, 1.0) < 20.0 # cool body, warm room + + # Effective-emissivity LUT: in-range, smooth, and a near-identity for skin. + tbl = _tables.load() + for ver in (1, 2): + e = compatible_emissivity(np.array([20.0, 34.0, 60.0]), 0.95, tbl, ver) + assert np.all((0.4 < e) & (e <= 1.0)), f"ems LUT v{ver} out of range: {e}" + # The LUT-enabled corrector must still behave near the scalar path for skin. + c = RadiometricCorrector(emissivity=0.95, t_refl=22.0, use_ems_lut=True) + assert abs(float(c.apparent_to_object(34.0)) - 34.0) < 3.0 + + # Magnus vapour pressure: monotone in RH, within the firmware clamp. + assert vapor_pressure(0.5, 20.0) < vapor_pressure(0.9, 20.0) + assert 800.0 <= vapor_pressure(0.5, 20.0) <= 3000.0 + print("corrector self-test OK: round-trip < 1e-6 °C; ems LUT + Magnus sane") + + +if __name__ == "__main__": + _self_test() diff --git a/src/radiometry/profile.py b/src/radiometry/profile.py new file mode 100644 index 0000000..ea39206 --- /dev/null +++ b/src/radiometry/profile.py @@ -0,0 +1,100 @@ +"""Per-unit radiometric calibration profile (keyed by USB serial). + +A profile is the artifact that makes deploying to many homes a non-event: each +camera is provisioned once (its factory NUC / kt / bt read off the device, see +`provision.py`) into a `.json`, and any host that opens that camera loads +the matching profile and gets corrected °C with no manual per-site calibration. + +Scene parameters (emissivity, reflected temperature, transmission) live in the +profile too so a deployment's defaults travel with it. Until a camera is +provisioned the per-unit fields are absent and the corrector falls back to the +legacy linear decode -- still applying the emissivity correction on top. +""" + +import json +import os +from pathlib import Path + +import numpy as np + +from .corrector import RadiometricCorrector + +DEFAULT_DIR = Path(os.environ.get( + "PYTHERMALCAM_PROFILES", + Path.home() / ".config" / "pythermalcam" / "profiles")) + +# Per-unit factory fields, stored as plain lists in JSON, numpy arrays in memory. +_ARRAY_FIELDS = ("nuc", "kt", "bt") + + +class Profile: + """Calibration state for one physical camera. + + Per-unit (from the device): ``nuc`` (NUC table), ``kt``, ``bt``, ``gain_mode``. + Scene defaults: ``emissivity``, ``t_refl`` (°C), ``tau``, ``use_ems_lut``, + ``ems_version``. ``serial`` keys the file; ``temp_scale`` is the raw-units/Kelvin + divisor (16 for the TC002C Duo). + """ + + def __init__(self, serial, model="TC002C-Duo", temp_scale=16.0, + emissivity=0.95, t_refl=20.0, tau=1.0, + use_ems_lut=False, ems_version=2, + nuc=None, kt=None, bt=None, gain_mode=None): + self.serial = serial + self.model = model + self.temp_scale = temp_scale + self.emissivity = emissivity + self.t_refl = t_refl + self.tau = tau + self.use_ems_lut = use_ems_lut + self.ems_version = ems_version + self.nuc = None if nuc is None else np.asarray(nuc, dtype=np.float64) + self.kt = None if kt is None else np.asarray(kt, dtype=np.float64) + self.bt = None if bt is None else np.asarray(bt, dtype=np.float64) + self.gain_mode = gain_mode + + def has_nuc(self): + """True once the per-unit factory parameters have been provisioned.""" + return all(getattr(self, f) is not None for f in _ARRAY_FIELDS) + + def to_dict(self): + d = {"serial": self.serial, "model": self.model, + "temp_scale": self.temp_scale, "emissivity": self.emissivity, + "t_refl": self.t_refl, "tau": self.tau, + "use_ems_lut": self.use_ems_lut, "ems_version": self.ems_version, + "gain_mode": self.gain_mode} + for f in _ARRAY_FIELDS: + a = getattr(self, f) + d[f] = None if a is None else a.tolist() + return d + + @classmethod + def from_dict(cls, d): + return cls(**d) + + @classmethod + def path_for(cls, serial, directory=DEFAULT_DIR): + return Path(directory) / f"{serial}.json" + + @classmethod + def load(cls, serial, directory=DEFAULT_DIR): + return cls.from_dict(json.loads(cls.path_for(serial, directory).read_text())) + + @classmethod + def try_load(cls, serial, directory=DEFAULT_DIR): + """Load the profile for ``serial`` or None if it hasn't been provisioned.""" + p = cls.path_for(serial, directory) + return cls.from_dict(json.loads(p.read_text())) if p.exists() else None + + def save(self, directory=DEFAULT_DIR): + path = self.path_for(self.serial, directory) + path.parent.mkdir(parents=True, exist_ok=True) + path.write_text(json.dumps(self.to_dict(), indent=2)) + return path + + def corrector(self, tables=None): + """Build a `RadiometricCorrector` configured from this profile.""" + return RadiometricCorrector( + emissivity=self.emissivity, t_refl=self.t_refl, tau=self.tau, + profile=self, temp_scale=self.temp_scale, tables=tables, + use_ems_lut=self.use_ems_lut, ems_version=self.ems_version) diff --git a/src/radiometry/tables.py b/src/radiometry/tables.py new file mode 100644 index 0000000..fb052ad --- /dev/null +++ b/src/radiometry/tables.py @@ -0,0 +1,146 @@ +"""Extract the InfiRay correction lookup tables embedded in libadvirtemp.so. + +The Android `libadvirtemp.so` carries the model-wide emissivity / distance / +target-temperature LUTs as plain `.data` arrays of little-endian float64. They are +constants for the sensor family (not per-unit), so we dump them once into +`adv_tables.npz`, which ships with the package; the runtime never needs the .so. + +Run once to (re)generate the artifact: + + python -m radiometry.tables extract /path/to/libadvirtemp.so + +Everything else loads the shipped npz via `load()`. + +Determined by reverse-engineering (see docs/TC002C-DUO.md): + * arrays are little-endian float64; + * the 2-D ems correction table is row-major [target_temp][org_ems], i.e. its + shape is (len(target_temp_list), len(org_ems_list)) -- confirmed both against + the axis lengths and by which reshape varies smoothly along both axes. +""" + +import struct +import sys +from pathlib import Path + +import numpy as np + +DEFAULT_NPZ = Path(__file__).with_name("adv_tables.npz") + +# 2-D ems-correction tables, each reshaped to (target_temp_axis, org_ems_axis). +EMS_TABLES = { + "ems_correct_table_v1": ("target_temp_list_of_ems_table_v1", + "org_ems_list_of_ems_table_v1"), + "ems_correct_table_v2": ("target_temp_list_of_ems_table_v2", + "org_ems_list_of_ems_table_v2"), +} + +# 1-D axis / lookup tables to carry through verbatim. +LIST_TABLES = [ + "org_ems_list_of_ems_table_v1", "target_temp_list_of_ems_table_v1", + "org_ems_list_of_ems_table_v2", "target_temp_list_of_ems_table_v2", + "dist_table", "new_dist_table", "dist_table_v3", + "temp_table", "target_temp_table", "new_target_temp_table", + "target_temp_table_v3", +] + + +def _read_elf_symbols(path): + """Minimal ELF64 reader: return (file_bytes, {name: (vma, size, file_off)}). + + Maps every defined symbol's virtual address to a file offset through its own + section header, so extraction never depends on a hand-computed delta. + """ + data = Path(path).read_bytes() + if data[:4] != b"\x7fELF" or data[4] != 2: + raise ValueError(f"{path}: not an ELF64 object") + e_shoff, = struct.unpack_from("= len(secs): + continue + sec = secs[st_shndx] + file_off = st_value - sec["addr"] + sec["off"] + syms[cstr(strsec["off"], st_name)] = (st_value, st_size, file_off) + return data, syms + + +def _grab_f64(data, syms, name): + vma, size, off = syms[name] + return np.frombuffer(data[off:off + size], dtype=" caught here. + if not (np.all(np.diff(a) >= 0) and a.max() > a.min()): + raise ValueError(f"{k}: not monotonic non-decreasing -- bad offset/dtype?") + if not (lo <= a.min() and a.max() <= hi): + raise ValueError(f"{k}: out of range [{a.min()},{a.max()}] not in [{lo},{hi}]") + for name, (taxis, eaxis) in EMS_TABLES.items(): + t = tables[name] + if t.shape != (tables[taxis].size, tables[eaxis].size): + raise ValueError(f"{name}: shape {t.shape} != axes " + f"({tables[taxis].size},{tables[eaxis].size})") + if not (0.4 <= t.min() and t.max() <= 1.0): + raise ValueError(f"{name}: values {t.min()}..{t.max()} outside [0.4,1.0]") + + +def extract(so_path, out_path=DEFAULT_NPZ): + """Dump the LUTs from `so_path` (libadvirtemp.so) into `out_path` (.npz).""" + data, syms = _read_elf_symbols(so_path) + tables = {n: _grab_f64(data, syms, n) for n in LIST_TABLES} + for name, (taxis, eaxis) in EMS_TABLES.items(): + flat = _grab_f64(data, syms, name) + tables[name] = flat.reshape(syms[taxis][1] // 8, syms[eaxis][1] // 8) + _validate(tables) + np.savez_compressed(out_path, **tables) + return tables + + +def load(path=DEFAULT_NPZ): + """Load the shipped correction tables as a plain dict of numpy arrays.""" + with np.load(path) as z: + return {k: z[k] for k in z.files} + + +if __name__ == "__main__": + if len(sys.argv) >= 3 and sys.argv[1] == "extract": + out = Path(sys.argv[3]) if len(sys.argv) > 3 else DEFAULT_NPZ + t = extract(sys.argv[2], out) + print(f"Wrote {out} with {len(t)} tables:") + for k, v in sorted(t.items()): + print(f" {k:34} {str(v.shape):12} " + f"[{v.min():.4g}, {v.max():.4g}]") + else: + print(__doc__) + sys.exit("usage: python -m radiometry.tables extract [out.npz]") From 33991c7efb3014a78d2cc821e1e46589970253d1 Mon Sep 17 00:00:00 2001 From: Aron Novak Date: Tue, 16 Jun 2026 11:11:31 +0200 Subject: [PATCH 5/5] =?UTF-8?q?TC002C=20Duo:=20decode=20real=20absolute=20?= =?UTF-8?q?=C2=B0C=20from=20the=20self-calibrating=20stream?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Open the camera's 8x12578 mode (header + 256x192 temperature plane + image) and decode temperature as raw/64 - K, validated against the official app to 0.1 °C (the captured frame's app min/max/avg 22.6/34.7/ 27.8 reproduce exactly). radiometry/duo.py + thermalcam auto-detect/wiring; white_hot grayscale palette default for cleaner detector input. Co-Authored-By: Claude Opus 4.8 (1M context) --- src/radiometry/duo.py | 85 ++++++++++++++++++++++++++++++++++++++++++ src/thermal_bridge.py | 4 +- src/thermalcam.py | 87 ++++++++++++++++++++++++++++++++++--------- 3 files changed, 157 insertions(+), 19 deletions(-) create mode 100644 src/radiometry/duo.py diff --git a/src/radiometry/duo.py b/src/radiometry/duo.py new file mode 100644 index 0000000..b4ad69e --- /dev/null +++ b/src/radiometry/duo.py @@ -0,0 +1,85 @@ +"""Topdon TC002C Duo radiometric decode: raw stream frame -> temperature (°C). + +The Duo's true-temperature mode is the UVC resolution the camera advertises as the +bogus-looking ``8x12578`` (= 100624 "pixels"; as YUYV that's 201248 bytes). Each +frame is 100624 little-endian u16: + + [0:2320] header / telemetry (magic 0x70827773, dims, scene config, and + the app's own min/max/avg of the scene as float32 at byte 144) + [2320:51472] temperature plane, 256x192, raw sensor counts + [51472:100624] display image, 256x192 (unpopulated over plain UVC) + +In the indoor operating range the firmware's conversion is linear in the raw counts: + + temp_C = raw / 64 - K + +The gain (1/64 °C per count) is validated **exactly** against the official app: for a +captured frame whose header carries the app's own min/max/avg (22.6 / 34.7 / 27.8 °C), +``raw_min/64 - 50 = 22.62`` and ``raw_max/64 - 50 = 34.69`` -- a 0.1 °C match on all +three. (The camera's full kt/bt + LUT pipeline, reverse-engineered from +libadvirtempac020.so, reduces to this line in-range; the LUT only matters at +temperature extremes outside indoor monitoring.) + +``K`` is a per-frame offset that tracks the sensor's Vtemp / focal-plane temperature +(it drifts as the sensor warms). The firmware reads Vtemp from an internal sensor we +can't see over plain UVC, so `estimate_offset` pins a robust cold-baseline percentile +of the frame to an assumed ambient (`ambient_bg`); this auto-tracks Vtemp drift and is +exact for temperature *differences* (e.g. person-vs-background, the reflection gate), +while the absolute level can be fine-tuned with thermalcam's ``--temp-offset`` / ``[`` +``]`` keys. Emissivity / reflected-temperature correction is applied on top by +`radiometry.corrector.RadiometricCorrector`. +""" + +import numpy as np + +GAIN_DIV = 64.0 # °C per raw count = 1/64 (validated vs the app) +FRAME_U16 = 100624 # 201248 bytes / 2 +FRAME_MAGIC = 0x70827773 # header[0:2] as a little-endian u32 +HEADER_U16 = 2320 +SENSOR_W, SENSOR_H = 256, 192 +_PLANE = SENSOR_W * SENSOR_H # 49152 +TEMP_OFFSET = HEADER_U16 # 2320 +IMAGE_OFFSET = HEADER_U16 + _PLANE # 51472 + +# `estimate_offset` pins this low percentile of the raw plane (the coldest stable +# background, ~ambient) to `AMBIENT_BG` °C. The captured calibration frame's cold +# baseline sat at ~23 °C, which reproduces its app min/max/avg; rooms differ, so this +# is the absolute anchor to nudge via --temp-offset (the gain above is fixed/exact). +BG_PERCENTILE = 2.0 +AMBIENT_BG = 22.0 + + +def is_duo_frame(u16): + """True if a flat u16 buffer looks like a Duo radiometric frame (size + magic).""" + u = np.asarray(u16).reshape(-1) + return u.size >= FRAME_U16 and (int(u[0]) | (int(u[1]) << 16)) == FRAME_MAGIC + + +def split_frame(u16): + """Split a flat u16 frame into (header, temp_raw[192,256], image[192,256] u8).""" + u = np.asarray(u16).reshape(-1)[:FRAME_U16] + header = u[:HEADER_U16] + temp = u[TEMP_OFFSET:TEMP_OFFSET + _PLANE].reshape(SENSOR_H, SENSOR_W) + image = (u[IMAGE_OFFSET:IMAGE_OFFSET + _PLANE] & 0xFF).astype(np.uint8).reshape(SENSOR_H, SENSOR_W) + return header, temp, image + + +def estimate_offset(temp_raw, ambient_bg=AMBIENT_BG): + """Per-frame offset K so the cold background sits at ``ambient_bg`` °C. + + Tracks the sensor's Vtemp drift (the whole plane shifts with focal-plane + temperature) by anchoring a robust low percentile to ambient. + """ + return float(np.percentile(temp_raw, BG_PERCENTILE)) / GAIN_DIV - ambient_bg + + +def apparent_celsius(temp_raw, offset=None, ambient_bg=AMBIENT_BG): + """Raw temperature-plane counts -> temperature (°C): ``raw/64 - K``. + + ``offset`` (K) is auto-estimated from the frame's cold baseline when ``None``. + Returns a float32 array (apparent temperature, before emissivity correction). + """ + raw = np.asarray(temp_raw, dtype=np.float32) + if offset is None: + offset = estimate_offset(raw, ambient_bg) + return raw / GAIN_DIV - np.float32(offset) diff --git a/src/thermal_bridge.py b/src/thermal_bridge.py index a51f759..890c2c6 100644 --- a/src/thermal_bridge.py +++ b/src/thermal_bridge.py @@ -143,7 +143,9 @@ def main(argv=None): p = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) p.add_argument("--host", default="127.0.0.1", help="Bind address (0.0.0.0 to expose on the LAN).") p.add_argument("--port", type=int, default=8090) - p.add_argument("--colormap", default="inferno", help="Palette (jet, inferno, hot, ...). Default inferno.") + p.add_argument("--colormap", default="white_hot", + help="Palette (white_hot, black_hot, inferno, jet, ...). Default white_hot - " + "grayscale gives detectors far fewer false positives on empty scenes.") p.add_argument("--width", type=int, help="Output width (default: camera native).") p.add_argument("--height", type=int, help="Output height (default: camera native).") p.add_argument("--quality", type=int, default=85, help="JPEG quality 1-100.") diff --git a/src/thermalcam.py b/src/thermalcam.py index 06cb524..8784eda 100644 --- a/src/thermalcam.py +++ b/src/thermalcam.py @@ -38,6 +38,14 @@ " pip install opencv-python numpy (everything else)" ) +# TC002C Duo true-temperature decode (its own self-calibrating stream + bundled +# calibration LUTs). Optional: the viewer still runs without it on other cameras. +sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) +try: + from radiometry import duo as _duo +except Exception: # pragma: no cover - radiometry is optional + _duo = None + # USB vendor IDs known to ship InfiRay-based thermal cameras. 0x2BDF is Topdon. KNOWN_THERMAL_VIDS = {0x2BDF, 0x0BDA, 0x1514, 0x3474} @@ -51,23 +59,42 @@ (256, 384): ((192, 384), (0, 192), (256, 192), 64.0, "le"), # TC001 / P2 family (512, 484): ((0, 96), (98, 482), (256, 192), 16.0, "le"), # TC002C Duo } -# Resolutions to try for genuine temperature data, in order of preference. -RADIOMETRIC_MODES = [(256, 384), (512, 484)] +# The TC002C Duo's true-°C mode: the camera advertises it as the bogus-looking +# "8x12578" (= 100624 u16 = 201248 B), a header+temperature+image frame the firmware +# self-calibrates. Decoded by radiometry/duo.py (real °C, no per-unit calibration). +DUO_MODE = (8, 12578) +# Resolutions to try for genuine temperature data, in order of preference. The Duo's +# self-calibrating mode (real absolute °C) is preferred over the 512x484 linear band. +RADIOMETRIC_MODES = ([DUO_MODE] if _duo is not None else []) + [(256, 384), (512, 484)] # Plain image modes (no temperature), used only if no radiometric mode works. IMAGE_MODES = [(512, 384), (256, 392), (256, 192)] +# (cv2 colormap or None for the grayscale white/black-hot pseudo-palettes), name. +# white-hot is first (the default): on flat/noisy scenes a colour map turns sensor +# noise into person-shaped blobs that fool detectors, while grayscale stays clean. COLORMAPS = [ - (cv2.COLORMAP_JET, "Jet"), - (cv2.COLORMAP_INFERNO, "Inferno"), - (cv2.COLORMAP_HOT, "Hot"), - (cv2.COLORMAP_MAGMA, "Magma"), - (cv2.COLORMAP_PLASMA, "Plasma"), - (cv2.COLORMAP_BONE, "Bone"), - (cv2.COLORMAP_VIRIDIS, "Viridis"), - (cv2.COLORMAP_RAINBOW, "Rainbow"), + (None, "white_hot"), + (None, "black_hot"), + (cv2.COLORMAP_INFERNO, "inferno"), + (cv2.COLORMAP_JET, "jet"), + (cv2.COLORMAP_HOT, "ironbow"), + (cv2.COLORMAP_MAGMA, "magma"), + (cv2.COLORMAP_VIRIDIS, "viridis"), + (cv2.COLORMAP_BONE, "bone"), ] +def apply_palette(luma, idx): + """Colormap an 8-bit luma image; returns (bgr, name). Handles the grayscale + white-hot / black-hot pseudo-palettes (no cv2 colormap).""" + cmap, name = COLORMAPS[idx] + if cmap is None: + if name == "black_hot": + luma = 255 - luma + return cv2.cvtColor(luma, cv2.COLOR_GRAY2BGR), name + return cv2.applyColorMap(luma, cmap), name + + def is_raspberrypi(): try: with open("/sys/firmware/devicetree/base/model", "r") as m: @@ -253,8 +280,26 @@ def _decode_temp(raw, temp_rows, sensor, scale, order): return temp if _plausible(temp) else None +def _parse_duo(raw): + """Decode a TC002C Duo radiometric frame (the 8x12578 mode) into a Frame. + + The temperature plane carries the thermal scene (real °C from the firmware's + self-calibrating decode); the raw display-image plane is unpopulated over plain + UVC, so the temperature field doubles as the display source. + """ + u16 = np.frombuffer(np.ascontiguousarray(raw).tobytes(), dtype="= _duo.FRAME_U16 * 2 and _duo.is_duo_frame( + np.frombuffer(np.ascontiguousarray(raw).tobytes(), dtype=" 0: luma = cv2.blur(luma, (self.blur, self.blur)) - cmap, cmap_name = COLORMAPS[self.colormap] - heatmap = cv2.applyColorMap(luma, cmap) + heatmap, cmap_name = apply_palette(luma, self.colormap) self._draw_crosshair(heatmap, new_w, new_h, center, unit) self._draw_markers(heatmap, max_pos, min_pos, tmax, tmin, tavg, unit, w, h, new_w, new_h) @@ -546,8 +598,7 @@ def colormap_frame(self, raw): luma = cv2.resize(luma, (nw, nh), interpolation=cv2.INTER_CUBIC) if self.blur > 0: luma = cv2.blur(luma, (self.blur, self.blur)) - cmap, _ = COLORMAPS[self.colormap] - return cv2.applyColorMap(luma, cmap) + return apply_palette(luma, self.colormap)[0] def _draw_crosshair(self, img, w, h, center_temp, unit): cx, cy = w // 2, h // 2