diff --git a/.github/workflows/release.yml b/.github/workflows/release.yml new file mode 100644 index 0000000..3dd6928 --- /dev/null +++ b/.github/workflows/release.yml @@ -0,0 +1,63 @@ +# Builds the person model and publishes it as a release. Run by hand +# (Actions -> release -> Run workflow) and only that way: release assets are +# never uploaded from a workstation or by any check that runs on a pull +# request or a push, so each one traces to the commit and run that built it. +name: release + +on: + workflow_dispatch: + inputs: + tag: + description: Release tag to create (e.g. v1). Must not exist yet. + required: true + +permissions: + contents: write + +concurrency: + group: release + cancel-in-progress: false + +jobs: + build: + runs-on: ubuntu-22.04 + timeout-minutes: 120 + steps: + - uses: actions/checkout@v4 + + - uses: actions/setup-python@v5 + with: + python-version: "3.11" + + - name: Dependencies + run: | + pip install -q -r requirements.txt + sudo apt-get update -qq + sudo apt-get install -y -qq gcc-arm-linux-gnueabihf qemu-user-static >/dev/null + + - uses: actions/cache@v4 + with: + path: .cache + key: coco-${{ hashFiles('coco.txt') }} + + - name: Build, check and measure + run: bash scripts/build.sh out + + - name: Publish + env: + GH_TOKEN: ${{ github.token }} + TAG: ${{ inputs.tag }} + run: | + { + echo "Built by run ${{ github.run_id }} from ${{ github.sha }}." + echo + echo "\`person.pdnet\` is the int8 person detector majestic's CPU engine loads;" + echo "\`person.json\`, installed beside it, is what its three object-detection" + echo "presets mean for this model. Both the host engine and its ARMv7 NEON" + echo "build were checked bit for bit against the reference arithmetic before" + echo "anything was measured. Measured on labelled crops of COCO val2017:" + echo + cat out/points.md + } > notes.md + gh release create "$TAG" --target "${{ github.sha }}" --title "$TAG" \ + --notes-file notes.md out/person.pdnet out/person.json out/SHA256SUMS out/points.md diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..83a4942 --- /dev/null +++ b/.gitignore @@ -0,0 +1,3 @@ +.cache/ +out/ +__pycache__/ diff --git a/README.md b/README.md index 25fabb5..7eb5658 100644 --- a/README.md +++ b/README.md @@ -1,2 +1,64 @@ # cpu-models -CPU-run detection models for OpenIPC cameras, built and published reproducibly + +Detection models that OpenIPC cameras run on their CPU, with majestic's own +int8 engine, and the release that publishes them. The first is the person +detector majestic's object detection uses on cameras without an NPU. + +## The person detector + +`models/person/model.onnx` is a statically quantized (QDQ) int8 network: +MobileNetV1-0.25 with an FPN-lite head, 66.7k parameters, a 320x192 YUV +input, and boxes at strides 8, 16 and 32. It was trained for majestic on COCO +2017 person boxes and fine-tuned on hard negatives mined from a lab camera's +own crops; while those mined crops were sorted, an Ultralytics YOLO model was +used to find which ones contained a person. + +The camera hands it a crop of the motion region at the scale it was trained +for, in NV21; it reports people with a score from 0 to 100. + +## Releases + +Published only by the `release` workflow (Actions -> release -> Run workflow, +with a new tag), never by hand. A release carries: + +- `person.pdnet`: the model as majestic's engine loads it, exported from the + ONNX by `tools/export_pdnet.py`; +- `person.json`: the score thresholds behind majestic's three presets, for this + model, installed beside it: + - `balanced`: the F1-optimal threshold; + - `fewerFalseAlarms`: the F0.5 optimum; + - `catchEverything`: the highest threshold keeping recall at 0.90, else the + F2 optimum; +- `points.md`: what each preset achieves; and `SHA256SUMS`. + +## Building + +```sh +pip install -r requirements.txt # Python 3.11 +sudo apt-get install gcc-arm-linux-gnueabihf qemu-user-static +bash scripts/build.sh [out] +``` + +`scripts/build.sh` is the whole build, the same here and in the workflow. It: + +1. fetches COCO val2017 and its annotations, hash-checked (`coco.txt`); +2. exports the blob; +3. builds the engine (`engine/pdnet.c`) for the host and for ARMv7 with NEON, + the camera's path; +4. checks both against `tools/pdq.py`, the model's arithmetic written out in + numpy. On 40 crops, every value of all six head maps must be identical, + and the ARM build, run under qemu-arm, must find exactly the host's boxes + on all 3000 crops; +5. measures the presets: `tools/crops.py` cuts 3000 labelled 320x192 crops + from COCO val2017 in the camera's geometry, 2000 around a person and 1000 + with nobody in them, and `tools/points.py` scores the engine's detections + on them. + +A build of one commit is byte-identical anywhere. + +`engine/` is the engine majestic carries; a change to it lands here and there +together. + +## License + +MIT. diff --git a/coco.txt b/coco.txt new file mode 100644 index 0000000..a4f7ccc --- /dev/null +++ b/coco.txt @@ -0,0 +1,4 @@ +# COCO val2017 and its annotations, for the presets' operating points, pinned +# by hash. file sha256 url +val2017.zip 4f7e2ccb2866ec5041993c9cf2a952bbed69647b115d0f74da7ce8f4bef82f05 http://images.cocodataset.org/zips/val2017.zip +annotations_trainval2017.zip 113a836d90195ee1f884e704da6304dfaaecff1f023f49b6ca93c4aaae470268 http://images.cocodataset.org/annotations/annotations_trainval2017.zip diff --git a/engine/pdnet.c b/engine/pdnet.c new file mode 100644 index 0000000..038b3a9 --- /dev/null +++ b/engine/pdnet.c @@ -0,0 +1,836 @@ +#define _POSIX_C_SOURCE 200809L +#include "pdnet.h" + +#include +#include +#include +#include +#include + +#if defined(__ARM_NEON) || defined(__ARM_NEON__) +#include +#define PDNET_NEON 1 +#endif + +#define PDNET_MAGIC "MAJPDNT1" +#define MAX_TENSORS 64 +#define MAX_OPS 64 +#define MAX_LEVELS 3 +#define NMS_IOU 0.55f + +#ifdef PDNET_PROFILE +#include +double pdnet_prof_us[MAX_OPS]; +void pdnet_profile_dump(const PdNet *n, int runs); +#endif + +enum { OP_CONV = 1, OP_RESIZE2X = 2, OP_ADD = 3 }; + +typedef struct { + int h, w, c; + int8_t *q; /* int8 activation, HWC */ + float *f; /* float head output, HWC */ +} Tensor; + +typedef struct { + int type; + int in, in2, out; + /* conv */ + int k, stride, pad, dw, cin, cout, float_out; + int32_t in_zp, out_zp; + const int32_t *ms; /* (mult, shift) pairs */ + const float *scale; /* float_out */ + const int32_t *bias; + const int8_t *w; + int32_t *bias_folded; /* 1x1: bias - in_zp * sum(w) */ + /* For the NEON kernels, prepared at load: weights widened to int16 with + * the output channel innermost (1x1: [cin][cout]; kxk: [k][k][cin][cout]; + * depthwise: [k][k][c]), and each channel's requantization as a + * multiplier and the negated right shift vrshlq_s64 wants. */ + int16_t *w16; + int32_t *mult; + int64_t *nshift; + int kpad; /* kxk: k*k*cin rounded up to 8, the im2col row */ + /* add */ + int32_t zp_a, zp_b, zp_o, m_a, s_a, m_b, s_b, m_o, s_o; +} Op; + +struct PdNet { + uint8_t *blob; + int in_w, in_h; + int8_t lut_y[256], lut_c[256]; + int n_t, n_ops, n_lv; + Tensor t[MAX_TENSORS]; + Op op[MAX_OPS]; + struct { int stride, cls, reg; } lv[MAX_LEVELS]; + int32_t *acc; + int16_t *col; /* im2col rows for the NEON kxk */ +}; + +/* ---- arithmetic ---------------------------------------------------------- */ + +/* acc * mult * 2^(shift - 31), rounded once, half up. */ +static inline int32_t requant(int32_t acc, int32_t mult, int32_t shift) { + const int n = 31 - shift; + const int64_t p = (int64_t)acc * mult; + return (int32_t)((p + ((int64_t)1 << (n - 1))) >> n); +} + +static inline int8_t clamp8(int32_t v) { + return (int8_t)(v < -128 ? -128 : v > 127 ? 127 : v); +} + +/* ---- loading --------------------------------------------------------------- */ + +typedef struct { + const uint8_t *p, *end; + bool bad; +} Rd; + +static const void *take(Rd *r, size_t n) { + if (r->bad || (size_t)(r->end - r->p) < n) { + r->bad = true; + return NULL; + } + const void *v = r->p; + r->p += n; + return v; +} + +static uint32_t u32(Rd *r) { + const uint8_t *b = take(r, 4); + return b ? (uint32_t)b[0] | (uint32_t)b[1] << 8 | (uint32_t)b[2] << 16 | + (uint32_t)b[3] << 24 + : 0; +} + +static int32_t i32(Rd *r) { return (int32_t)u32(r); } + +static float f32(Rd *r) { + uint32_t v = u32(r); + float f; + memcpy(&f, &v, 4); + return f; +} + +static bool tensor_ok(const PdNet *n, uint32_t i) { return i < (uint32_t)n->n_t; } + +static bool parse(PdNet *n, Rd *r) { + const void *m = take(r, 8); + if (!m || memcmp(m, PDNET_MAGIC, 8) || u32(r) != 1) + return false; + n->in_w = (int)u32(r); + n->in_h = (int)u32(r); + const float in_scale = f32(r); + const int32_t in_zp = i32(r); + if (r->bad || n->in_w <= 0 || n->in_w > 4096 || n->in_h <= 0 || n->in_h > 4096 || + !(in_scale > 0)) + return false; + /* The frame's bytes to the network's int8 input, as the float model's + * preprocessing defines it: Y from BT.709 limited range to -1..1, U and + * V from their offset to about -1..1. Two tables of 256. */ + const float ky = (float)(255.0 / 219.0 / 128.0), kc = (float)(255.0 / 224.0 / 128.0); + for (int v = 0; v < 256; v++) { + const float y = ((float)v - 16.0f) * ky - 1.0f; + const float c = ((float)v - 128.0f) * kc; + n->lut_y[v] = clamp8((int32_t)rintf(y / in_scale) + in_zp); + n->lut_c[v] = clamp8((int32_t)rintf(c / in_scale) + in_zp); + } + n->n_t = (int)u32(r); + if (r->bad || n->n_t <= 0 || n->n_t > MAX_TENSORS) + return false; + for (int i = 0; i < n->n_t; i++) { + n->t[i].h = (int)u32(r); + n->t[i].w = (int)u32(r); + n->t[i].c = (int)u32(r); + if (r->bad || n->t[i].h <= 0 || n->t[i].h > 4096 || n->t[i].w <= 0 || + n->t[i].w > 4096 || n->t[i].c <= 0 || n->t[i].c > 4096) + return false; + } + n->n_ops = (int)u32(r); + if (r->bad || n->n_ops <= 0 || n->n_ops > MAX_OPS) + return false; + for (int i = 0; i < n->n_ops; i++) { + Op *o = &n->op[i]; + o->type = (int)u32(r); + if (o->type == OP_CONV) { + o->in = (int)u32(r); + o->out = (int)u32(r); + o->k = (int)u32(r); + o->stride = (int)u32(r); + o->pad = (int)u32(r); + o->dw = (int)u32(r); + o->cin = (int)u32(r); + o->cout = (int)u32(r); + o->float_out = (int)u32(r); + o->in_zp = i32(r); + o->out_zp = i32(r); + if (r->bad || !tensor_ok(n, o->in) || !tensor_ok(n, o->out) || o->k < 1 || + o->k > 7 || o->stride < 1 || o->stride > 4 || o->pad > o->k) + return false; + const Tensor *ti = &n->t[o->in], *to = &n->t[o->out]; + if (ti->c != o->cin || to->c != o->cout || (o->dw && o->cin != o->cout) || + to->h != (ti->h + 2 * o->pad - o->k) / o->stride + 1 || + to->w != (ti->w + 2 * o->pad - o->k) / o->stride + 1) + return false; + if (o->float_out) + o->scale = take(r, (size_t)o->cout * 4); + else + o->ms = take(r, (size_t)o->cout * 8); + o->bias = take(r, (size_t)o->cout * 4); + const size_t wn = o->dw ? (size_t)o->k * o->k * o->cout + : (size_t)o->cout * o->k * o->k * o->cin; + o->w = take(r, wn); + take(r, (4 - wn % 4) % 4); + if (r->bad) + return false; + if (!o->float_out) + for (int c = 0; c < o->cout; c++) { + const int32_t sh = o->ms[2 * c + 1]; + if (sh > 30 || sh < -31) + return false; + } + } else if (o->type == OP_RESIZE2X) { + o->in = (int)u32(r); + o->out = (int)u32(r); + if (r->bad || !tensor_ok(n, o->in) || !tensor_ok(n, o->out)) + return false; + const Tensor *ti = &n->t[o->in], *to = &n->t[o->out]; + if (to->h != 2 * ti->h || to->w != 2 * ti->w || to->c != ti->c) + return false; + } else if (o->type == OP_ADD) { + o->in = (int)u32(r); + o->in2 = (int)u32(r); + o->out = (int)u32(r); + o->zp_a = i32(r); o->zp_b = i32(r); o->zp_o = i32(r); + o->m_a = i32(r); o->s_a = i32(r); + o->m_b = i32(r); o->s_b = i32(r); + o->m_o = i32(r); o->s_o = i32(r); + if (r->bad || !tensor_ok(n, o->in) || !tensor_ok(n, o->in2) || + !tensor_ok(n, o->out)) + return false; + const Tensor *a = &n->t[o->in], *b = &n->t[o->in2], *to = &n->t[o->out]; + if (a->h != b->h || a->w != b->w || a->c != b->c || to->h != a->h || + to->w != a->w || to->c != a->c) + return false; + } else { + return false; + } + } + n->n_lv = (int)u32(r); + if (r->bad || n->n_lv <= 0 || n->n_lv > MAX_LEVELS) + return false; + for (int i = 0; i < n->n_lv; i++) { + n->lv[i].stride = (int)u32(r); + n->lv[i].cls = (int)u32(r); + n->lv[i].reg = (int)u32(r); + if (r->bad || !tensor_ok(n, n->lv[i].cls) || !tensor_ok(n, n->lv[i].reg) || + n->t[n->lv[i].cls].c != 1 || n->t[n->lv[i].reg].c != 4) + return false; + } + return !r->bad && n->t[0].h == n->in_h && n->t[0].w == n->in_w && n->t[0].c == 3; +} + +/* Widened, channel-innermost weights and vector requantization constants + * for a convolution's NEON kernel. Harmless without NEON. */ +static bool prepare(Op *o) { + if (o->type != OP_CONV) + return true; + const int k = o->k, ci = o->cin, co = o->cout; + const size_t wn = o->dw ? (size_t)k * k * co : (size_t)co * k * k * ci; + if (!(o->w16 = malloc(wn * sizeof(int16_t)))) + return false; + if (o->dw) { + for (size_t i = 0; i < wn; i++) + o->w16[i] = o->w[i]; + } else { + /* [k][k][cin] rows, cout innermost; a kxk's rows padded with zeros to + * a multiple of 8, the im2col row its kernel reads 8 at a time. */ + o->kpad = (k * k * ci + 7) & ~7; + free(o->w16); + if (!(o->w16 = calloc((size_t)o->kpad * co, sizeof(int16_t)))) + return false; + for (int c = 0; c < co; c++) + for (int t = 0; t < k * k; t++) + for (int i = 0; i < ci; i++) + o->w16[((size_t)t * ci + i) * co + c] = + o->w[((size_t)c * k * k + t) * ci + i]; + } + if (o->float_out) + return true; + if (!(o->mult = malloc((size_t)co * sizeof(int32_t))) || + !(o->nshift = malloc((size_t)co * sizeof(int64_t)))) + return false; + for (int c = 0; c < co; c++) { + o->mult[c] = o->ms[2 * c]; + o->nshift[c] = -(int64_t)(31 - o->ms[2 * c + 1]); + } + return true; +} + +PdNet *pdnet_open(const char *path) { + FILE *f = fopen(path, "rb"); + if (!f) + return NULL; + PdNet *n = calloc(1, sizeof(*n)); + long len = -1; + if (n && !fseek(f, 0, SEEK_END) && (len = ftell(f)) > 0 && len < (16 << 20) && + !fseek(f, 0, SEEK_SET) && (n->blob = malloc((size_t)len)) && + fread(n->blob, 1, (size_t)len, f) == (size_t)len) { + Rd r = {n->blob, n->blob + len, false}; + if (parse(n, &r)) { + fclose(f); + f = NULL; + /* Buffers: every tensor its own, which for this network is a few + * hundred KB and spares a liveness analysis. */ + bool ok = true; + for (int i = 0; ok && i < n->n_t; i++) { + const size_t e = (size_t)n->t[i].h * n->t[i].w * n->t[i].c; + bool is_float = false; + for (int j = 0; j < n->n_ops; j++) + if (n->op[j].type == OP_CONV && n->op[j].out == i && n->op[j].float_out) + is_float = true; + if (is_float) + ok = (n->t[i].f = malloc(e * sizeof(float))) != NULL; + else + ok = (n->t[i].q = malloc(e)) != NULL; + } + int max_c = 0; + for (int j = 0; j < n->n_ops; j++) + if (n->op[j].type == OP_CONV && n->op[j].cout > max_c) + max_c = n->op[j].cout; + ok = ok && (n->acc = malloc((size_t)max_c * sizeof(int32_t))); + for (int j = 0; ok && j < n->n_ops; j++) { + Op *o = &n->op[j]; + if (o->type != OP_CONV || o->k != 1 || o->dw) + continue; + if (!(o->bias_folded = malloc((size_t)o->cout * sizeof(int32_t)))) { + ok = false; + break; + } + for (int c = 0; c < o->cout; c++) { + int32_t s = 0; + for (int i = 0; i < o->cin; i++) + s += o->w[(size_t)c * o->cin + i]; + o->bias_folded[c] = o->bias[c] - o->in_zp * s; + } + } + size_t col = 0; + for (int j = 0; ok && j < n->n_ops; j++) { + ok = prepare(&n->op[j]); + const Op *o = &n->op[j]; + if (ok && o->type == OP_CONV && !o->dw && o->k > 1) { + const Tensor *to = &n->t[o->out]; + const size_t e = (size_t)to->w * o->kpad; + if (e > col) + col = e; + } + } + ok = ok && (!col || (n->col = malloc(col * sizeof(int16_t)))); + if (ok) + return n; + } + } + if (f) + fclose(f); + pdnet_close(n); + return NULL; +} + +void pdnet_close(PdNet *n) { + if (!n) + return; + for (int i = 0; i < n->n_t; i++) { + free(n->t[i].q); + free(n->t[i].f); + } + for (int j = 0; j < n->n_ops; j++) { + free(n->op[j].bias_folded); + free(n->op[j].w16); + free(n->op[j].mult); + free(n->op[j].nshift); + } + free(n->acc); + free(n->col); + free(n->blob); + free(n); +} + +int pdnet_width(const PdNet *n) { return n->in_w; } +int pdnet_height(const PdNet *n) { return n->in_h; } + +/* ---- the layers -------------------------------------------------------------- */ + +static inline void emit(const Op *o, int32_t acc, int c, int8_t *q, float *f) { + if (o->float_out) + *f = (float)((double)acc * (double)o->scale[c]); + else + *q = clamp8(requant(acc, o->ms[2 * c], o->ms[2 * c + 1]) + o->out_zp); +} + +/* Sum of a[i] * w[i] over n int8 pairs. */ +static inline int32_t dot_i8(const int8_t *a, const int8_t *w, int n) { + int i = 0; + int32_t s = 0; +#ifdef PDNET_NEON + int32x4_t c = vdupq_n_s32(0); + for (; i + 16 <= n; i += 16) { + const int8x16_t av = vld1q_s8(a + i), wv = vld1q_s8(w + i); + c = vpadalq_s16(c, vmull_s8(vget_low_s8(av), vget_low_s8(wv))); + c = vpadalq_s16(c, vmull_s8(vget_high_s8(av), vget_high_s8(wv))); + } + for (; i + 8 <= n; i += 8) + c = vpadalq_s16(c, vmull_s8(vld1_s8(a + i), vld1_s8(w + i))); + s = vgetq_lane_s32(c, 0) + vgetq_lane_s32(c, 1) + vgetq_lane_s32(c, 2) + + vgetq_lane_s32(c, 3); +#endif + for (; i < n; i++) + s += a[i] * w[i]; + return s; +} + +#ifdef PDNET_NEON +/* Eight channels' accumulators to int8: acc * mult rounded once by the + * channel's shift (vrshlq_s64 adds the half before shifting, as requant() + * does), plus the zero point, saturated. Bit for bit requant() + clamp8(). */ +static inline int8x8_t requant8(const Op *o, int c, int32x4_t lo, int32x4_t hi) { + const int32x4_t m0 = vld1q_s32(o->mult + c), m1 = vld1q_s32(o->mult + c + 4); + const int64_t *ns = o->nshift + c; + int64x2_t p0 = vrshlq_s64(vmull_s32(vget_low_s32(lo), vget_low_s32(m0)), vld1q_s64(ns)); + int64x2_t p1 = vrshlq_s64(vmull_s32(vget_high_s32(lo), vget_high_s32(m0)), vld1q_s64(ns + 2)); + int64x2_t p2 = vrshlq_s64(vmull_s32(vget_low_s32(hi), vget_low_s32(m1)), vld1q_s64(ns + 4)); + int64x2_t p3 = vrshlq_s64(vmull_s32(vget_high_s32(hi), vget_high_s32(m1)), vld1q_s64(ns + 6)); + const int32x4_t zp = vdupq_n_s32(o->out_zp); + const int32x4_t r0 = vaddq_s32(vcombine_s32(vqmovn_s64(p0), vqmovn_s64(p1)), zp); + const int32x4_t r1 = vaddq_s32(vcombine_s32(vqmovn_s64(p2), vqmovn_s64(p3)), zp); + return vqmovn_s16(vcombine_s16(vqmovn_s32(r0), vqmovn_s32(r1))); +} + +/* The multiply both 1x1 and kxk come down to: npix rows of K inputs (K a + * multiple of 8) against a [K][cout] int16 matrix, cout a multiple of 8. + * Four pixels share every weight load; inputs come in 8 at a time and are + * used lane by lane (vmlal_lane_s16), so nothing crosses from the NEON + * registers to the core ones in the inner loop. IN is int8_t (a 1x1's + * activations, zero point folded into its bias) or int16_t (an im2col row, + * zero point already taken out). */ +#define GEMM_LOAD8_int8_t(p) vmovl_s8(vld1_s8(p)) +#define GEMM_LOAD8_int16_t(p) vld1q_s16(p) +#define DEFINE_GEMM(IN) \ + static void gemm_##IN(const Op *o, const IN *a, int K, const int32_t *bias, \ + int8_t *out, int npix) { \ + const int cout = o->cout; \ + int m = 0; \ + for (; m + 4 <= npix; m += 4) { \ + const IN *a0 = a + (size_t)m * K, *a1 = a0 + K, *a2 = a1 + K, *a3 = a2 + K; \ + for (int c = 0; c < cout; c += 8) { \ + const int32x4_t b0 = vld1q_s32(bias + c), b1 = vld1q_s32(bias + c + 4); \ + int32x4_t l0 = b0, h0 = b1, l1 = b0, h1 = b1; \ + int32x4_t l2 = b0, h2 = b1, l3 = b0, h3 = b1; \ + const int16_t *w = o->w16 + c; \ + for (int i = 0; i < K; i += 8) { \ + const int16x8_t x0 = GEMM_LOAD8_##IN(a0 + i); \ + const int16x8_t x1 = GEMM_LOAD8_##IN(a1 + i); \ + const int16x8_t x2 = GEMM_LOAD8_##IN(a2 + i); \ + const int16x8_t x3 = GEMM_LOAD8_##IN(a3 + i); \ + GEMM_STEP(0, vget_low_s16) GEMM_STEP(1, vget_low_s16) \ + GEMM_STEP(2, vget_low_s16) GEMM_STEP(3, vget_low_s16) \ + GEMM_STEP(0, vget_high_s16) GEMM_STEP(1, vget_high_s16) \ + GEMM_STEP(2, vget_high_s16) GEMM_STEP(3, vget_high_s16) \ + } \ + vst1_s8(out + (size_t)(m + 0) * cout + c, requant8(o, c, l0, h0)); \ + vst1_s8(out + (size_t)(m + 1) * cout + c, requant8(o, c, l1, h1)); \ + vst1_s8(out + (size_t)(m + 2) * cout + c, requant8(o, c, l2, h2)); \ + vst1_s8(out + (size_t)(m + 3) * cout + c, requant8(o, c, l3, h3)); \ + } \ + } \ + for (; m < npix; m++) { \ + const IN *a0 = a + (size_t)m * K; \ + for (int c = 0; c < cout; c += 8) { \ + int32x4_t lo = vld1q_s32(bias + c), hi = vld1q_s32(bias + c + 4); \ + const int16_t *w = o->w16 + c; \ + for (int i = 0; i < K; i++, w += cout) { \ + const int16x8_t wv = vld1q_s16(w); \ + lo = vmlal_n_s16(lo, vget_low_s16(wv), a0[i]); \ + hi = vmlal_n_s16(hi, vget_high_s16(wv), a0[i]); \ + } \ + vst1_s8(out + (size_t)m * cout + c, requant8(o, c, lo, hi)); \ + } \ + } \ + } +/* One half (4 inputs) of an 8-input step: four weight rows, each against + * its input lane, for all four pixels. */ +#define GEMM_LANE(P, HALF, L) \ + { const int16x8_t wv = vld1q_s16(w); w += cout; \ + l0 = vmlal_lane_s16(l0, vget_low_s16(wv), HALF(x0), L); \ + h0 = vmlal_lane_s16(h0, vget_high_s16(wv), HALF(x0), L); \ + l1 = vmlal_lane_s16(l1, vget_low_s16(wv), HALF(x1), L); \ + h1 = vmlal_lane_s16(h1, vget_high_s16(wv), HALF(x1), L); \ + l2 = vmlal_lane_s16(l2, vget_low_s16(wv), HALF(x2), L); \ + h2 = vmlal_lane_s16(h2, vget_high_s16(wv), HALF(x2), L); \ + l3 = vmlal_lane_s16(l3, vget_low_s16(wv), HALF(x3), L); \ + h3 = vmlal_lane_s16(h3, vget_high_s16(wv), HALF(x3), L); } +/* GEMM_STEP(lane, half): the step order above walks lanes 0..3 of the low + * half, then 0..3 of the high half, matching the weight rows i..i+7. */ +#define GEMM_STEP(L, HALF) GEMM_LANE(0, HALF, L) +DEFINE_GEMM(int8_t) +DEFINE_GEMM(int16_t) + +static void conv1x1_neon(const Op *o, const Tensor *ti, Tensor *to) { + gemm_int8_t(o, ti->q, o->cin, o->bias_folded, to->q, to->h * to->w); +} + +/* kxk (the stem) as im2col and the same multiply: each output pixel's k*k*cin + * inputs, zero point taken out, padding as 0, rows padded to a multiple of 8. */ +static void conv_kxk_neon(const Op *o, const Tensor *ti, Tensor *to, int16_t *col) { + const int k = o->k, st = o->stride, pad = o->pad, cin = o->cin, K = o->kpad; + /* One output row at a time, so the rows stay in L1 between being written + * and being read: the first layer's whole im2col is a megabyte. */ + for (int oy = 0; oy < to->h; oy++) { + int16_t *r = col; + for (int ox = 0; ox < to->w; ox++, r += K) { + int16_t *q = r; + for (int ky = 0; ky < k; ky++) { + const int iy = oy * st + ky - pad; + for (int kx = 0; kx < k; kx++) { + const int ix = ox * st + kx - pad; + if (iy < 0 || iy >= ti->h || ix < 0 || ix >= ti->w) { + for (int i = 0; i < cin; i++) + *q++ = 0; + continue; + } + const int8_t *ip = ti->q + ((size_t)iy * ti->w + ix) * cin; + for (int i = 0; i < cin; i++) + *q++ = (int16_t)(ip[i] - o->in_zp); + } + } + for (; q < r + K; q++) + *q = 0; + } + gemm_int16_t(o, col, K, o->bias, to->q + (size_t)oy * to->w * o->cout, to->w); + } +} + +/* Depthwise: eight channels at a time, the chunk's nine weight vectors held + * in registers across every pixel, and the accumulators across the taps. + * Interior pixels skip the bounds tests; the border takes them. */ +static void depthwise_neon(const Op *o, const Tensor *ti, Tensor *to) { + const int k = o->k, st = o->stride, pad = o->pad, c = o->cout; + const int16x8_t zp = vdupq_n_s16((int16_t)o->in_zp); + const size_t rs = (size_t)ti->w * c; /* input row stride */ + for (int ch = 0; ch < c; ch += 8) { + int16x8_t w[9]; + if (k == 3) + for (int t = 0; t < 9; t++) + w[t] = vld1q_s16(o->w16 + (size_t)t * c + ch); + const int32x4_t b0 = vld1q_s32(o->bias + ch), b1 = vld1q_s32(o->bias + ch + 4); + for (int oy = 0; oy < to->h; oy++) { + const int y0 = oy * st - pad; + for (int ox = 0; ox < to->w; ox++) { + const int x0 = ox * st - pad; + int32x4_t lo = b0, hi = b1; + if (k == 3 && y0 >= 0 && y0 + 2 < ti->h && x0 >= 0 && x0 + 2 < ti->w) { + const int8_t *p = ti->q + (size_t)y0 * rs + (size_t)x0 * c + ch; +#define DW_TAP(T, OFF) \ + { const int16x8_t a = vsubq_s16(vmovl_s8(vld1_s8(p + (OFF))), zp); \ + lo = vmlal_s16(lo, vget_low_s16(a), vget_low_s16(w[T])); \ + hi = vmlal_s16(hi, vget_high_s16(a), vget_high_s16(w[T])); } + DW_TAP(0, 0) DW_TAP(1, c) DW_TAP(2, 2 * c) + DW_TAP(3, rs) DW_TAP(4, rs + c) DW_TAP(5, rs + 2 * c) + DW_TAP(6, 2 * rs) DW_TAP(7, 2 * rs + c) DW_TAP(8, 2 * rs + 2 * c) +#undef DW_TAP + } else { + for (int ky = 0; ky < k; ky++) { + const int iy = y0 + ky; + if (iy < 0 || iy >= ti->h) + continue; + for (int kx = 0; kx < k; kx++) { + const int ix = x0 + kx; + if (ix < 0 || ix >= ti->w) + continue; + const int16x8_t a = vsubq_s16( + vmovl_s8(vld1_s8(ti->q + (size_t)iy * rs + (size_t)ix * c + ch)), zp); + const int16x8_t wv = vld1q_s16(o->w16 + (size_t)(ky * k + kx) * c + ch); + lo = vmlal_s16(lo, vget_low_s16(a), vget_low_s16(wv)); + hi = vmlal_s16(hi, vget_high_s16(a), vget_high_s16(wv)); + } + } + } + vst1_s8(to->q + ((size_t)oy * to->w + ox) * c + ch, requant8(o, ch, lo, hi)); + } + } + } +} +#endif + +/* 1x1: every output pixel is a row of cin bytes against cout weight rows. + * The zero point is folded into the bias at load. */ +static void conv1x1(const Op *o, const Tensor *ti, Tensor *to) { + const int hw = to->h * to->w, cin = o->cin, cout = o->cout; + for (int m = 0; m < hw; m++) { + const int8_t *a = ti->q + (size_t)m * cin; + for (int c = 0; c < cout; c++) { + const int32_t acc = o->bias_folded[c] + dot_i8(a, o->w + (size_t)c * cin, cin); + emit(o, acc, c, to->q ? to->q + (size_t)m * cout + c : NULL, + to->f ? to->f + (size_t)m * cout + c : NULL); + } + } +} + +/* kxk, all channels mixed: the stem. Taps outside the picture are padding, + * which contributes (zp - zp) * w = 0 and is skipped. */ +static void conv_kxk(const Op *o, const Tensor *ti, Tensor *to, int32_t *acc) { + const int k = o->k, st = o->stride, pad = o->pad, cin = o->cin, cout = o->cout; + for (int oy = 0; oy < to->h; oy++) + for (int ox = 0; ox < to->w; ox++) { + for (int c = 0; c < cout; c++) + acc[c] = o->bias[c]; + for (int ky = 0; ky < k; ky++) { + const int iy = oy * st + ky - pad; + if (iy < 0 || iy >= ti->h) + continue; + for (int kx = 0; kx < k; kx++) { + const int ix = ox * st + kx - pad; + if (ix < 0 || ix >= ti->w) + continue; + const int8_t *ip = ti->q + ((size_t)iy * ti->w + ix) * cin; + for (int c = 0; c < cout; c++) { + const int8_t *wp = o->w + (((size_t)c * k + ky) * k + kx) * cin; + int32_t s = 0; + for (int i = 0; i < cin; i++) + s += ((int32_t)ip[i] - o->in_zp) * wp[i]; + acc[c] += s; + } + } + } + const size_t base = ((size_t)oy * to->w + ox) * cout; + for (int c = 0; c < cout; c++) + emit(o, acc[c], c, to->q ? to->q + base + c : NULL, + to->f ? to->f + base + c : NULL); + } +} + +/* Depthwise kxk: channels innermost, so each tap is one contiguous run. */ +static void depthwise(const Op *o, const Tensor *ti, Tensor *to, int32_t *acc) { + const int k = o->k, st = o->stride, pad = o->pad, c = o->cout; + for (int oy = 0; oy < to->h; oy++) + for (int ox = 0; ox < to->w; ox++) { + for (int ch = 0; ch < c; ch++) + acc[ch] = o->bias[ch]; + for (int ky = 0; ky < k; ky++) { + const int iy = oy * st + ky - pad; + if (iy < 0 || iy >= ti->h) + continue; + for (int kx = 0; kx < k; kx++) { + const int ix = ox * st + kx - pad; + if (ix < 0 || ix >= ti->w) + continue; + const int8_t *ip = ti->q + ((size_t)iy * ti->w + ix) * c; + const int8_t *wp = o->w + ((size_t)ky * k + kx) * c; + int ch = 0; +#ifdef PDNET_NEON + const int16x8_t zp = vdupq_n_s16((int16_t)o->in_zp); + for (; ch + 8 <= c; ch += 8) { + const int16x8_t a = vsubq_s16(vmovl_s8(vld1_s8(ip + ch)), zp); + const int16x8_t w = vmovl_s8(vld1_s8(wp + ch)); + int32x4_t lo = vld1q_s32(acc + ch), hi = vld1q_s32(acc + ch + 4); + lo = vmlal_s16(lo, vget_low_s16(a), vget_low_s16(w)); + hi = vmlal_s16(hi, vget_high_s16(a), vget_high_s16(w)); + vst1q_s32(acc + ch, lo); + vst1q_s32(acc + ch + 4, hi); + } +#endif + for (; ch < c; ch++) + acc[ch] += ((int32_t)ip[ch] - o->in_zp) * wp[ch]; + } + } + const size_t base = ((size_t)oy * to->w + ox) * c; + for (int ch = 0; ch < c; ch++) + emit(o, acc[ch], ch, to->q ? to->q + base + ch : NULL, + to->f ? to->f + base + ch : NULL); + } +} + +static void resize2x(const Tensor *ti, Tensor *to) { + const int c = ti->c; + for (int y = 0; y < to->h; y++) + for (int x = 0; x < to->w; x++) + memcpy(to->q + ((size_t)y * to->w + x) * c, + ti->q + ((size_t)(y / 2) * ti->w + x / 2) * c, (size_t)c); +} + +static void add(const Op *o, const Tensor *a, const Tensor *b, Tensor *to) { + const size_t e = (size_t)to->h * to->w * to->c; + size_t i = 0; +#ifdef PDNET_NEON + /* The same three roundings, two lanes of 64 bits at a time. */ + const int64x2_t na = vdupq_n_s64(-(int64_t)(31 - o->s_a)); + const int64x2_t nb = vdupq_n_s64(-(int64_t)(31 - o->s_b)); + const int64x2_t no = vdupq_n_s64(-(int64_t)(31 - o->s_o)); + const int32x2_t ma = vdup_n_s32(o->m_a), mb = vdup_n_s32(o->m_b), mo = vdup_n_s32(o->m_o); + const int32x4_t za = vdupq_n_s32(o->zp_a), zb = vdupq_n_s32(o->zp_b), zo = vdupq_n_s32(o->zp_o); +#define ADD_RQ(x, m, nsh) \ + vcombine_s32(vqmovn_s64(vrshlq_s64(vmull_s32(vget_low_s32(x), m), nsh)), \ + vqmovn_s64(vrshlq_s64(vmull_s32(vget_high_s32(x), m), nsh))) + for (; i + 8 <= e; i += 8) { + const int16x8_t va = vmovl_s8(vld1_s8(a->q + i)), vb = vmovl_s8(vld1_s8(b->q + i)); + int32x4_t out[2]; + for (int h = 0; h < 2; h++) { + const int16x4_t ha = h ? vget_high_s16(va) : vget_low_s16(va); + const int16x4_t hb = h ? vget_high_s16(vb) : vget_low_s16(vb); + const int32x4_t xa = vshlq_n_s32(vsubq_s32(vmovl_s16(ha), za), 20); + const int32x4_t xb = vshlq_n_s32(vsubq_s32(vmovl_s16(hb), zb), 20); + const int32x4_t sum = vaddq_s32(ADD_RQ(xa, ma, na), ADD_RQ(xb, mb, nb)); + out[h] = vaddq_s32(ADD_RQ(sum, mo, no), zo); + } + vst1_s8(to->q + i, vqmovn_s16(vcombine_s16(vqmovn_s32(out[0]), vqmovn_s32(out[1])))); + } +#undef ADD_RQ +#endif + for (; i < e; i++) { + const int32_t ra = requant(((int32_t)a->q[i] - o->zp_a) * (1 << 20), o->m_a, o->s_a); + const int32_t rb = requant(((int32_t)b->q[i] - o->zp_b) * (1 << 20), o->m_b, o->s_b); + to->q[i] = clamp8(requant(ra + rb, o->m_o, o->s_o) + o->zp_o); + } +} + +/* ---- the frame in, the people out ------------------------------------------- */ + +static float sigmoidf(float x) { return 1.0f / (1.0f + expf(-x)); } + +static float iou(const PdBox *a, const PdBox *b) { + const int x0 = a->x > b->x ? a->x : b->x, y0 = a->y > b->y ? a->y : b->y; + const int x1 = a->x + a->w < b->x + b->w ? a->x + a->w : b->x + b->w; + const int y1 = a->y + a->h < b->y + b->h ? a->y + a->h : b->y + b->h; + if (x1 <= x0 || y1 <= y0) + return 0; + const float i = (float)(x1 - x0) * (y1 - y0); + return i / ((float)a->w * a->h + (float)b->w * b->h - i); +} + +static int by_score(const void *a, const void *b) { + return ((const PdBox *)b)->score - ((const PdBox *)a)->score; +} + +int pdnet_detect( + PdNet *n, const uint8_t *y, int ystride, const uint8_t *vu, int vustride, + PdBox *out, int max, int score_min) { + if (!n || !y || !vu || !out || max <= 0) + return -1; + int8_t *in = n->t[0].q; + for (int j = 0; j < n->in_h; j++) { + const uint8_t *yr = y + (size_t)j * ystride, *cr = vu + (size_t)(j / 2) * vustride; + int8_t *p = in + (size_t)j * n->in_w * 3; + for (int i = 0; i < n->in_w; i++) { + const int ci = i & ~1; + p[3 * i + 0] = n->lut_y[yr[i]]; + p[3 * i + 1] = n->lut_c[cr[ci + 1]]; /* U; NV21 puts V first */ + p[3 * i + 2] = n->lut_c[cr[ci]]; /* V */ + } + } +#ifdef PDNET_PROFILE + struct timespec t0, t1; +#endif + for (int j = 0; j < n->n_ops; j++) { + const Op *o = &n->op[j]; + Tensor *to = &n->t[o->out]; +#ifdef PDNET_PROFILE + clock_gettime(CLOCK_MONOTONIC, &t0); +#endif + if (o->type == OP_CONV) { + const Tensor *ti = &n->t[o->in]; +#ifdef PDNET_NEON + /* Eight output channels a time; the heads (1 and 4 channels, and + * float) take the plain path, which costs them nothing. */ + if (!o->float_out && o->cout % 8 == 0 && + (o->dw || o->k > 1 || o->cin % 8 == 0)) { + if (o->dw) + depthwise_neon(o, ti, to); + else if (o->k == 1) + conv1x1_neon(o, ti, to); + else + conv_kxk_neon(o, ti, to, n->col); + } else +#endif + if (o->dw) + depthwise(o, ti, to, n->acc); + else if (o->k == 1) + conv1x1(o, ti, to); + else + conv_kxk(o, ti, to, n->acc); + } else if (o->type == OP_RESIZE2X) { + resize2x(&n->t[o->in], to); + } else { + add(o, &n->t[o->in], &n->t[o->in2], to); + } +#ifdef PDNET_PROFILE + clock_gettime(CLOCK_MONOTONIC, &t1); + pdnet_prof_us[j] += (t1.tv_sec - t0.tv_sec) * 1e6 + (t1.tv_nsec - t0.tv_nsec) / 1e3; +#endif + } + + /* Decode as the TVM shim did: an anchor-free head, a box per cell whose + * score clears the floor, distances in stride units, then NMS. */ + static PdBox cand[4096]; + int nc = 0; + for (int l = 0; l < n->n_lv; l++) { + const Tensor *tc = &n->t[n->lv[l].cls], *tr = &n->t[n->lv[l].reg]; + const int st = n->lv[l].stride; + for (int cy = 0; cy < tc->h; cy++) + for (int cx = 0; cx < tc->w; cx++) { + const int i = cy * tc->w + cx; + const int s = (int)(sigmoidf(tc->f[i]) * 100.0f + 0.5f); + if (s < score_min || nc == (int)(sizeof(cand) / sizeof(cand[0]))) + continue; + const float *r = tr->f + (size_t)i * 4; + const float ccx = (cx + 0.5f) * st, ccy = (cy + 0.5f) * st; + const int x0 = (int)(ccx - fmaxf(0.0f, r[0]) * st); + const int y0 = (int)(ccy - fmaxf(0.0f, r[1]) * st); + const int x1 = (int)(ccx + fmaxf(0.0f, r[2]) * st); + const int y1 = (int)(ccy + fmaxf(0.0f, r[3]) * st); + if (x1 <= x0 || y1 <= y0) + continue; + cand[nc++] = (PdBox){x0, y0, x1 - x0, y1 - y0, s}; + } + } + qsort(cand, (size_t)nc, sizeof(cand[0]), by_score); + int kept = 0; + for (int i = 0; i < nc && kept < max; i++) { + bool dup = false; + for (int j = 0; j < kept && !dup; j++) + dup = iou(&cand[i], &out[j]) > NMS_IOU; + if (!dup) + out[kept++] = cand[i]; + } + return kept; +} + +const float *pdnet_head(const PdNet *n, int level, int reg, int *h, int *w) { + if (!n || level < 0 || level >= n->n_lv) + return NULL; + const Tensor *t = &n->t[reg ? n->lv[level].reg : n->lv[level].cls]; + *h = t->h; + *w = t->w; + return t->f; +} + +#ifdef PDNET_PROFILE +void pdnet_profile_dump(const PdNet *n, int runs) { + double total = 0; + for (int j = 0; j < n->n_ops; j++) + total += pdnet_prof_us[j]; + for (int j = 0; j < n->n_ops; j++) { + const Op *o = &n->op[j]; + const Tensor *to = &n->t[o->out]; + const char *kind = o->type == OP_RESIZE2X ? "resize" : o->type == OP_ADD ? "add" + : o->dw ? "dw" : o->k == 1 ? "1x1" : "kxk"; + long macs = 0; + if (o->type == OP_CONV) + macs = (long)to->h * to->w * o->cout * o->k * o->k * (o->dw ? 1 : o->cin); + fprintf(stderr, "op %2d %-6s k%d s%d %3d->%3d out %3dx%3d %7.2f ms %5.1f%% %6.1f MMAC/s\n", + j, kind, o->k, o->stride, o->cin, o->cout, to->h, to->w, + pdnet_prof_us[j] / runs / 1e3, 100 * pdnet_prof_us[j] / total, + macs ? macs / (pdnet_prof_us[j] / runs) : 0.0); + } + fprintf(stderr, "total %.2f ms per run\n", total / runs / 1e3); +} +#endif diff --git a/engine/pdnet.h b/engine/pdnet.h new file mode 100644 index 0000000..7840fb2 --- /dev/null +++ b/engine/pdnet.h @@ -0,0 +1,34 @@ +#pragma once + +/* pdnet: the CPU person detector's int8 network, run directly. + * + * Loads the blob export_pdnet.py writes from the quantized model and runs it + * on an NV21 frame of the network's input size: int8 activations, per-channel + * int8 weights, int32 accumulators, one rounding per requantization, and the + * detection heads taken to float straight from their accumulators. */ + +#include +#include + +typedef struct PdNet PdNet; + +typedef struct { + int x, y, w, h; + int score; /* 0..100 */ +} PdBox; + +PdNet *pdnet_open(const char *path); +void pdnet_close(PdNet *n); +int pdnet_width(const PdNet *n); +int pdnet_height(const PdNet *n); + +/* Detect people on one NV21 frame of exactly pdnet_width() x pdnet_height(): + * boxes in input pixels, scores 0..100, at most `max`, best first, after + * suppression. Returns the count, or -1. */ +int pdnet_detect( + PdNet *n, const uint8_t *y, int ystride, const uint8_t *vu, int vustride, + PdBox *out, int max, int score_min); + +/* The raw head outputs of the last pdnet_detect(), for checking against the + * reference: level 0..2, cls (1 channel) or reg (4), HWC float. */ +const float *pdnet_head(const PdNet *n, int level, int reg, int *h, int *w); diff --git a/models/person/model.onnx b/models/person/model.onnx new file mode 100644 index 0000000..320dcc4 Binary files /dev/null and b/models/person/model.onnx differ diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..26e7bca --- /dev/null +++ b/requirements.txt @@ -0,0 +1,3 @@ +numpy==2.4.6 +onnx==1.23.2 +opencv-python-headless==5.0.0.93 diff --git a/scripts/build.sh b/scripts/build.sh new file mode 100755 index 0000000..7b777ad --- /dev/null +++ b/scripts/build.sh @@ -0,0 +1,59 @@ +#!/bin/bash +# build.sh [out-dir] +# +# The whole build, the same here and in .github/workflows/release.yml: +# 1. fetch COCO val2017 and its annotations, hash-checked (coco.txt) +# 2. export the int8 model to the blob majestic's engine loads (person.pdnet) +# 3. build the engine's tools for this host, and for ARMv7 with NEON (the +# camera's path), run the ARM ones under qemu-arm +# 4. check both against the reference arithmetic on 40 crops, bit for bit: +# nothing below is measured on an engine that does not compute the model +# 5. measure the presets on the labelled crops (points.py) and write +# person.json, points.md and SHA256SUMS +# +# Needs Python 3.11 with requirements.txt, gcc, gcc-arm-linux-gnueabihf and +# qemu-user-static. +set -euo pipefail +src=$(cd "$(dirname "$0")/.." && pwd) +out=$(readlink -m "${1:-$src/out}") +cache=$src/.cache +mkdir -p "$out" "$cache" + +bash "$src/scripts/fetch.sh" "$src/coco.txt" "$cache" +# Extracted fresh every build, from the archives just hash-checked: files left +# in a reused output directory are not what was pinned. +coco=$out/coco +rm -rf "$coco" "$out/eval" +mkdir -p "$coco" +unzip -qo -d "$coco" "$cache/val2017.zip" +unzip -qo -d "$coco" "$cache/annotations_trainval2017.zip" annotations/instances_val2017.json + +python3 "$src/tools/export_pdnet.py" "$src/models/person/model.onnx" "$out/person.pdnet" + +cflags=(-O2 -std=c11 -Wall -Wextra -I"$src/engine") +gcc "${cflags[@]}" -o "$out/dets" "$src/tools/dets.c" "$src/engine/pdnet.c" -lm +gcc "${cflags[@]}" -o "$out/dump_heads" "$src/tools/dump_heads.c" "$src/engine/pdnet.c" -lm +# _FILE_OFFSET_BITS: a 32-bit readdir cannot return the 64-bit directory +# offsets some host filesystems hand out, and finds nothing. +arm=(arm-linux-gnueabihf-gcc "${cflags[@]}" -D_FILE_OFFSET_BITS=64 -static -march=armv7-a -mfpu=neon -mfloat-abi=hard) +"${arm[@]}" -o "$out/dump_heads.arm" "$src/tools/dump_heads.c" "$src/engine/pdnet.c" -lm +"${arm[@]}" -o "$out/dets.arm" "$src/tools/dets.c" "$src/engine/pdnet.c" -lm + +python3 "$src/tools/crops.py" "$coco" "$out/eval" +mapfile -t sample < <(ls "$out"/eval/crops/*.nv21 | sort | awk 'NR % 75 == 1' | head -40) +for run in "$out/dump_heads" "qemu-arm-static $out/dump_heads.arm"; do + rm -f "$out"/eval/crops/*.heads + $run "$out/person.pdnet" "${sample[@]}" + echo "== exactness: ${run%% *}" + python3 "$src/tools/check_exact.py" "$src/models/person/model.onnx" "${sample[@]}" +done +rm -f "$out"/eval/crops/*.heads + +"$out/dets" "$out/person.pdnet" "$out/eval/crops" > "$out/dets.txt" +# The camera's path must find the same boxes, not just the same heads. +qemu-arm-static "$out/dets.arm" "$out/person.pdnet" "$out/eval/crops" > "$out/dets.arm.txt" +cmp "$out/dets.txt" "$out/dets.arm.txt" || { echo "ARM and host detections differ"; exit 1; } + +python3 "$src/tools/points.py" "$out/eval/labels.json" "$out/dets.txt" \ + "$out/person.json" "$out/points.md" +cd "$out" && sha256sum person.pdnet person.json > SHA256SUMS && cat SHA256SUMS diff --git a/scripts/fetch.sh b/scripts/fetch.sh new file mode 100755 index 0000000..8a4413b --- /dev/null +++ b/scripts/fetch.sh @@ -0,0 +1,22 @@ +#!/bin/bash +# fetch.sh : every "file sha256 url" line of , into +# , hash-checked. A cached file whose hash does not match is +# fetched again, and a download lands in the cache only once its hash +# matches, so an interrupted one can never be mistaken for the file. +set -euo pipefail +mkdir -p "$2" +grep -v '^#' "$1" | while read -r file sha url; do + [ -n "$file" ] || continue + dst="$2/$file" + if [ -f "$dst" ] && echo "$sha $dst" | sha256sum -c --status 2>/dev/null; then + continue + fi + tmp="$dst.part" + if ! curl -sfL --retry 3 -o "$tmp" "$url" || + ! echo "$sha $tmp" | sha256sum -c --status; then + rm -f "$tmp" + echo "fetch.sh: $file from $url: download failed or sha256 mismatch" >&2 + exit 1 + fi + mv "$tmp" "$dst" +done diff --git a/tools/check_exact.py b/tools/check_exact.py new file mode 100644 index 0000000..d6d637d --- /dev/null +++ b/tools/check_exact.py @@ -0,0 +1,46 @@ +# check_exact.py ... +# +# The engine must compute what the quantized model defines. For each crop, +# .heads (written by dump_heads, from pdnet.c) must equal pdq.py's +# reference arithmetic exactly, every value of all six head maps. Exits +# non-zero on any difference, so a release stops before it measures or +# publishes anything. +import sys + +import numpy as np + +sys.path.insert(0, __file__.rsplit("/", 1)[0]) +import pdq # noqa: E402 + +W, H = 320, 192 + + +def chw_from_nv21(y, vu): + """The float model's input, as majestic's preprocessing defines it.""" + ky = np.float32(255 / 219 / 128) + kc = np.float32(255 / 224 / 128) + Y = (y.astype(np.float32) - np.float32(16)) * ky - np.float32(1) + V = (vu[:, 0::2].astype(np.float32) - np.float32(128)) * kc + U = (vu[:, 1::2].astype(np.float32) - np.float32(128)) * kc + up = lambda p: np.repeat(np.repeat(p, 2, 0), 2, 1) + return np.stack([Y, up(U), up(V)]).astype(np.float32) + + +g = pdq.Graph(sys.argv[1]) +crops = sys.argv[2:] +if not crops: + sys.exit("check_exact: no crops") +bad = 0 +for f in crops: + raw = np.fromfile(f, np.uint8) + ref = pdq.run(g, chw_from_nv21(raw[:W * H].reshape(H, W), raw[W * H:].reshape(H // 2, W)), + head="float") + got = np.fromfile(f + ".heads", np.float32) + want = np.concatenate([r.transpose(1, 2, 0).reshape(-1) for r in ref]).astype(np.float32) + if got.shape != want.shape or not np.array_equal(got, want): + n = int(np.sum(got != want)) if got.shape == want.shape else -1 + print(f"check_exact: {f}: {n} values differ") + bad += 1 +if bad: + sys.exit(f"check_exact: {bad} of {len(crops)} crops differ from the reference") +print(f"check_exact: all {len(crops)} crops identical to the reference") diff --git a/tools/crops.py b/tools/crops.py new file mode 100644 index 0000000..e5fbaef --- /dev/null +++ b/tools/crops.py @@ -0,0 +1,107 @@ +# crops.py [n_person] [n_background] +# +# The labelled set the person model's presets are measured on, in the +# geometry the camera hands the model: 320x192 crops around a person at the +# scale the detector was trained for (the person ~47% of the crop's height, as +# majestic places a crop), off-centre the way a motion region is, and crops of +# person-free COCO val2017 images for the false alarms. Each crop is NV21, +# BT.709 limited range, as the VPSS delivers it; labels.json lists the people +# inside each in crop pixels, marking crowds and people less than half inside +# the crop to be ignored. Seeded: the same inputs give the same crops. +import json +import os +import random +import sys + +import cv2 +import numpy as np + +W, H = 320, 192 +FILL = 0.47 + + +def nv21(rgb): + r, g, b = rgb[..., 0], rgb[..., 1], rgb[..., 2] + y = 16 + (0.2126 * r + 0.7152 * g + 0.0722 * b) * 219 / 255 + cb = 128 + (-0.1146 * r - 0.3854 * g + 0.5 * b) * 224 / 255 + cr = 128 + (0.5 * r - 0.4542 * g - 0.0458 * b) * 224 / 255 + sub = lambda c: c.reshape(H // 2, 2, W // 2, 2).mean((1, 3)) + vu = np.stack([sub(cr), sub(cb)], -1).reshape(H // 2, W) + return (np.clip(y, 0, 255).round().astype(np.uint8), + np.clip(vu, 0, 255).round().astype(np.uint8)) + + +def main(coco, out, n_person=2000, n_bg=1000): + n_person, n_bg = int(n_person), int(n_bg) + rnd = random.Random(20261011) + ann = json.load(open(f"{coco}/annotations/instances_val2017.json")) + imgs = {i["id"]: i for i in ann["images"]} + people = {} + for a in ann["annotations"]: + if a["category_id"] == 1: + people.setdefault(a["image_id"], []).append(a) + # A directory holding crops from another run would mix them in. + if os.path.isdir(f"{out}/crops") and os.listdir(f"{out}/crops"): + sys.exit(f"crops.py: {out}/crops is not empty") + os.makedirs(f"{out}/crops", exist_ok=True) + labels = {} + + def emit(name, img, x0, y0, cw, ch, anns): + x0, y0, cw, ch = (int(round(v)) for v in (x0, y0, cw, ch)) + crop = img[y0:y0 + ch, x0:x0 + cw] + rgb = cv2.resize(cv2.cvtColor(crop, cv2.COLOR_BGR2RGB), (W, H), + interpolation=cv2.INTER_AREA).astype(np.float32) + y, vu = nv21(rgb) + with open(f"{out}/crops/{name}.nv21", "wb") as f: + f.write(y.tobytes() + vu.tobytes()) + sx, sy = W / cw, H / ch + boxes = [] + for a in anns: + bx, by, bw, bh = a["bbox"] + ix0, iy0 = max(bx, x0), max(by, y0) + ix1, iy1 = min(bx + bw, x0 + cw), min(by + bh, y0 + ch) + if ix1 <= ix0 or iy1 <= iy0: + continue + visible = (ix1 - ix0) * (iy1 - iy0) / max(bw * bh, 1e-9) + boxes.append({"box": [(ix0 - x0) * sx, (iy0 - y0) * sy, + (ix1 - ix0) * sx, (iy1 - iy0) * sy], + "ignore": bool(a["iscrowd"]) or visible < 0.5}) + labels[name] = boxes + + cands = [(iid, a) for iid, anns in sorted(people.items()) for a in anns + if not a["iscrowd"] and a["bbox"][3] >= 24] + rnd.shuffle(cands) + made = 0 + for iid, a in cands: + if made >= n_person: + break + im = imgs[iid] + bx, by, bw, bh = a["bbox"] + ch = bh / FILL + cw = ch * W / H + if cw > im["width"] or ch > im["height"]: + continue + cx = bx + bw / 2 + rnd.uniform(-cw / 3, cw / 3) + cy = by + bh / 2 + rnd.uniform(-ch / 6, ch / 6) + x0 = min(max(cx - cw / 2, 0), im["width"] - cw) + y0 = min(max(cy - ch / 2, 0), im["height"] - ch) + img = cv2.imread(f"{coco}/val2017/{im['file_name']}") + emit(f"p{made:05d}_{iid}", img, x0, y0, cw, ch, people[iid]) + made += 1 + empty = sorted(i for i in imgs if i not in people) + rnd.shuffle(empty) + for k, iid in enumerate(empty[:n_bg]): + im = imgs[iid] + ch = rnd.uniform(0.35, 1.0) * min(im["height"], im["width"] * H / W) + cw = ch * W / H + x0 = rnd.uniform(0, im["width"] - cw) + y0 = rnd.uniform(0, im["height"] - ch) + img = cv2.imread(f"{coco}/val2017/{im['file_name']}") + emit(f"b{k:05d}_{iid}", img, x0, y0, cw, ch, []) + with open(f"{out}/labels.json", "w") as f: + json.dump(labels, f, sort_keys=True) + print(f"crops.py: {made} person crops, {min(n_bg, len(empty))} person-free") + + +if __name__ == "__main__": + main(*sys.argv[1:]) diff --git a/tools/dets.c b/tools/dets.c new file mode 100644 index 0000000..f53601a --- /dev/null +++ b/tools/dets.c @@ -0,0 +1,66 @@ +#define _POSIX_C_SOURCE 200809L +/* dets : every box the engine finds on every crop, as + * " x y w h score" (input pixels, score 0..100), at a score floor of + * 1, so points.py can sweep every threshold. The crops in name order. */ +#include +#include +#include +#include + +#include "pdnet.h" + +static int cmp(const void *a, const void *b) { + return strcmp(*(char *const *)a, *(char *const *)b); +} + +int main(int argc, char **argv) { + if (argc != 3) { + fprintf(stderr, "usage: %s \n", argv[0]); + return 2; + } + PdNet *n = pdnet_open(argv[1]); + if (!n) { + fprintf(stderr, "dets: cannot open %s\n", argv[1]); + return 1; + } + const int w = pdnet_width(n), h = pdnet_height(n); + const size_t frame = (size_t)w * h * 3 / 2; + DIR *d = opendir(argv[2]); + if (!d) + return 1; + char **names = NULL; + size_t nn = 0, cap = 0; + struct dirent *e; + while ((e = readdir(d))) { + const size_t l = strlen(e->d_name); + if (l <= 5 || strcmp(e->d_name + l - 5, ".nv21")) + continue; + if (nn == cap && !(names = realloc(names, (cap = cap ? 2 * cap : 1024) * sizeof(*names)))) + return 1; + names[nn++] = strdup(e->d_name); + } + closedir(d); + qsort(names, nn, sizeof(*names), cmp); + unsigned char *buf = malloc(frame); + PdBox box[256]; + for (size_t i = 0; i < nn; i++) { + char path[4096]; + snprintf(path, sizeof path, "%s/%s", argv[2], names[i]); + FILE *f = fopen(path, "rb"); + if (!f || fread(buf, 1, frame, f) != frame) { + fprintf(stderr, "dets: cannot read %s\n", path); + return 1; + } + fclose(f); + const int k = pdnet_detect(n, buf, w, buf + (size_t)w * h, w, box, 256, 1); + if (k < 0) + return 1; + names[i][strlen(names[i]) - 5] = 0; + for (int j = 0; j < k; j++) + printf("%s %d %d %d %d %d\n", names[i], box[j].x, box[j].y, box[j].w, + box[j].h, box[j].score); + } + fprintf(stderr, "dets: %zu crops\n", nn); + pdnet_close(n); + return 0; +} diff --git a/tools/dump_heads.c b/tools/dump_heads.c new file mode 100644 index 0000000..43250b8 --- /dev/null +++ b/tools/dump_heads.c @@ -0,0 +1,28 @@ +/* dump_heads ...: each crop's six head maps, raw float32, + * to .heads, for comparison with pdq.py. */ +#include +#include +#include "pdnet.h" +int main(int argc, char **argv) { + PdNet *n = pdnet_open(argv[1]); + if (!n) { fprintf(stderr, "cannot open %s\n", argv[1]); return 1; } + const int w = pdnet_width(n), h = pdnet_height(n); + unsigned char *buf = malloc((size_t)w * h * 3 / 2); + PdBox box[64]; + for (int a = 2; a < argc; a++) { + FILE *f = fopen(argv[a], "rb"); + if (!f || fread(buf, 1, (size_t)w * h * 3 / 2, f) != (size_t)w * h * 3 / 2) return 1; + fclose(f); + pdnet_detect(n, buf, w, buf + (size_t)w * h, w, box, 64, 101); + char out[1024]; snprintf(out, sizeof out, "%s.heads", argv[a]); + FILE *o = fopen(out, "wb"); + for (int l = 0; l < 3; l++) + for (int r = 0; r < 2; r++) { + int hh, ww; const float *p = pdnet_head(n, l, r, &hh, &ww); + fwrite(p, sizeof(float), (size_t)hh * ww * (r ? 4 : 1), o); + } + fclose(o); + } + pdnet_close(n); + return 0; +} diff --git a/tools/export_pdnet.py b/tools/export_pdnet.py new file mode 100644 index 0000000..638491e --- /dev/null +++ b/tools/export_pdnet.py @@ -0,0 +1,134 @@ +# export_pdnet.py +# +# The person detector's QDQ graph as the blob pdnet.c runs: the same walk, +# constants and arithmetic as pdq.run(head="float"), which is the reference +# the C engine is checked against bit for bit. +# +# Layout, little-endian, every section 4-byte aligned: +# "MAJPDNT1" u32 version=1 +# u32 in_w in_h f32 in_scale i32 in_zp +# u32 n_tensors, then per tensor u32 h w c +# u32 n_ops, then per op u32 type, and: +# CONV(1): u32 in out k stride pad depthwise cin cout float_out +# i32 in_zp out_zp +# float_out ? f32 scale[cout] : i32 (mult, shift)[cout] +# i32 bias[cout] +# i8 w[...] then pad to 4 -- [cout][cin] for 1x1, +# [cout][k][k][cin] for kxk, [k][k][c] for depthwise +# RESIZE2X(2): u32 in out (nearest, same quantization) +# ADD(3): u32 a b out i32 zp_a zp_b zp_out +# i32 mult_a shift_a mult_b shift_b mult_out shift_out +# u32 n_levels, then per level u32 stride cls_tensor reg_tensor +import struct +import sys + +import numpy as np +import onnx + +sys.path.insert(0, __file__.rsplit("/", 1)[0]) +import pdq # noqa: E402 + +CONV, RESIZE2X, ADD = 1, 2, 3 +ADD_SHIFT = 20 + + +def main(src, dst): + g = pdq.Graph(src) + shapes = {} # tensor name -> (h, w, c) + ids = {} + out = bytearray() + + def tid(name, shape): + if name not in ids: + ids[name] = len(ids) + shapes[name] = shape + return ids[name] + + in_w, in_h = 320, 192 + qi = next(n for n in g.nodes if n.op_type == "QuantizeLinear" and n.input[0] == g.input) + in_scale, in_zp = g.qparams[qi.output[0]] + tid(qi.output[0], (in_h, in_w, 3)) + ops = [] + final = {} + for n in g.nodes: + if n.op_type == "Conv": + src_t, s_in, zp_in = g.dq_source(n.input[0]) + wq, s_w, _ = g.dq_source(n.input[1]) + w = g.init[wq].astype(np.int8) + bias = np.zeros(w.shape[0], np.int32) + if len(n.input) > 2: + bq, _, _ = g.dq_source(n.input[2]) + bias = g.init[bq].astype(np.int32) + a = {x.name: onnx.helper.get_attribute_value(x) for x in n.attribute} + k, st, grp = a["kernel_shape"][0], a["strides"][0], a.get("group", 1) + pad = a.get("pads", [0, 0, 0, 0])[0] + ih, iw, ic = shapes[src_t] + oh, ow = (ih + 2 * pad - k) // st + 1, (iw + 2 * pad - k) // st + 1 + co = w.shape[0] + dw = grp > 1 + real = float(s_in) * s_w.astype(np.float64) + q = next(c for c in g.nodes if n.output[0] in c.input) + dq = next(c for c in g.nodes if q.output[0] in c.input) + is_final = any(o in g.outputs for o in dq.output) + o = tid(q.output[0], (oh, ow, co)) + if is_final: + final[dq.output[0]] = o + zp_out = 0 + else: + s_out, zp_out = g.qparams[q.output[0]] + rec = struct.pack("<10I2i", CONV, ids[src_t], o, k, st, pad, int(dw), ic, co, + int(is_final), int(zp_in), int(zp_out)) + if is_final: + rec += real.astype(np.float32).tobytes() + else: + ms = [pdq.quantize_multiplier(r / s_out) for r in real] + rec += np.array(ms, np.int32).reshape(-1).tobytes() + rec += bias.tobytes() + if dw: + wb = w[:, 0].transpose(1, 2, 0) # k, k, c + elif k == 1: + wb = w[:, :, 0, 0] # co, ci + else: + wb = w.transpose(0, 2, 3, 1) # co, k, k, ci + wb = np.ascontiguousarray(wb).tobytes() + rec += wb + bytes(-len(wb) % 4) + ops.append(rec) + elif n.op_type == "Resize": + src_t, _, _ = g.dq_source(n.input[0]) + q = next(c for c in g.nodes if n.output[0] in c.input) + h, w_, c = shapes[src_t] + ops.append(struct.pack("<3I", RESIZE2X, ids[src_t], tid(q.output[0], (2 * h, 2 * w_, c)))) + elif n.op_type == "Add": + (ta, s_a, zp_a), (tb, s_b, zp_b) = (g.dq_source(i) for i in n.input) + q = next(c for c in g.nodes if n.output[0] in c.input) + s_out, zp_out = g.qparams[q.output[0]] + twice = 2 * max(float(s_a), float(s_b)) + ma = pdq.quantize_multiplier(float(s_a) / twice) + mb = pdq.quantize_multiplier(float(s_b) / twice) + mo = pdq.quantize_multiplier(twice / ((1 << ADD_SHIFT) * s_out)) + assert shapes[ta] == shapes[tb] + o = tid(q.output[0], shapes[ta]) + ops.append(struct.pack("<4I9i", ADD, ids[ta], ids[tb], o, int(zp_a), int(zp_b), + int(zp_out), *ma, *mb, *mo)) + out += b"MAJPDNT1" + struct.pack(" 0 as (int32 mult, shift): m == mult * 2^(shift - 31).""" + if m == 0: + return 0, 0 + mant, exp = np.frexp(m) + q = int(round(mant * (1 << 31))) + if q == 1 << 31: + q //= 2 + exp += 1 + return q, int(exp) + + +def requant(acc, mult, shift): + """acc * mult * 2^(shift - 31), rounded once (half up). + + Not TFLite's MultiplyByQuantizedMultiplier, which rounds twice -- once in + the doubling high multiply, again in the power-of-two divide -- and so + lands one off the exact result for a quarter of negative values at some + shifts. One 64-bit product and one rounding shift agree with exact + rounding except at exact halves; ARMv7 NEON has the 32x32->64 multiply + this needs (vmull_s32).""" + n = (31 - np.asarray(shift, np.int64)) + prod = acc.astype(np.int64) * np.asarray(mult, np.int64) + return (prod + (np.int64(1) << (n - 1))) >> n + + +class Graph: + def __init__(self, path): + m = onnx.load(path) + self.init = {t.name: nh.to_array(t) for t in m.graph.initializer} + for n in m.graph.node: + if n.op_type == "Constant": + self.init[n.output[0]] = nh.to_array(n.attribute[0].t) + self.nodes = [n for n in m.graph.node if n.op_type != "Constant"] + self.prod = {o: n for n in self.nodes for o in n.output} + self.outputs = [o.name for o in m.graph.output] + self.input = m.graph.input[0].name + self.qparams = {} # quantized tensor name -> (scale, zp) + for n in self.nodes: + if n.op_type == "QuantizeLinear": + self.qparams[n.output[0]] = (float(self.init[n.input[1]]), + int(self.init[n.input[2]])) + + def dq_source(self, name): + """A DequantizeLinear's quantized input and its (scale, zp).""" + n = self.prod[name] + assert n.op_type == "DequantizeLinear", (name, n.op_type) + src = n.input[0] + s = self.init[n.input[1]] + zp = self.init[n.input[2]] if len(n.input) > 2 else np.zeros_like(s, np.int32) + return src, s, zp + + +def run(g, x_chw, head="float", keep=None, rq="fixed"): + """x_chw: float32 [3, H, W] as the float model takes it.""" + t = {} # tensor name -> int8 NHWC (H, W, C) for quantized activations + qi = next(n for n in g.nodes if n.op_type == "QuantizeLinear" and n.input[0] == g.input) + s, zp = g.qparams[qi.output[0]] + t[qi.output[0]] = np.clip(np.round(x_chw.transpose(1, 2, 0) / s) + zp, -128, 127).astype(np.int32) + outs = {} + for n in g.nodes: + if n.op_type == "Conv": + src, s_in, zp_in = g.dq_source(n.input[0]) + wq, s_w, zp_w = g.dq_source(n.input[1]) + w = g.init[wq].astype(np.int32) # [Co, Ci/g, k, k] + assert not np.any(zp_w) + bias = np.zeros(w.shape[0], np.int64) + if len(n.input) > 2: + bq, s_b, _ = g.dq_source(n.input[2]) + bias = g.init[bq].astype(np.int64) + a = {x.name: onnx.helper.get_attribute_value(x) for x in n.attribute} + k, st, grp = a["kernel_shape"][0], a["strides"][0], a.get("group", 1) + pads = a.get("pads", [0, 0, 0, 0]) + x = t[src] - int(zp_in) + H, W, C = x.shape + xp = np.zeros((H + pads[0] + pads[2], W + pads[1] + pads[3], C), np.int32) + xp[pads[0]:pads[0] + H, pads[1]:pads[1] + W] = x + Ho = (xp.shape[0] - k) // st + 1 + Wo = (xp.shape[1] - k) // st + 1 + Co = w.shape[0] + acc = np.zeros((Ho, Wo, Co), np.int64) + cpg_in, cpg_out = C // grp, Co // grp + for ky in range(k): + for kx in range(k): + patch = xp[ky:ky + st * Ho:st, kx:kx + st * Wo:st] # Ho, Wo, C + if grp == 1: + acc += patch.astype(np.int64) @ w[:, :, ky, kx].T.astype(np.int64) + elif cpg_in == 1 and cpg_out == 1: + acc += patch.astype(np.int64) * w[:, 0, ky, kx].astype(np.int64) + else: + raise NotImplementedError("grouped conv") + acc += bias + real_scale = float(s_in) * s_w.astype(np.float64) # per output channel + q = next(c for c in g.nodes if n.output[0] in c.input) + final = any(o in g.outputs for o in + next(c for c in g.nodes if q.output[0] in c.input).output) + if final and head == "float": + # the scale as the blob stores it, float32, so the C engine + # can match this exactly + outs[q.output[0]] = acc * real_scale.astype(np.float32).astype(np.float64) + continue + s_out, zp_out = g.qparams[q.output[0]] + if rq == "float": # onnxruntime's way: float, round half to even + y = np.round((acc * real_scale).astype(np.float32) / np.float32(s_out)).astype(np.int64) + zp_out + else: + mults = [quantize_multiplier(r / s_out) for r in real_scale] + mult = np.array([m for m, _ in mults], np.int64) + shift = np.array([e for _, e in mults], np.int64) + y = requant(acc, mult, shift) + zp_out + t[q.output[0]] = np.clip(y, -128, 127).astype(np.int32) + elif n.op_type == "Resize": + src, s_in, zp_in = g.dq_source(n.input[0]) + q = next(c for c in g.nodes if n.output[0] in c.input) + s_out, zp_out = g.qparams[q.output[0]] + assert abs(float(s_in) - s_out) < 1e-9 and int(zp_in) == zp_out + t[q.output[0]] = t[src].repeat(2, 0).repeat(2, 1) # nearest, x2 + elif n.op_type == "Add": + (sa, s_a, zp_a), (sb, s_b, zp_b) = (g.dq_source(i) for i in n.input) + q = next(c for c in g.nodes if n.output[0] in c.input) + s_out, zp_out = g.qparams[q.output[0]] + # TFLite's quantized add: both sides to a shared 2^20 grid + ls = 20 + twice = 2 * max(float(s_a), float(s_b)) + ma, ea = quantize_multiplier(float(s_a) / twice) + mb, eb = quantize_multiplier(float(s_b) / twice) + mo, eo = quantize_multiplier(twice / ((1 << ls) * s_out)) + A = (t[sa] - int(zp_a)).astype(np.int64) << ls + B = (t[sb] - int(zp_b)).astype(np.int64) << ls + y = requant(requant(A, ma, ea) + requant(B, mb, eb), mo, eo) + zp_out + t[q.output[0]] = np.clip(y, -128, 127).astype(np.int32) + if keep is not None: + keep.update(t) + res = [] + for o in g.outputs: + dq = g.prod[o] + if dq.input[0] in outs: + v = outs[dq.input[0]] + else: + s, zp = g.qparams[dq.input[0]] + v = (t[dq.input[0]] - zp) * s + res.append(v.transpose(2, 0, 1).astype(np.float32)) # C, H, W + return res diff --git a/tools/points.py b/tools/points.py new file mode 100644 index 0000000..cd836fc --- /dev/null +++ b/tools/points.py @@ -0,0 +1,100 @@ +# points.py +# +# What majestic's three object-detection presets mean for the person model: +# the score threshold of each, measured on the labelled crops (crops.py) from +# the engine's own detections (dets). A box is right when it overlaps a person +# by IoU 0.5; people marked ignore are neither a hit nor a miss. Matching is +# greedy by score, so one pass gives exact precision and recall at every +# threshold. +# +# balanced the F1-optimal threshold +# fewerFalseAlarms the F0.5 optimum (precision weighted) +# catchEverything the highest threshold with recall >= 0.90, else the +# F2 optimum +# +# Writes {"presets": {preset: {"human": threshold}}} for majestic to read +# beside the model, and a table of what each preset achieves, including the +# false alarms per person-free crop. +import json +import sys +from collections import defaultdict + +import numpy as np + +PRESETS = ["fewerFalseAlarms", "balanced", "catchEverything"] + + +def iou(a, b): + w = min(a[0] + a[2], b[0] + b[2]) - max(a[0], b[0]) + h = min(a[1] + a[3], b[1] + b[3]) - max(a[1], b[1]) + if w <= 0 or h <= 0: + return 0.0 + i = w * h + return i / (a[2] * a[3] + b[2] * b[3] - i) + + +def main(labels_path, dets_path, out_json, out_md): + labels = json.load(open(labels_path)) + dets = defaultdict(list) + for line in open(dets_path): + p = line.split() + if len(p) == 6: + dets[p[0]].append((int(p[5]), [float(v) for v in p[1:5]])) + npos = sum(1 for g in labels.values() for b in g if not b["ignore"]) + empty = [c for c, g in labels.items() if not g] + rec, bg = [], [] + for c, g in labels.items(): + used = [False] * len(g) + for s, box in sorted(dets.get(c, []), key=lambda d: -d[0]): + best, bi = 0.0, -1 + for i, b in enumerate(g): + if used[i] and not b["ignore"]: + continue + o = iou(box, b["box"]) + if o > best: + best, bi = o, i + if best >= 0.5: + if g[bi]["ignore"]: + continue + used[bi] = True + rec.append((s, True)) + else: + rec.append((s, False)) + if not g: + bg.append(s) + s = np.array([r[0] for r in rec]) + hit = np.array([r[1] for r in rec], bool) + bg = np.array(bg) + ts = np.arange(1, 100) + tp = np.array([(hit & (s >= t)).sum() for t in ts], float) + fp = np.array([(~hit & (s >= t)).sum() for t in ts], float) + some = tp + fp > 0 + p = np.divide(tp, tp + fp, out=np.zeros_like(tp), where=some) + r = tp / max(npos, 1) + f = lambda b: (1 + b * b) * p * r / np.maximum(b * b * p + r, 1e-9) + ok = np.nonzero(r >= 0.90)[0] + idx = {"fewerFalseAlarms": int(f(0.5).argmax()), "balanced": int(f(1).argmax()), + "catchEverything": int(ok[-1]) if len(ok) else int(f(2).argmax())} + presets, rows = {}, [] + for name in PRESETS: + i = idx[name] + t = ts[i] / 100 + presets[name] = {"human": round(float(t), 2)} + rows.append(f"| {name} | {t:.2f} | {p[i]:.3f} | {r[i]:.3f} | " + f"{2 * p[i] * r[i] / max(p[i] + r[i], 1e-9):.3f} | " + f"{(bg >= ts[i]).sum() / max(len(empty), 1):.3f} |") + with open(out_json, "w") as o: + json.dump({"presets": presets}, o, indent=1, sort_keys=True) + table = "\n".join([ + f"{len(labels)} crops of COCO val2017 ({len(empty)} with nobody in them), " + f"{npos} people counted.", + "", + "| preset | threshold | precision | recall | F1 | false alarms per person-free crop |", + "|---|---|---|---|---|---|"] + rows) + "\n" + with open(out_md, "w") as o: + o.write(table) + print(table) + + +if __name__ == "__main__": + main(*sys.argv[1:])